only implement if a use case is presented.
neg. binomial tests use the following parameters, which are currently hard-coded (to default) or not supported:

• latent.vars
• min.cells.feature
• ...

• Get input from the team.
• cap heatmaps like in enrichment analysis e.g.,
  • simple cap ie LFC>cap_th -> set to cap_th
  • log? (log again the LFC?)
  • scale the heatmaps?

Consider: it is a breaking change as it requires a new config field...

• add raw p-values, but still colored according to config on adjusted p-values
• indicate below that absolute LFC being used by e.g., |log2FC| > 1
• remove duplicate legend entries (probably fix by adding something like drop=FALSE during plotting)

like here: https://github.com/epigen/dea_limma/blob/4b31eaf3de1beb4199245084b48ee466955a4520/workflow/scripts/volcanos.R

• upgrade and test with Seurat version 4.4.0
• include this function if assay is SCT (can't hurt)
  • https://satijalab.org/seurat/reference/prepsctfindmarkers
  • motivated here

Hello,
thanks for the effort to turn the Seurat pipeline into a snakemake.
The following line in the Snakefile gives me an error since the function srcdir is not defined:

SDIR = os.path.realpath(os.path.dirname(srcdir("Snakefile")))

Had to comment it out to get to to run.
cheers

To get a better overview of the number of differential features in each direction(!), up and down, the barplot could go literally up and down from zero. "down" features need a negative sign (-).
plot_df are the DEA statistics.

# plot config
direction_col <- list("up"="red", "down"="blue")

# plot
ggplot(plot_df, aes(x=groups, y=n_features, fill=direction)) + 
        geom_bar(stat="identity", position="identity") +
        scale_fill_manual(values=direction_col) +
        custom_theme +
        theme(legend.position = "none",
             axis.text.x = element_text(angle = 45, vjust = 1, hjust=1)
             )

• {analysis}/plots/volcano/{feature_list}/{group}_{raw|adj}.png
• {analysis}/plots/heatmap/{feature_list}.png

feature_lists: ALL, ...

If only two cell groups are provided, there is only one column in the pheatmap plot (heatmap.R) and hclust causes an error.

This fixes it:

clusterCols <- ncol(dea_lfc) > 1

make heatmap

lfc_heatmap <- as.ggplot(pheatmap(dea_lfc,
show_rownames=F,
show_colnames=T,
cluster_cols = clusterCols,
fontsize = 5,
angle_col = 45,
treeheight_row = 25,
treeheight_col = 10,

annotation_col = annot,

           breaks=seq(-max(abs(dea_lfc)), max(abs(dea_lfc)), length.out=200),
           color=colorRampPalette(c("blue", "white", "red"))(200),
           annotation_names_col = F,
           silent = TRUE
          ))

• update current docs (work through current README)
  • rulegraph
  • remove patchwork from env.yamls
  • add package version numbers
• incorporate changes to latest features
  • restructuring of results and plots
• add latest features
  • feature lists
  • raw & adj p value volcano plots
  • statistics barplot
  • parallelization
  • example data set
    • provide rds object via Zenodo (3.2GB)
    • LFC heatmap of celltypemarkers
• update config.yaml

like here: https://github.com/epigen/dea_limma/blob/4b31eaf3de1beb4199245084b48ee466955a4520/config/config.yaml#L63

• configurable feature lists
• volcano plots for each feature list
• LFC heatmap for each feature list
• remove LFC matrix CSV generation (redundant, bc duplication of data)

only implement if a use case is presented.
Implement automatic DEA comparison of all pairs within a metadata label.
Annotation column control could be "PAIRWISE" or "PAIRS"

https://satijalab.org/seurat/archive/v3.0/future_vignette.html

setup Lee2020NatGenet as example data.

same functionality as implemented in dea_limma module.

this is important to be used eg for correction in enrichment analyses (eg as background gene set)