Differential gene test

This function utilizes the Seurat package to perform a differential expression (DE) test on gene expression data. Users have the flexibility to specify custom cell groups, marker types, and various options for DE analysis.

Usage

RunDEtest(
  srt,
  group_by = NULL,
  group1 = NULL,
  group2 = NULL,
  cells1 = NULL,
  cells2 = NULL,
  features = NULL,
  markers_type = c("all", "paired", "conserved", "disturbed"),
  grouping.var = NULL,
  meta.method = c("maximump", "minimump", "wilkinsonp", "meanp", "sump", "votep"),
  test.use = "wilcox",
  only.pos = TRUE,
  fc.threshold = 1.5,
  base = 2,
  pseudocount.use = 1,
  mean.fxn = NULL,
  min.pct = 0.1,
  min.diff.pct = -Inf,
  max.cells.per.ident = Inf,
  latent.vars = NULL,
  min.cells.feature = 3,
  min.cells.group = 3,
  norm.method = "LogNormalize",
  p.adjust.method = "bonferroni",
  layer = "data",
  assay = NULL,
  seed = 11,
  verbose = TRUE,
  cores = 1,
  ...
)

Arguments

srt

A Seurat object.

group_by

A grouping variable in the dataset to define the groups or conditions for the differential test. If not provided, the function uses the "active.ident" variable in the Seurat object.

group1

A vector of cell IDs or a character vector specifying the cells that belong to the first group. If both group_by and group1 are provided, group1 takes precedence.

group2

A vector of cell IDs or a character vector specifying the cells that belong to the second group. This parameter is only used when group_by or group1 is provided.

cells1

A vector of cell IDs specifying the cells that belong to group1. If provided, group1 is ignored.

cells2

A vector of cell IDs specifying the cells that belong to group2. This parameter is only used when cells1 is provided.

features

A vector of gene names specifying the features to consider for the differential test. If not provided, all features in the dataset are considered.

markers_type

A character value specifying the type of markers to find. Possible values are "all", "paired", "conserved", and "disturbed".

grouping.var

A character value specifying the grouping variable for finding conserved or disturbed markers. This parameter is only used when markers_type is "conserved" or "disturbed".

meta.method

A character value specifying the method to use for combining p-values in the conserved markers test. Possible values are "maximump", "minimump", "wilkinsonp", "meanp", "sump", and "votep".

test.use

Denotes which test to use. Available options are:

"wilcox" : Identifies differentially expressed genes between two groups of cells using a Wilcoxon Rank Sum test (default); will use a fast implementation by Presto if installed
"wilcox_limma" : Identifies differentially expressed genes between two groups of cells using the limma implementation of the Wilcoxon Rank Sum test; set this option to reproduce results from Seurat v4
"bimod" : Likelihood-ratio test for single cell gene expression, (McDavid et al., Bioinformatics, 2013)
"roc" : Identifies 'markers' of gene expression using ROC analysis. For each gene, evaluates (using AUC) a classifier built on that gene alone, to classify between two groups of cells. An AUC value of 1 means that expression values for this gene alone can perfectly classify the two groupings (i.e. Each of the cells in cells.1 exhibit a higher level than each of the cells in cells.2). An AUC value of 0 also means there is perfect classification, but in the other direction. A value of 0.5 implies that the gene has no predictive power to classify the two groups. Returns a 'predictive power' (abs(AUC-0.5) * 2) ranked matrix of putative differentially expressed genes.
"t" : Identify differentially expressed genes between two groups of cells using the Student's t-test.
"negbinom" : Identifies differentially expressed genes between two groups of cells using a negative binomial generalized linear model. Use only for UMI-based datasets
"poisson" : Identifies differentially expressed genes between two groups of cells using a poisson generalized linear model. Use only for UMI-based datasets
"LR" : Uses a logistic regression framework to determine differentially expressed genes. Constructs a logistic regression model predicting group membership based on each feature individually and compares this to a null model with a likelihood ratio test.
"MAST" : Identifies differentially expressed genes between two groups of cells using a hurdle model tailored to scRNA-seq data. Utilizes the MAST package to run the DE testing.
"DESeq2" : Identifies differentially expressed genes between two groups of cells based on a model using DESeq2 which uses a negative binomial distribution (Love et al, Genome Biology, 2014).This test does not support pre-filtering of genes based on average difference (or percent detection rate) between cell groups. However, genes may be pre-filtered based on their minimum detection rate (min.pct) across both cell groups. To use this method, please install DESeq2, using the instructions at https://bioconductor.org/packages/release/bioc/html/DESeq2.html

only.pos

Only return positive markers (FALSE by default)

fc.threshold

A numeric value used to filter genes for testing based on their average fold change between/among the two groups. Default is 1.5.

base

The base with respect to which logarithms are computed.

pseudocount.use

Pseudocount to add to averaged expression values when calculating logFC. 1 by default.

mean.fxn

Function to use for fold change or average difference calculation. The default depends on the the value of fc.slot:

"counts" : difference in the log of the mean counts, with pseudocount.
"data" : difference in the log of the average exponentiated data, with pseudocount. This adjusts for differences in sequencing depth between cells, and assumes that "data" has been log-normalized.
"scale.data" : difference in the means of scale.data.

min.pct

only test genes that are detected in a minimum fraction of min.pct cells in either of the two populations. Meant to speed up the function by not testing genes that are very infrequently expressed. Default is 0.01

min.diff.pct

only test genes that show a minimum difference in the fraction of detection between the two groups. Set to -Inf by default

max.cells.per.ident

Down sample each identity class to a max number. Default is no downsampling. Not activated by default (set to Inf)

latent.vars

Variables to test, used only when test.use is one of 'LR', 'negbinom', 'poisson', or 'MAST'

min.cells.feature

Minimum number of cells expressing the feature in at least one of the two groups, currently only used for poisson and negative binomial tests

min.cells.group

Minimum number of cells in one of the groups

norm.method

Normalization method for fold change calculation when layer is 'data'. Default is "LogNormalize".

p.adjust.method

A character value specifying the method to use for adjusting p-values. Default is "bonferroni".

layer

The layer used.

assay

Assay to use in differential expression testing

seed

An integer value specifying the seed. Default is 11.

verbose

Print a progress bar once expression testing begins

cores

The number of cores to use for parallelization with foreach::foreach. Default is 1.

...

Additional arguments to pass to the Seurat::FindMarkers function.

Examples

data(pancreas_sub)
pancreas_sub <- standard_scop(pancreas_sub)
#> ℹ [2025-09-20 13:33:51] Start standard scop workflow...
#> ℹ [2025-09-20 13:33:52] Checking a list of <Seurat> object...
#> ! [2025-09-20 13:33:52] Data 1/1 of the `srt_list` is "unknown"
#> ℹ [2025-09-20 13:33:52] Perform `NormalizeData()` with `normalization.method = 'LogNormalize'` on the data 1/1 of the `srt_list`...
#> ℹ [2025-09-20 13:33:54] Perform `Seurat::FindVariableFeatures()` on the data 1/1 of the `srt_list`...
#> ℹ [2025-09-20 13:33:55] Use the separate HVF from srt_list
#> ℹ [2025-09-20 13:33:55] Number of available HVF: 2000
#> ℹ [2025-09-20 13:33:55] Finished check
#> ℹ [2025-09-20 13:33:55] Perform `Seurat::ScaleData()`
#> Warning: Different features in new layer data than already exists for scale.data
#> ℹ [2025-09-20 13:33:55] Perform pca linear dimension reduction
#> StandardPC_ 1 
#> Positive:  Aplp1, Cpe, Gnas, Fam183b, Map1b, Hmgn3, Pcsk1n, Chga, Tuba1a, Bex2 
#> 	   Syt13, Isl1, 1700086L19Rik, Pax6, Chgb, Scgn, Rbp4, Scg3, Gch1, Camk2n1 
#> 	   Cryba2, Pcsk2, Pyy, Tspan7, Mafb, Hist3h2ba, Dbpht2, Abcc8, Rap1b, Slc38a5 
#> Negative:  Spp1, Anxa2, Sparc, Dbi, 1700011H14Rik, Wfdc2, Gsta3, Adamts1, Clu, Mgst1 
#> 	   Bicc1, Ldha, Vim, Cldn3, Cyr61, Rps2, Mt1, Ptn, Phgdh, Nudt19 
#> 	   Smtnl2, Smco4, Habp2, Mt2, Col18a1, Rpl12, Galk1, Cldn10, Acot1, Ccnd1 
#> StandardPC_ 2 
#> Positive:  Rbp4, Tagln2, Tuba1b, Fkbp2, Pyy, Pcsk2, Iapp, Tmem27, Meis2, Tubb4b 
#> 	   Pcsk1n, Dbpht2, Rap1b, Dynll1, Tubb2a, Sdf2l1, Scgn, 1700086L19Rik, Scg2, Abcc8 
#> 	   Atp1b1, Hspa5, Fam183b, Papss2, Slc38a5, Scg3, Mageh1, Tspan7, Ppp1r1a, Ociad2 
#> Negative:  Neurog3, Btbd17, Gadd45a, Ppp1r14a, Neurod2, Sox4, Smarcd2, Mdk, Pax4, Btg2 
#> 	   Sult2b1, Hes6, Grasp, Igfbpl1, Gpx2, Cbfa2t3, Foxa3, Shf, Mfng, Tmsb4x 
#> 	   Amotl2, Gdpd1, Cdc14b, Epb42, Rcor2, Cotl1, Upk3bl, Rbfox3, Cldn6, Cer1 
#> StandardPC_ 3 
#> Positive:  Nusap1, Top2a, Birc5, Aurkb, Cdca8, Pbk, Mki67, Tpx2, Plk1, Ccnb1 
#> 	   2810417H13Rik, Incenp, Cenpf, Ccna2, Prc1, Racgap1, Cdk1, Aurka, Cdca3, Hmmr 
#> 	   Spc24, Kif23, Sgol1, Cenpe, Cdc20, Hist1h1b, Cdca2, Mxd3, Kif22, Ska1 
#> Negative:  Anxa5, Pdzk1ip1, Acot1, Tpm1, Anxa2, Dcdc2a, Capg, Sparc, Ttr, Pamr1 
#> 	   Clu, Cxcl12, Ndrg2, Hnf1aos1, Gas6, Gsta3, Krt18, Ces1d, Atp1b1, Muc1 
#> 	   Hhex, Acadm, Spp1, Enpp2, Bcl2l14, Sat1, Smtnl2, 1700011H14Rik, Tgm2, Fam159a 
#> StandardPC_ 4 
#> Positive:  Glud1, Tm4sf4, Akr1c19, Cldn4, Runx1t1, Fev, Pou3f4, Gm43861, Pgrmc1, Arx 
#> 	   Cd200, Lrpprc, Hmgn3, Ppp1r14c, Pam, Etv1, Tsc22d1, Slc25a5, Akap17b, Pgf 
#> 	   Fam43a, Emb, Jun, Krt8, Dnajc12, Mid1ip1, Ids, Rgs17, Uchl1, Alcam 
#> Negative:  Ins2, Ins1, Ppp1r1a, Nnat, Calr, Sytl4, Sdf2l1, Iapp, Pdia6, Mapt 
#> 	   G6pc2, C2cd4b, Npy, Gng12, P2ry1, Ero1lb, Adra2a, Papss2, Arhgap36, Fam151a 
#> 	   Dlk1, Creld2, Gip, Tmem215, Gm27033, Cntfr, Prss53, C2cd4a, Lyve1, Ociad2 
#> StandardPC_ 5 
#> Positive:  Pdx1, Nkx6-1, Npepl1, Cldn4, Cryba2, Fev, Jun, Chgb, Gng12, Adra2a 
#> 	   Mnx1, Sytl4, Pdk3, Gm27033, Nnat, Chga, Ins2, 1110012L19Rik, Enho, Krt7 
#> 	   Mlxipl, Tmsb10, Flrt1, Pax4, Tubb3, Prrg2, Gars, Frzb, BC023829, Gm2694 
#> Negative:  Irx2, Irx1, Gcg, Ctxn2, Tmem27, Ctsz, Tmsb15l, Nap1l5, Pou6f2, Gria2 
#> 	   Ghrl, Peg10, Smarca1, Arx, Lrpap1, Rgs4, Ttr, Gast, Tmsb15b2, Serpina1b 
#> 	   Slc16a10, Wnk3, Ly6e, Auts2, Sct, Arg1, Dusp10, Sphkap, Dock11, Edn3 
#> ℹ [2025-09-20 13:33:56] Perform `Seurat::FindClusters()` with louvain and `cluster_resolution` = 0.6
#> ℹ [2025-09-20 13:33:57] Reorder clusters...
#> ! [2025-09-20 13:33:57] Using `Seurat::AggregateExpression()` to calculate pseudo-bulk data for <Assay5>
#> ℹ [2025-09-20 13:33:57] Perform umap nonlinear dimension reduction
#> ℹ [2025-09-20 13:33:57] Non-linear dimensionality reduction (umap) using (Standardpca) dims (1-50) as input
#> ℹ [2025-09-20 13:33:57] UMAP will return its model
#> ℹ [2025-09-20 13:34:01] Non-linear dimensionality reduction (umap) using (Standardpca) dims (1-50) as input
#> ℹ [2025-09-20 13:34:01] UMAP will return its model
#> ✔ [2025-09-20 13:34:05] Run scop standard workflow done
pancreas_sub <- RunDEtest(
  pancreas_sub,
  group_by = "SubCellType"
)
#> ℹ [2025-09-20 13:34:05] immunogenomics/presto installed successfully
#> ℹ [2025-09-20 13:34:06] Data type is log-normalized
#> ℹ [2025-09-20 13:34:06] Start differential expression test
#> ℹ [2025-09-20 13:34:06] Find all markers(wilcox) among 8 groups...
#> ℹ [2025-09-20 13:34:06] Using 1 core
#> ⠙ [2025-09-20 13:34:06] Running [1/8] ETA:  1s
#> ⠹ [2025-09-20 13:34:06] Running [4/8] ETA:  1s
#> ✔ [2025-09-20 13:34:06] Completed 8 tasks in 1.3s
#> 
#> ℹ [2025-09-20 13:34:07] Building results
#> ✔ [2025-09-20 13:34:07] Differential expression test completed
AllMarkers <- dplyr::filter(
  pancreas_sub@tools$DEtest_SubCellType$AllMarkers_wilcox,
  p_val_adj < 0.05 & avg_log2FC > 1
)
table(AllMarkers$group1)
#> 
#>        Ductal  Ngn3-high-EP          Beta   Ngn3-low-EP Pre-endocrine 
#>          1050           366           459           336           225 
#>         Alpha       Epsilon         Delta 
#>           272           169           148 
ht1 <- GroupHeatmap(
  pancreas_sub,
  features = AllMarkers$gene,
  feature_split = AllMarkers$group1,
  group.by = "SubCellType"
)
#> 'magick' package is suggested to install to give better rasterization.
#> 
#> Set `ht_opt$message = FALSE` to turn off this message.
#> ℹ [2025-09-20 13:34:10] The size of the heatmap is fixed because certain elements are not scalable.
#> ℹ [2025-09-20 13:34:10] The width and height of the heatmap are determined by the size of the current viewport.
#> ℹ [2025-09-20 13:34:10] If you want to have more control over the size, you can manually set the parameters 'width' and 'height'.

ht1$plot


TopMarkers <- AllMarkers |>
  dplyr::group_by(gene) |>
  dplyr::top_n(1, avg_log2FC) |>
  dplyr::group_by(group1) |>
  dplyr::top_n(3, avg_log2FC)
ht2 <- GroupHeatmap(
  pancreas_sub,
  features = TopMarkers$gene,
  feature_split = TopMarkers$group1,
  group.by = "SubCellType",
  show_row_names = TRUE
)
#> 'magick' package is suggested to install to give better rasterization.
#> 
#> Set `ht_opt$message = FALSE` to turn off this message.
#> ℹ [2025-09-20 13:34:13] The size of the heatmap is fixed because certain elements are not scalable.
#> ℹ [2025-09-20 13:34:13] The width and height of the heatmap are determined by the size of the current viewport.
#> ℹ [2025-09-20 13:34:13] If you want to have more control over the size, you can manually set the parameters 'width' and 'height'.

ht2$plot


pancreas_sub <- RunDEtest(
  pancreas_sub,
  group_by = "SubCellType",
  markers_type = "paired"
)
#> ℹ [2025-09-20 13:34:13] immunogenomics/presto installed successfully
#> ℹ [2025-09-20 13:34:14] Data type is log-normalized
#> ℹ [2025-09-20 13:34:14] Start differential expression test
#> ℹ [2025-09-20 13:34:14] Find paired markers(wilcox) among 8 groups...
#> ℹ [2025-09-20 13:34:14] Using 1 core
#> ⠙ [2025-09-20 13:34:14] Running [1/56] ETA:  6s
#> ⠹ [2025-09-20 13:34:14] Running [13/56] ETA:  5s
#> ⠸ [2025-09-20 13:34:14] Running [44/56] ETA:  1s
#> ✔ [2025-09-20 13:34:14] Completed 56 tasks in 5.4s
#> 
#> ℹ [2025-09-20 13:34:19] Building results
#> ✔ [2025-09-20 13:34:19] Differential expression test completed
PairedMarkers <- dplyr::filter(
  pancreas_sub@tools$DEtest_SubCellType$PairedMarkers_wilcox,
  p_val_adj < 0.05 & avg_log2FC > 1
)
ht3 <- GroupHeatmap(
  pancreas_sub,
  features = PairedMarkers$gene,
  feature_split = PairedMarkers$group1,
  group.by = "SubCellType"
)
#> 'magick' package is suggested to install to give better rasterization.
#> 
#> Set `ht_opt$message = FALSE` to turn off this message.
#> ℹ [2025-09-20 13:35:05] The size of the heatmap is fixed because certain elements are not scalable.
#> ℹ [2025-09-20 13:35:05] The width and height of the heatmap are determined by the size of the current viewport.
#> ℹ [2025-09-20 13:35:05] If you want to have more control over the size, you can manually set the parameters 'width' and 'height'.

ht3$plot


data(panc8_sub)
panc8_sub <- integration_scop(
  panc8_sub,
  batch = "tech",
  integration_method = "Seurat"
)
#> ℹ [2025-09-20 13:35:08] Run Seurat integration...
#> ℹ [2025-09-20 13:35:08] Spliting `srt_merge` into `srt_list` by column "tech"...
#> ℹ [2025-09-20 13:35:09] Checking a list of <Seurat> object...
#> ! [2025-09-20 13:35:09] Data 1/5 of the `srt_list` is "unknown"
#> ℹ [2025-09-20 13:35:09] Perform `NormalizeData()` with `normalization.method = 'LogNormalize'` on the data 1/5 of the `srt_list`...
#> ℹ [2025-09-20 13:35:10] Perform `Seurat::FindVariableFeatures()` on the data 1/5 of the `srt_list`...
#> ! [2025-09-20 13:35:11] Data 2/5 of the `srt_list` is "unknown"
#> ℹ [2025-09-20 13:35:11] Perform `NormalizeData()` with `normalization.method = 'LogNormalize'` on the data 2/5 of the `srt_list`...
#> ℹ [2025-09-20 13:35:12] Perform `Seurat::FindVariableFeatures()` on the data 2/5 of the `srt_list`...
#> ! [2025-09-20 13:35:13] Data 3/5 of the `srt_list` is "unknown"
#> ℹ [2025-09-20 13:35:13] Perform `NormalizeData()` with `normalization.method = 'LogNormalize'` on the data 3/5 of the `srt_list`...
#> ℹ [2025-09-20 13:35:14] Perform `Seurat::FindVariableFeatures()` on the data 3/5 of the `srt_list`...
#> ! [2025-09-20 13:35:15] Data 4/5 of the `srt_list` is "unknown"
#> ℹ [2025-09-20 13:35:15] Perform `NormalizeData()` with `normalization.method = 'LogNormalize'` on the data 4/5 of the `srt_list`...
#> ℹ [2025-09-20 13:35:16] Perform `Seurat::FindVariableFeatures()` on the data 4/5 of the `srt_list`...
#> ! [2025-09-20 13:35:17] Data 5/5 of the `srt_list` is "unknown"
#> ℹ [2025-09-20 13:35:17] Perform `NormalizeData()` with `normalization.method = 'LogNormalize'` on the data 5/5 of the `srt_list`...
#> ℹ [2025-09-20 13:35:18] Perform `Seurat::FindVariableFeatures()` on the data 5/5 of the `srt_list`...
#> ℹ [2025-09-20 13:35:19] Use the separate HVF from srt_list
#> ℹ [2025-09-20 13:35:19] Number of available HVF: 2000
#> ℹ [2025-09-20 13:35:20] Finished check
#> ℹ [2025-09-20 13:35:22] Perform FindIntegrationAnchors
#> Warning: Different features in new layer data than already exists for scale.data
#> Warning: Different features in new layer data than already exists for scale.data
#> Warning: Different features in new layer data than already exists for scale.data
#> Warning: Different features in new layer data than already exists for scale.data
#> Warning: Different features in new layer data than already exists for scale.data
#> ℹ [2025-09-20 13:35:53] Perform integration(Seurat)
#> Warning: Layer counts isn't present in the assay object; returning NULL
#> Warning: Different cells in new layer data than already exists for scale.data
#> Warning: Layer counts isn't present in the assay object; returning NULL
#> Warning: Different cells in new layer data than already exists for scale.data
#> Warning: Layer counts isn't present in the assay object; returning NULL
#> Warning: Different cells in new layer data than already exists for scale.data
#> Warning: Layer counts isn't present in the assay object; returning NULL
#> ℹ [2025-09-20 13:36:05] Perform ScaleData
#> ℹ [2025-09-20 13:36:06] Perform linear dimension reduction ("pca")
#> ℹ [2025-09-20 13:36:07] Perform FindClusters ("louvain")
#> ℹ [2025-09-20 13:36:07] Reorder clusters...
#> ! [2025-09-20 13:36:07] Using `Seurat::AverageExpression()` to calculate pseudo-bulk data for <Assay>
#> ℹ [2025-09-20 13:36:07] Perform nonlinear dimension reduction (umap)
#> ℹ [2025-09-20 13:36:07] Non-linear dimensionality reduction (umap) using (Seuratpca) dims (1-12) as input
#> ℹ [2025-09-20 13:36:12] Non-linear dimensionality reduction (umap) using (Seuratpca) dims (1-12) as input
#> ✔ [2025-09-20 13:36:18] Run Seurat integration done
CellDimPlot(
  panc8_sub,
  group.by = c("celltype", "tech")
)
#> Warning: No shared levels found between `names(values)` of the manual scale and the
#> data's fill values.
#> Warning: No shared levels found between `names(values)` of the manual scale and the
#> data's fill values.


panc8_sub <- RunDEtest(
  srt = panc8_sub,
  group_by = "celltype",
  grouping.var = "tech",
  markers_type = "conserved"
)
#> ℹ [2025-09-20 13:36:19] immunogenomics/presto installed successfully
#> ℹ [2025-09-20 13:36:19] Data type is log-normalized
#> ℹ [2025-09-20 13:36:19] Start differential expression test
#> ℹ [2025-09-20 13:36:19] Find conserved markers(wilcox) among 13 groups...
#> ℹ [2025-09-20 13:36:19] Using 1 core
#> ⠙ [2025-09-20 13:36:19] Running [1/13] ETA: 10s
#> ⠹ [2025-09-20 13:36:19] Running [3/13] ETA:  9s
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s
#> ! [2025-09-20 13:36:25] Fluidigmc1 has fewer than 3 cells. Skipping fluidigmc1
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:25] Celseq has fewer than 3 cells. Skipping celseq
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:25] Celseq2 has fewer than 3 cells. Skipping celseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:25] Fluidigmc1 has fewer than 3 cells. Skipping fluidigmc1
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:25] Smartseq2 has fewer than 3 cells. Skipping smartseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Only a single group was tested
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Celseq has fewer than 3 cells. Skipping celseq
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Celseq2 has fewer than 3 cells. Skipping celseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Fluidigmc1 has fewer than 3 cells. Skipping fluidigmc1
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Smartseq2 has fewer than 3 cells. Skipping smartseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Only a single group was tested
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Celseq has fewer than 3 cells. Skipping celseq
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Celseq2 has fewer than 3 cells. Skipping celseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Fluidigmc1 has fewer than 3 cells. Skipping fluidigmc1
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Smartseq2 has fewer than 3 cells. Skipping smartseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Only a single group was tested
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:26] Celseq has fewer than 3 cells. Skipping celseq
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ! [2025-09-20 13:36:27] Smartseq2 has fewer than 3 cells. Skipping smartseq2
#> ⠸ [2025-09-20 13:36:19] Running [6/13] ETA:  6s

#> ✔ [2025-09-20 13:36:19] Completed 13 tasks in 7.7s
#> 
#> ℹ [2025-09-20 13:36:27] Building results
#> Error in names(output_list) <- x: 'names' attribute [13] must be the same length as the vector [12]
ConservedMarkers1 <- dplyr::filter(
  panc8_sub@tools$DEtest_celltype$ConservedMarkers_wilcox,
  p_val_adj < 0.05 & avg_log2FC > 1
)
#> Error in UseMethod("filter"): no applicable method for 'filter' applied to an object of class "NULL"
ht4 <- GroupHeatmap(
  panc8_sub,
  layer = "data",
  features = ConservedMarkers1$gene,
  feature_split = ConservedMarkers1$group1,
  group.by = "tech",
  split.by = "celltype",
  within_groups = TRUE
)
#> Error: object 'ConservedMarkers1' not found
ht4$plot
#> Error: object 'ht4' not found

panc8_sub <- RunDEtest(
  srt = panc8_sub,
  group_by = "tech",
  grouping.var = "celltype",
  markers_type = "conserved"
)
#> ℹ [2025-09-20 13:36:27] immunogenomics/presto installed successfully
#> ℹ [2025-09-20 13:36:28] Data type is log-normalized
#> ℹ [2025-09-20 13:36:28] Start differential expression test
#> ℹ [2025-09-20 13:36:28] Find conserved markers(wilcox) among 5 groups...
#> ℹ [2025-09-20 13:36:28] Using 1 core
#> ! [2025-09-20 13:36:28] Schwann has fewer than 3 cells. Skipping schwann
#> ! [2025-09-20 13:36:28] Epsilon has fewer than 3 cells. Skipping epsilon
#> ! [2025-09-20 13:36:28] Mast has fewer than 3 cells. Skipping mast
#> ! [2025-09-20 13:36:28] Macrophage has fewer than 3 cells. Skipping macrophage
#> ! [2025-09-20 13:36:28] Quiescent-stellate has fewer than 3 cells. Skipping quiescent-stellate
#> ! [2025-09-20 13:36:28] Endothelial has fewer than 3 cells. Skipping endothelial
#> ⠙ [2025-09-20 13:36:28] Running [1/5] ETA:  4s
#> ! [2025-09-20 13:36:29] Schwann has fewer than 3 cells. Skipping schwann
#> ⠙ [2025-09-20 13:36:28] Running [1/5] ETA:  4s

#> ! [2025-09-20 13:36:29] Epsilon has fewer than 3 cells. Skipping epsilon
#> ⠙ [2025-09-20 13:36:28] Running [1/5] ETA:  4s

#> ! [2025-09-20 13:36:29] Mast has fewer than 3 cells. Skipping mast
#> ⠙ [2025-09-20 13:36:28] Running [1/5] ETA:  4s

#> ! [2025-09-20 13:36:29] Macrophage has fewer than 3 cells. Skipping macrophage
#> ⠙ [2025-09-20 13:36:28] Running [1/5] ETA:  4s

#> ! [2025-09-20 13:36:29] Quiescent-stellate has fewer than 3 cells. Skipping quiescent-stellate
#> ⠙ [2025-09-20 13:36:28] Running [1/5] ETA:  4s

#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s
#> ! [2025-09-20 13:36:30] Schwann has fewer than 3 cells. Skipping schwann
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:30] Epsilon has fewer than 3 cells. Skipping epsilon
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:30] Mast has fewer than 3 cells. Skipping mast
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:30] Macrophage has fewer than 3 cells. Skipping macrophage
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:30] Quiescent-stellate has fewer than 3 cells. Skipping quiescent-stellate
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:30] Endothelial has fewer than 3 cells. Skipping endothelial
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:31] Schwann has fewer than 3 cells. Skipping schwann
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:31] Epsilon has fewer than 3 cells. Skipping epsilon
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:31] Mast has fewer than 3 cells. Skipping mast
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:31] Macrophage has fewer than 3 cells. Skipping macrophage
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:31] Quiescent-stellate has fewer than 3 cells. Skipping quiescent-stellate
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:32] Activated-stellate has fewer than 3 cells. Skipping activated-stellate
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:33] Schwann has fewer than 3 cells. Skipping schwann
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:33] Epsilon has fewer than 3 cells. Skipping epsilon
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:33] Mast has fewer than 3 cells. Skipping mast
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ! [2025-09-20 13:36:33] Quiescent-stellate has fewer than 3 cells. Skipping quiescent-stellate
#> ⠹ [2025-09-20 13:36:28] Running [2/5] ETA:  4s

#> ✔ [2025-09-20 13:36:28] Completed 5 tasks in 5.9s
#> 
#> ℹ [2025-09-20 13:36:34] Building results
#> ✔ [2025-09-20 13:36:34] Differential expression test completed
ConservedMarkers2 <- dplyr::filter(
  panc8_sub@tools$DEtest_tech$ConservedMarkers_wilcox,
  p_val_adj < 0.05 & avg_log2FC > 1
)
table(ConservedMarkers2$group1)
#> 
#>     celseq    celseq2  smartseq2 fluidigmc1     indrop 
#>         89        319        266        890         32 
ht4 <- GroupHeatmap(
  srt = panc8_sub,
  layer = "data",
  features = ConservedMarkers2$gene,
  feature_split = ConservedMarkers2$group1,
  group.by = "tech",
  split.by = "celltype"
)
#> 'magick' package is suggested to install to give better rasterization.
#> 
#> Set `ht_opt$message = FALSE` to turn off this message.
#> ℹ [2025-09-20 13:36:36] The size of the heatmap is fixed because certain elements are not scalable.
#> ℹ [2025-09-20 13:36:36] The width and height of the heatmap are determined by the size of the current viewport.
#> ℹ [2025-09-20 13:36:36] If you want to have more control over the size, you can manually set the parameters 'width' and 'height'.

ht4$plot


panc8_sub <- RunDEtest(
  srt = panc8_sub,
  group_by = "celltype",
  grouping.var = "tech",
  markers_type = "disturbed"
)
#> ℹ [2025-09-20 13:36:42] immunogenomics/presto installed successfully
#> ℹ [2025-09-20 13:36:42] Data type is log-normalized
#> ℹ [2025-09-20 13:36:42] Start differential expression test
#> ℹ [2025-09-20 13:36:42] Find disturbed markers(wilcox) among 13 groups...
#> ℹ [2025-09-20 13:36:42] Using 1 core
#> ℹ [2025-09-20 13:36:42] immunogenomics/presto installed successfully
#> ℹ [2025-09-20 13:36:43] Data type is log-normalized
#> ℹ [2025-09-20 13:36:43] Find all markers(wilcox) among 5 groups...
#> ⠙ [2025-09-20 13:36:42] Running [1/13] ETA: 18s
#> ℹ [2025-09-20 13:36:44] immunogenomics/presto installed successfully
#> ⠙ [2025-09-20 13:36:42] Running [1/13] ETA: 18s

#> ℹ [2025-09-20 13:36:45] Data type is log-normalized
#> ⠙ [2025-09-20 13:36:42] Running [1/13] ETA: 18s

#> ℹ [2025-09-20 13:36:45] Find all markers(wilcox) among 5 groups...
#> ⠙ [2025-09-20 13:36:42] Running [1/13] ETA: 18s

#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s
#> ℹ [2025-09-20 13:36:45] immunogenomics/presto installed successfully
#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s

#> ℹ [2025-09-20 13:36:46] Data type is log-normalized
#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s

#> ℹ [2025-09-20 13:36:46] Find all markers(wilcox) among 5 groups...
#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s

#> ℹ [2025-09-20 13:36:47] immunogenomics/presto installed successfully
#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s

#> ℹ [2025-09-20 13:36:48] Data type is log-normalized
#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s

#> ℹ [2025-09-20 13:36:48] Find all markers(wilcox) among 5 groups...
#> ⠹ [2025-09-20 13:36:42] Running [2/13] ETA: 16s

#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s
#> ℹ [2025-09-20 13:36:49] immunogenomics/presto installed successfully
#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s

#> ℹ [2025-09-20 13:36:49] Data type is log-normalized
#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s

#> ℹ [2025-09-20 13:36:49] Find all markers(wilcox) among 5 groups...
#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s

#> ℹ [2025-09-20 13:36:50] immunogenomics/presto installed successfully
#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s

#> ℹ [2025-09-20 13:36:51] Data type is log-normalized
#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s

#> ℹ [2025-09-20 13:36:51] Find all markers(wilcox) among 5 groups...
#> ⠸ [2025-09-20 13:36:42] Running [4/13] ETA: 14s

#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s
#> ℹ [2025-09-20 13:36:52] immunogenomics/presto installed successfully
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:53] Data type is log-normalized
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:53] Find all markers(wilcox) among 5 groups...
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:53] immunogenomics/presto installed successfully
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:54] Data type is log-normalized
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:54] Find all markers(wilcox) among 3 groups...
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:54] immunogenomics/presto installed successfully
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:55] Data type is log-normalized
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ℹ [2025-09-20 13:36:55] Find all markers(wilcox) among 4 groups...
#> ⠼ [2025-09-20 13:36:42] Running [6/13] ETA: 11s

#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s
#> ℹ [2025-09-20 13:36:55] immunogenomics/presto installed successfully
#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s

#> ℹ [2025-09-20 13:36:56] Data type is log-normalized
#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s

#> ℹ [2025-09-20 13:36:56] Find all markers(wilcox) among 2 groups...
#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s

#> ℹ [2025-09-20 13:36:56] immunogenomics/presto installed successfully
#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s

#> ℹ [2025-09-20 13:36:57] Data type is log-normalized
#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s

#> ℹ [2025-09-20 13:36:57] Find all markers(wilcox) among 4 groups...
#> ⠴ [2025-09-20 13:36:42] Running [9/13] ETA:  6s

#> ✔ [2025-09-20 13:36:42] Completed 13 tasks in 14.7s
#> 
#> ℹ [2025-09-20 13:36:57] Building results
#> Error in names(output_list) <- x: 'names' attribute [13] must be the same length as the vector [11]
DisturbedMarkers <- dplyr::filter(
  panc8_sub@tools$DEtest_celltype$DisturbedMarkers_wilcox,
  p_val_adj < 0.05 & avg_log2FC > 1 & var1 == "smartseq2"
)
#> Error in UseMethod("filter"): no applicable method for 'filter' applied to an object of class "NULL"
table(DisturbedMarkers$group1)
#> Error: object 'DisturbedMarkers' not found
ht5 <- GroupHeatmap(
  srt = panc8_sub,
  layer = "data",
  features = DisturbedMarkers$gene,
  feature_split = DisturbedMarkers$group1,
  group.by = "celltype",
  split.by = "tech"
)
#> Error: object 'DisturbedMarkers' not found
ht5$plot
#> Error: object 'ht5' not found

gene_specific <- names(which(table(DisturbedMarkers$gene) == 1))
#> Error in h(simpleError(msg, call)): error in evaluating the argument 'x' in selecting a method for function 'which': object 'DisturbedMarkers' not found
DisturbedMarkers_specific <- DisturbedMarkers[
  DisturbedMarkers$gene %in% gene_specific,
]
#> Error: object 'DisturbedMarkers' not found
table(DisturbedMarkers_specific$group1)
#> Error: object 'DisturbedMarkers_specific' not found
ht6 <- GroupHeatmap(
  srt = panc8_sub,
  layer = "data",
  features = DisturbedMarkers_specific$gene,
  feature_split = DisturbedMarkers_specific$group1,
  group.by = "celltype",
  split.by = "tech"
)
#> Error: object 'DisturbedMarkers_specific' not found
ht6$plot
#> Error: object 'ht6' not found

ht7 <- GroupHeatmap(
  srt = panc8_sub,
  layer = "data",
  aggregate_fun = function(x) mean(expm1(x)) + 1,
  features = DisturbedMarkers_specific$gene,
  feature_split = DisturbedMarkers_specific$group1,
  group.by = "celltype",
  grouping.var = "tech",
  numerator = "smartseq2"
)
#> Error: object 'DisturbedMarkers_specific' not found
ht7$plot
#> Error: object 'ht7' not found

Usage

Arguments

See also

Examples