DE peaks/genes between two conditions (non-treated and treated cells)

[ccruizm]

Good day!

I have a related question to issue #215. When comparing two conditions (e.g., non-treated and treated cells) and after merging (or integration), should I use the RNA assay (using the inferred gene activities based on ATAC) or the ATAC (peaks) assay itself?

I have tried this:

DefaultAssay(combined) <- 'peaks'

da_peaks <- FindMarkers(
  object = combined,
  ident.1 = "treated",
    group.by = 'treatment',
    subset.ident = 'HSC',
  min.pct = 0.2,
  test.use = 'LR',
  latent.vars = 'peak_region_fragments'
)

and also following the 'Integrating stimulated vs. control PBMC datasets to learn cell-type-specific responses' vignette from RNA using the RNA assay from my ATAC data:

hsc <- subset(combined, idents = "HSC")
Idents(hsc) <- "treatment"
avg.hsc <- log1p(AverageExpression(hsc, verbose = FALSE)$RNA)
avg.hsc$gene <- rownames(avg.hsc)
but I am not sure which approach is the best for this case.

Thanks in advance for your help!

[timoast]

It depends if you want to find differentially accessible peaks or differentially accessible genes between the treatment groups. Usually we're interested in peaks and the gene activities are used mainly for annotating clusters, either manually by the accessibility of marker genes or automatically by integration with a matched annotated scRNA-seq dataset.

hsc <- subset(combined, idents = "HSC") Idents(hsc) <- "treatment" avg.hsc <- log1p(AverageExpression(hsc, verbose = FALSE)$RNA) avg.hsc$gene <- rownames(avg.hsc)

I don't think there's any value in running this for scATAC data. This was mainly a visualization used in the Seurat scRNA-seq immune alignment vignette.