Merge pull request satijalab#665 from satijalab/fix/fix_sct_findmarkers

Fix for incorrect logFC calculation on SCTassay

DESCRIPTION

@@ -1,5 +1,5 @@
 Package: Seurat
-Version: 4.2.0
+Version: 4.2.0.9001
 Date: 2022-09-21
 Title: Tools for Single Cell Genomics
 Description: A toolkit for quality control, analysis, and exploration of single cell RNA sequencing data. 'Seurat' aims to enable users to identify and interpret sources of heterogeneity from single cell transcriptomic measurements, and to integrate diverse types of single cell data. See Satija R, Farrell J, Gennert D, et al (2015) <doi:10.1038/nbt.3192>, Macosko E, Basu A, Satija R, et al (2015) <doi:10.1016/j.cell.2015.05.002>, Stuart T, Butler A, et al (2019) <doi:10.1016/j.cell.2019.05.031>, and Hao, Hao, et al (2020) <doi:10.1101/2020.10.12.335331> for more details.

NEWS.md

@@ -1,3 +1,8 @@
+# Unreleased
+
+## Changes
+- Fix bug in `FindMarkers()` when using `SCT` assay ([#6856](https://github.com/satijalab/seurat/issues/6586))
+
 # Seurat 4.2.0 (2022-09-21)
 
 ## Changes

R/differential_expression.R

@@ -752,6 +752,12 @@ FindMarkers.SCTAssay <- function(
     'scale.data' = GetAssayData(object = object, slot = "counts"),
     numeric()
   )
+  if (is.null(x = mean.fxn)){
+    mean.fxn <- function(x) {
+      return(log(x = rowMeans(x = expm1(x = x)) + pseudocount.use, base = base))
+
+    }
+  }
   fc.results <- FoldChange(
     object = object,
     slot = data.slot,
@@ -762,7 +768,7 @@ FindMarkers.SCTAssay <- function(
     mean.fxn = mean.fxn,
     fc.name = fc.name,
     base = base
-  )
+    )
   de.results <- FindMarkers(
     object = data.use,
     slot = data.slot,

tests/testthat/test_differential_expression.R

@@ -46,11 +46,11 @@ test_that("Default settings work as expected", {
 
   # SCT normalization
   expect_equal(results.sct[1, "p_val"], 6.225491e-11)
-  expect_equal(results.sct[1, "avg_logFC"], -0.6768721, tolerance = 1e-6)
+  expect_equal(results.sct[1, "avg_logFC"], -1.081321, tolerance = 1e-6)
   expect_equal(results.sct[1, "pct.1"], 0.111)
   expect_equal(results.sct[1, "pct.2"], 0.96)
   expect_equal(results.sct[1, "p_val_adj"], 1.369608e-08)
-  expect_equal(nrow(x = results.sct), 92)
+  expect_equal(nrow(x = results.sct), 158)
   expect_equal(rownames(x = results.sct)[1], "TYMP")
 })
 
@@ -96,8 +96,8 @@ test_that("setting pseudocount.use works", {
   expect_equal(results[1, "avg_logFC"], -2.630395, tolerance = 1e-6)
   expect_equal(nrow(x = results.clr), 182)
   expect_equal(results.clr[1, "avg_logFC"], -2.317338, tolerance = 1e-6)
-  expect_equal(nrow(results.sct), 185)
-  expect_equal(results.sct[1, "avg_logFC"], -1.845681, tolerance = 1e-6)
+  expect_equal(nrow(results.sct), 194)
+  expect_equal(results.sct[1, "avg_logFC"], -2.253920, tolerance = 1e-6)
 })
 
 results <- suppressWarnings(FindMarkers(object = pbmc_small, ident.1 = 0, ident.2 = 1, verbose = FALSE, base = exp(1), mean.fxn = rowMeans))
@@ -308,12 +308,12 @@ test_that("FindAllMarkers works as expected", {
   expect_equal(rownames(x = results.clr)[1], "HLA-DPB1")
 
   # SCT normalization
-  expect_equal(results.sct[1, "p_val"], 6.225491e-11)
-  expect_equal(results.sct[1, "avg_log2FC"], -1.265307, tolerance = 1e-6)
+  expect_equal(results.sct[1, "p_val"], 4.25861e-12)
+  expect_equal(results.sct[1, "avg_log2FC"], -2.70188, tolerance = 1e-6)
   expect_equal(results.sct[1, "pct.1"], 0.167)
   expect_equal(results.sct[1, "pct.2"], 0.909)
-  expect_equal(results.sct[1, "p_val_adj"], 1.369608e-08)
-  expect_equal(nrow(x = results.sct), 201)
+  expect_equal(results.sct[1, "p_val_adj"], 9.368941e-10)
+  expect_equal(nrow(x = results.sct), 210)
   expect_equal(rownames(x = results.sct)[1], "HLA-DPB1")
 
   # pseudocount.use = 0.1