FindMarkers with MAST using SCT assay, erroneous result with slots = data ?

[maxc2020]

I'm getting variable results with FindMarkers and looking for input on how to approach it, or if it is a bug to share it.

A SeuratObject was constructed with a variety of samples normalized with SCT v2. Now I am looking at Differential Gene expression across celltypes or conditions, with multiple samples in the condition. When I run comparison with FindMarkers and MAST using RNA assay (slots as counts or data), MAST using SCT assay (slots as data), or Wilcoxan test with RNA or SCT (slots as counts or data), a cell-type specific gene for Celltype1 vs Celltype2 has a log2_FC of ~5. When I run FindMarkers and MAST using SCT assay (slots as counts), I get a log2_FC of ~300. When I run my data in MAST package (e.g., convert to SCE format), the log2_FC is ~5. Using AverageExpression function, the expression of this gene in SCT assay is 60.519564 for Celltype1 and 1.126311 for Celltype2 in counts slot and 60.519564 and 1.126311 in data slots (and grossly similar in RNA assay). So it would appear that FindMarker with test.use = MAST has an error when assay = SCT and slot = Counts. Or is there another interpretation here?

Examples of the code lines to get DE:

FindMarkers(Obj, ident.1 = c("Celltype1"), ident.2 = c("Celltype2"), min.pct = 0.25, min.diff.pct = 0.25, only.pos = F, assay = "SCT", slot = "counts", recorrect_umi = F, test.use = "MAST", verbose = T)
FindMarkers(Obj, ident.1 = c("Celltype1"), ident.2 = c("Celltype2"), min.pct = 0.25, min.diff.pct = 0.25, only.pos = F, assay = "SCT", slot = "data", recorrect_umi = F, test.use = "MAST", verbose = T)
FindMarkers(Obj, ident.1 = c("Celltype1"), ident.2 = c("Celltype2"), min.pct = 0.25, min.diff.pct = 0.25, only.pos = F, assay = "RNA", slot = "counts", test.use = "MAST", verbose = T)

SessionInfo

R version 4.2.2 (2022-10-31)
Platform: aarch64-apple-darwin20 (64-bit)
Running under: macOS Ventura 13.2.1

Matrix products: default
LAPACK: /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] dittoSeq_1.10.0             parallelly_1.35.0           future_1.32.0               openxlsx_4.2.5.2           
 [5] Matrix_1.5-3                ggplot2_3.4.1               dplyr_1.1.1                 patchwork_1.1.2            
 [9] cowplot_1.1.1               MAST_1.25.1                 harmony_0.1.1               Rcpp_1.0.10                
[13] batchelor_1.14.1            SingleCellExperiment_1.20.1 SummarizedExperiment_1.28.0 Biobase_2.58.0             
[17] GenomicRanges_1.50.2        GenomeInfoDb_1.34.9         IRanges_2.32.0              S4Vectors_0.36.2           
[21] BiocGenerics_0.44.0         MatrixGenerics_1.10.0       matrixStats_0.63.0          SeuratWrappers_0.3.1       
[25] sctransform_0.3.5.9002      Seurat_4.9.9.9041           SeuratObject_4.9.9.9081     sp_1.6-0                   

loaded via a namespace (and not attached):
  [1] utf8_1.2.3                spatstat.explore_3.1-0    reticulate_1.28           R.utils_2.12.2           
  [5] tidyselect_1.2.0          htmlwidgets_1.6.2         grid_4.2.2                BiocParallel_1.32.6      
  [9] Rtsne_0.16                devtools_2.4.5            munsell_0.5.0             ScaledMatrix_1.6.0       
 [13] codetools_0.2-19          ica_1.0-3                 miniUI_0.1.1.1            withr_2.5.0              
 [17] spatstat.random_3.1-4     colorspace_2.1-0          progressr_0.13.0          rstudioapi_0.14          
 [21] ROCR_1.0-11               tensor_1.5                listenv_0.9.0             GenomeInfoDbData_1.2.9   
 [25] polyclip_1.10-4           pheatmap_1.0.12           rprojroot_2.0.3           vctrs_0.6.1              
 [29] generics_0.1.3            R6_2.5.1                  rsvd_1.0.5                bitops_1.0-7             
 [33] spatstat.utils_3.0-2      cachem_1.0.7              DelayedArray_0.24.0       promises_1.2.0.1         
 [37] scales_1.2.1              gtable_0.3.3              beachmat_2.14.0           globals_0.16.2           
 [41] processx_3.8.0            goftest_1.2-3             spam_2.9-1                rlang_1.1.0              
 [45] splines_4.2.2             lazyeval_0.2.2            spatstat.geom_3.1-0       BiocManager_1.30.20      
 [49] reshape2_1.4.4            abind_1.4-5               httpuv_1.6.9              usethis_2.1.6            
 [53] tools_4.2.2               ellipsis_0.3.2            RColorBrewer_1.1-3        sessioninfo_1.2.2        
 [57] ggridges_0.5.4            plyr_1.8.8                progress_1.2.2            sparseMatrixStats_1.10.0 
 [61] zlibbioc_1.44.0           purrr_1.0.1               RCurl_1.98-1.12           ps_1.7.3                 
 [65] prettyunits_1.1.1         deldir_1.0-6              pbapply_1.7-0             urlchecker_1.0.1         
 [69] zoo_1.8-11                ggrepel_0.9.3             cluster_2.1.4             fs_1.6.1                 
 [73] magrittr_2.0.3            data.table_1.14.8         RSpectra_0.16-1           scattermore_0.8          
 [77] ResidualMatrix_1.8.0      lmtest_0.9-40             RANN_2.6.1                fitdistrplus_1.1-8       
 [81] pkgload_1.3.2             hms_1.1.3                 mime_0.12                 xtable_1.8-4             
 [85] fastDummies_1.6.3         gridExtra_2.3             compiler_4.2.2            tibble_3.2.1             
 [89] KernSmooth_2.23-20        crayon_1.5.2              R.oo_1.25.0               htmltools_0.5.5          
 [93] later_1.3.0               tidyr_1.3.0               MASS_7.3-58.3             cli_3.6.1                
 [97] R.methodsS3_1.8.2         parallel_4.2.2            dotCall64_1.0-2           igraph_1.4.1             
[101] pkgconfig_2.0.3           plotly_4.10.1             scuttle_1.8.4             spatstat.sparse_3.0-1    
[105] XVector_0.38.0            stringr_1.5.0             callr_3.7.3               digest_0.6.31            
[109] RcppAnnoy_0.0.20          spatstat.data_3.0-1       leiden_0.4.3              uwot_0.1.14              
[113] DelayedMatrixStats_1.20.0 curl_5.0.0                shiny_1.7.4               lifecycle_1.0.3          
[117] nlme_3.1-162              jsonlite_1.8.4            BiocNeighbors_1.16.0      limma_3.54.2             
[121] desc_1.4.2                viridisLite_0.4.1         fansi_1.0.4               pillar_1.9.0             
[125] lattice_0.20-45           fastmap_1.1.1             httr_1.4.5                pkgbuild_1.4.0           
[129] survival_3.5-5            glue_1.6.2                remotes_2.4.2             zip_2.2.2                
[133] png_0.1-8                 stringi_1.7.12            profvis_0.3.7             RcppHNSW_0.4.1           
[137] BiocSingular_1.14.0       memoise_2.0.1             irlba_2.3.5.1             future.apply_1.10.0      

[lepriolc]

Hi @maxc2020,

I am not part of Seutat team but I have also noticed this kind of error in average log2FC computation made by FindMarkers function. It appears that the function to calculate average expression is not correct when using an SCT assay and one the following tests : negbinom, poisson, MAST. I described the error in more details in the issue #7095 and pulled a request to fix this error (#7115).

Christophe

[maxc2020]

Thanks for the note Cristophe and sorry I didn't catch your issue submission.

[saketkc]

Hi @lepriolc, @maxc2020
I have just pushed 992c1a9 to develop branch which should fix this and related issues. Can you pull from the develop branch and re-open this issue if you are still facing issues?

devtools::install_github("satijalab/seurat", ref="develop")