satijalab / seurat

R toolkit for single cell genomics
http://www.satijalab.org/seurat
Other
2.24k stars 902 forks source link

What slot does DEnrichRPlot use within the assay? #9141

Closed cristanchoa closed 21 hours ago

cristanchoa commented 1 month ago

I was running DEnrichRPlot using my SCT assay after the PrepSCTFindMarkers. The data didn't quite make sense, so I reran it under the "RNA" assay. The data was drastically different (honestly made more sense). Since I couldn't quite tell from the code, which slot was used and the "RNA" data hadn't been scaled since the object was integrated with SCT, I scaled the RNA assay. Now I reran DEnrichRPlot in RNA and the results are different again (and very similar to the SCT results). To clarify, the reason don't quite quite make much sense is that markers for a totally different cell type and otherwise sparsely expressed genes are dominating the analysis (verified by Violin and feature plots). Also, it is not the "NormalizeData" that changes the DEnrichRPlot output, it is just when I "ScaleData."

In FindMarkers the default assay is "data" and not "scale.data" (and documented nicely elsewhere the reasoning). Is DEnrichRPlot using "scale.data" and if so is there a way to change it?

samuel-marsh commented 1 month ago

Hi,

Not member of dev team but hopefully can be helpful. DEenrichRPlot uses the default layer for whatever DE test is being used through FindMarkers. See source code here:

https://github.com/satijalab/seurat/blob/1549dcb3075eaeac01c925c4b4bb73c73450fc50/R/mixscape.R#L186-L194

Best, Sam

cristanchoa commented 1 month ago

Hi Sam, Thanks for answering. Then perhaps there is a bug of some kind somewhere, because the data does change with the addition of the "scale.data" for the assay which does not appear like it should based on the FindMarkers() default (https://satijalab.org/seurat/reference/findmarkers).

Although I would certainly welcome if there is another explanation and appreciate the help! Ana

@samuel-marsh, if you would please able to re-open the issue until the question about the change in results may be addressed.

samuel-marsh commented 1 month ago

Hi Ana,

I'll happily reopen the issue for others to chime in but are you able tp provide reproducible example of results changing after running ScaleData. I'm unable to replicate with simple example below:

Best, Sam

library(tidyverse)
library(Seurat)
#> Loading required package: SeuratObject
#> Loading required package: sp
#> 'SeuratObject' was built under R 4.4.0 but the current version is
#> 4.4.1; it is recomended that you reinstall 'SeuratObject' as the ABI
#> for R may have changed
#> 
#> Attaching package: 'SeuratObject'
#> The following objects are masked from 'package:base':
#> 
#>     intersect, t
library(enrichR)
#> Welcome to enrichR
#> Checking connection ...
#> Enrichr ... Connection is Live!
#> FlyEnrichr ... Connection is Live!
#> WormEnrichr ... Connection is Live!
#> YeastEnrichr ... Connection is Live!
#> FishEnrichr ... Connection is Live!
#> OxEnrichr ... Connection is Live!
library(patchwork)

pbmc <- pbmc3k.SeuratData::pbmc3k.final
pbmc <- suppressMessages(UpdateSeuratObject(pbmc))
#> Warning: Assay RNA changing from Assay to Assay
#> Warning: Graph RNA_nn changing from Graph to Graph
#> Warning: Graph RNA_snn changing from Graph to Graph
#> Warning: DimReduc pca changing from DimReduc to DimReduc
#> Warning: DimReduc umap changing from DimReduc to DimReduc
#> Warning: Adding a command log without an assay associated with it
#> Adding a command log without an assay associated with it

pbmc[["idents"]] <- Idents(pbmc)

pbmc_labels <- pbmc@meta.data %>% 
  select(idents)

pbmc2 <- CreateSeuratObject(counts = pbmc@assays$RNA$counts, min.cells = 0, min.features = 0)

pbmc2 <- AddMetaData(pbmc2, metadata = pbmc_labels)

Idents(pbmc2) <- "idents"

pbmc2 <- NormalizeData(pbmc2)
#> Normalizing layer: counts

p1 <- DEenrichRPlot(object = pbmc2, ident.1 = "B", ident.2 = "NK", enrich.database = "GO_Molecular_Function_2023", max.genes = 100)
#> Uploading data to Enrichr... Done.
#>   Querying GO_Molecular_Function_2023... Done.
#> Parsing results... Done.
#> Uploading data to Enrichr... Done.
#>   Querying GO_Molecular_Function_2023... Done.
#> Parsing results... Done.

pbmc2 <- ScaleData(object = pbmc2, features = Features(pbmc2))
#> Centering and scaling data matrix

p2 <- DEenrichRPlot(object = pbmc2, ident.1 = "B", ident.2 = "NK", enrich.database = "GO_Molecular_Function_2023", max.genes = 100)
#> Uploading data to Enrichr... Done.
#>   Querying GO_Molecular_Function_2023... Done.
#> Parsing results... Done.
#> Uploading data to Enrichr... Done.
#>   Querying GO_Molecular_Function_2023... Done.
#> Parsing results... Done.

p1 / p2

Created on 2024-07-25 with reprex v2.1.1

Session info ``` r sessioninfo::session_info() #> ─ Session info ─────────────────────────────────────────────────────────────── #> setting value #> version R version 4.4.1 (2024-06-14) #> os macOS Big Sur 11.7.10 #> system x86_64, darwin20 #> ui X11 #> language (EN) #> collate en_US.UTF-8 #> ctype en_US.UTF-8 #> tz America/New_York #> date 2024-07-25 #> pandoc 3.1.11 @ /Applications/RStudio.app/Contents/Resources/app/quarto/bin/tools/x86_64/ (via rmarkdown) #> #> ─ Packages ─────────────────────────────────────────────────────────────────── #> package * version date (UTC) lib source #> abind 1.4-5 2016-07-21 [1] CRAN (R 4.4.0) #> cli 3.6.3 2024-06-21 [1] CRAN (R 4.4.0) #> cluster 2.1.6 2023-12-01 [1] CRAN (R 4.4.0) #> codetools 0.2-20 2024-03-31 [1] CRAN (R 4.4.1) #> colorspace 2.1-0 2023-01-23 [1] CRAN (R 4.4.0) #> cowplot 1.1.3 2024-01-22 [1] CRAN (R 4.4.0) #> curl 5.2.1 2024-03-01 [1] CRAN (R 4.4.0) #> data.table 1.15.4 2024-03-30 [1] CRAN (R 4.4.0) #> deldir 2.0-4 2024-02-28 [1] CRAN (R 4.4.0) #> digest 0.6.36 2024-06-23 [1] CRAN (R 4.4.0) #> dotCall64 1.1-1 2023-11-28 [1] CRAN (R 4.4.0) #> dplyr * 1.1.4 2023-11-17 [1] CRAN (R 4.4.0) #> enrichR * 3.2 2023-04-14 [1] CRAN (R 4.4.0) #> evaluate 0.24.0 2024-06-10 [1] CRAN (R 4.4.0) #> fansi 1.0.6 2023-12-08 [1] CRAN (R 4.4.0) #> farver 2.1.2 2024-05-13 [1] CRAN (R 4.4.0) #> fastDummies 1.7.3 2023-07-06 [1] CRAN (R 4.4.0) #> fastmap 1.2.0 2024-05-15 [1] CRAN (R 4.4.0) #> fitdistrplus 1.2-1 2024-07-12 [1] CRAN (R 4.4.0) #> forcats * 1.0.0 2023-01-29 [1] CRAN (R 4.4.0) #> fs 1.6.4 2024-04-25 [1] CRAN (R 4.4.0) #> future 1.33.2 2024-03-26 [1] CRAN (R 4.4.0) #> future.apply 1.11.2 2024-03-28 [1] CRAN (R 4.4.0) #> generics 0.1.3 2022-07-05 [1] CRAN (R 4.4.0) #> ggplot2 * 3.5.1 2024-04-23 [1] CRAN (R 4.4.0) #> ggrepel 0.9.5 2024-01-10 [1] CRAN (R 4.4.0) #> ggridges 0.5.6 2024-01-23 [1] CRAN (R 4.4.0) #> globals 0.16.3 2024-03-08 [1] CRAN (R 4.4.0) #> glue 1.7.0 2024-01-09 [1] CRAN (R 4.4.0) #> goftest 1.2-3 2021-10-07 [1] CRAN (R 4.4.0) #> gridExtra 2.3 2017-09-09 [1] CRAN (R 4.4.0) #> gtable 0.3.5 2024-04-22 [1] CRAN (R 4.4.0) #> highr 0.11 2024-05-26 [1] CRAN (R 4.4.0) #> hms 1.1.3 2023-03-21 [1] CRAN (R 4.4.0) #> htmltools 0.5.8.1 2024-04-04 [1] CRAN (R 4.4.0) #> htmlwidgets 1.6.4 2023-12-06 [1] CRAN (R 4.4.0) #> httpuv 1.6.15 2024-03-26 [1] CRAN (R 4.4.0) #> httr 1.4.7 2023-08-15 [1] CRAN (R 4.4.0) #> ica 1.0-3 2022-07-08 [1] CRAN (R 4.4.0) #> igraph 2.0.3 2024-03-13 [1] CRAN (R 4.4.0) #> irlba 2.3.5.1 2022-10-03 [1] CRAN (R 4.4.0) #> jsonlite 1.8.8 2023-12-04 [1] CRAN (R 4.4.0) #> KernSmooth 2.23-24 2024-05-17 [1] CRAN (R 4.4.1) #> knitr 1.48 2024-07-07 [1] CRAN (R 4.4.0) #> labeling 0.4.3 2023-08-29 [1] CRAN (R 4.4.0) #> later 1.3.2 2023-12-06 [1] CRAN (R 4.4.0) #> lattice 0.22-6 2024-03-20 [1] CRAN (R 4.4.1) #> lazyeval 0.2.2 2019-03-15 [1] CRAN (R 4.4.0) #> leiden 0.4.3.1 2023-11-17 [1] CRAN (R 4.4.0) #> lifecycle 1.0.4 2023-11-07 [1] CRAN (R 4.4.0) #> limma 3.60.4 2024-07-17 [1] Bioconductor 3.19 (R 4.4.1) #> listenv 0.9.1 2024-01-29 [1] CRAN (R 4.4.0) #> lmtest 0.9-40 2022-03-21 [1] CRAN (R 4.4.0) #> lubridate * 1.9.3 2023-09-27 [1] CRAN (R 4.4.0) #> magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.4.0) #> MASS 7.3-61 2024-06-13 [1] CRAN (R 4.4.0) #> Matrix 1.7-0 2024-04-26 [1] CRAN (R 4.4.1) #> matrixStats 1.3.0 2024-04-11 [1] CRAN (R 4.4.0) #> mime 0.12 2021-09-28 [1] CRAN (R 4.4.0) #> miniUI 0.1.1.1 2018-05-18 [1] CRAN (R 4.4.0) #> munsell 0.5.1 2024-04-01 [1] CRAN (R 4.4.0) #> nlme 3.1-165 2024-06-06 [1] CRAN (R 4.4.0) #> parallelly 1.37.1 2024-02-29 [1] CRAN (R 4.4.0) #> patchwork * 1.2.0 2024-01-08 [1] CRAN (R 4.4.0) #> pbapply 1.7-2 2023-06-27 [1] CRAN (R 4.4.0) #> pbmc3k.SeuratData 3.1.4 2024-07-24 [1] local #> pillar 1.9.0 2023-03-22 [1] CRAN (R 4.4.0) #> pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.4.0) #> plotly 4.10.4 2024-01-13 [1] CRAN (R 4.4.0) #> plyr 1.8.9 2023-10-02 [1] CRAN (R 4.4.0) #> png 0.1-8 2022-11-29 [1] CRAN (R 4.4.0) #> polyclip 1.10-7 2024-07-23 [1] CRAN (R 4.4.0) #> presto 1.0.0 2024-07-24 [1] Github (immunogenomics/presto@7636b3d) #> progressr 0.14.0 2023-08-10 [1] CRAN (R 4.4.0) #> promises 1.3.0 2024-04-05 [1] CRAN (R 4.4.0) #> purrr * 1.0.2 2023-08-10 [1] CRAN (R 4.4.0) #> R6 2.5.1 2021-08-19 [1] CRAN (R 4.4.0) #> RANN 2.6.1 2019-01-08 [1] CRAN (R 4.4.0) #> RColorBrewer 1.1-3 2022-04-03 [1] CRAN (R 4.4.0) #> Rcpp 1.0.13 2024-07-17 [1] CRAN (R 4.4.0) #> RcppAnnoy 0.0.22 2024-01-23 [1] CRAN (R 4.4.0) #> RcppHNSW 0.6.0 2024-02-04 [1] CRAN (R 4.4.0) #> readr * 2.1.5 2024-01-10 [1] CRAN (R 4.4.0) #> reprex 2.1.1 2024-07-06 [1] CRAN (R 4.4.0) #> reshape2 1.4.4 2020-04-09 [1] CRAN (R 4.4.0) #> reticulate 1.38.0 2024-06-19 [1] CRAN (R 4.4.0) #> rjson 0.2.21 2022-01-09 [1] CRAN (R 4.4.0) #> rlang 1.1.4 2024-06-04 [1] CRAN (R 4.4.0) #> rmarkdown 2.27 2024-05-17 [1] CRAN (R 4.4.0) #> ROCR 1.0-11 2020-05-02 [1] CRAN (R 4.4.0) #> RSpectra 0.16-2 2024-07-18 [1] CRAN (R 4.4.0) #> rstudioapi 0.16.0 2024-03-24 [1] CRAN (R 4.4.0) #> Rtsne 0.17 2023-12-07 [1] CRAN (R 4.4.0) #> scales 1.3.0 2023-11-28 [1] CRAN (R 4.4.0) #> scattermore 1.2 2023-06-12 [1] CRAN (R 4.4.0) #> sctransform 0.4.1 2023-10-19 [1] CRAN (R 4.4.0) #> sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.4.0) #> Seurat * 5.1.0 2024-05-10 [1] CRAN (R 4.4.0) #> SeuratObject * 5.0.2 2024-05-08 [1] CRAN (R 4.4.0) #> shiny 1.8.1.1 2024-04-02 [1] CRAN (R 4.4.0) #> sp * 2.1-4 2024-04-30 [1] CRAN (R 4.4.0) #> spam 2.10-0 2023-10-23 [1] CRAN (R 4.4.0) #> spatstat.data 3.1-2 2024-06-21 [1] CRAN (R 4.4.0) #> spatstat.explore 3.3-1 2024-07-15 [1] CRAN (R 4.4.0) #> spatstat.geom 3.3-2 2024-07-15 [1] CRAN (R 4.4.0) #> spatstat.random 3.3-1 2024-07-15 [1] CRAN (R 4.4.0) #> spatstat.sparse 3.1-0 2024-06-21 [1] CRAN (R 4.4.0) #> spatstat.univar 3.0-0 2024-06-28 [1] CRAN (R 4.4.0) #> spatstat.utils 3.0-5 2024-06-17 [1] CRAN (R 4.4.0) #> statmod 1.5.0 2023-01-06 [1] CRAN (R 4.4.0) #> stringi 1.8.4 2024-05-06 [1] CRAN (R 4.4.0) #> stringr * 1.5.1 2023-11-14 [1] CRAN (R 4.4.0) #> survival 3.7-0 2024-06-05 [1] CRAN (R 4.4.0) #> tensor 1.5 2012-05-05 [1] CRAN (R 4.4.0) #> tibble * 3.2.1 2023-03-20 [1] CRAN (R 4.4.0) #> tidyr * 1.3.1 2024-01-24 [1] CRAN (R 4.4.0) #> tidyselect 1.2.1 2024-03-11 [1] CRAN (R 4.4.0) #> tidyverse * 2.0.0 2023-02-22 [1] CRAN (R 4.4.0) #> timechange 0.3.0 2024-01-18 [1] CRAN (R 4.4.0) #> tzdb 0.4.0 2023-05-12 [1] CRAN (R 4.4.0) #> utf8 1.2.4 2023-10-22 [1] CRAN (R 4.4.0) #> uwot 0.2.2 2024-04-21 [1] CRAN (R 4.4.0) #> vctrs 0.6.5 2023-12-01 [1] CRAN (R 4.4.0) #> viridisLite 0.4.2 2023-05-02 [1] CRAN (R 4.4.0) #> withr 3.0.0 2024-01-16 [1] CRAN (R 4.4.0) #> WriteXLS 6.7.0 2024-07-20 [1] CRAN (R 4.4.0) #> xfun 0.46 2024-07-18 [1] CRAN (R 4.4.0) #> xml2 1.3.6 2023-12-04 [1] CRAN (R 4.4.0) #> xtable 1.8-4 2019-04-21 [1] CRAN (R 4.4.0) #> yaml 2.3.9 2024-07-05 [1] CRAN (R 4.4.0) #> zoo 1.8-12 2023-04-13 [1] CRAN (R 4.4.0) #> #> [1] /Library/Frameworks/R.framework/Versions/4.4-x86_64/Resources/library #> #> ────────────────────────────────────────────────────────────────────────────── ```
mhkowalski commented 1 month ago

I agree with Sam that DEenrichRPlot will have the same behavior as FindMarkers, so the results should not depend on whether scaled data is present or not. If you could provide additional detail or a reproducible example of what is going on, we'd be happy to try to help further.

cristanchoa commented 1 month ago

Hi all, my apologies for the delay in response. I think it might be something on our end because the servers have been hanging in weird places suddenly. Once I get it up and running again, if I can reproduce it, I'll post again. If I can't reproduce it but I figure out why, I'll post again for others to learn from. Thanks for your time!