I have added a new assay to my seurat object and I want to test for differential expression of all these features by a condition in a metadata column (which for this example can take the values of "condition1" or "condition2"). My command looks like:
dea <- FindMarkers(obj, ident.1 = "condition1", ident.2 = "condition2",
group.by = "condition", assay = "another_assay", features = rownames(obj@assays$another_assay),
test.use = "wilcox", min.pct = 0, logfc.threshold = 0, base = 2)
This does not return the test results for all the features in the assay :-/ I know the default for features = is all genes in the assay but i was not getting all the features, so i specifically put the rownames of the assay to specify it, but i still get less than I would expect.
I have also tried logfc.threshold = -Inf, I get the same output.
I have tried setting min.diff.pct = 0, I get even fewer features.
In addition to this, I have a question regarding the output I do get. The documentation says:
avg_logFC: log fold-chage of the average expression between the two groups.
For a given feature, I'm getting an avg_log2FC of 8.5 . Just to verify, I did:
feature1 <- data.frame(V1=GetAssayData(obj, assay = "another_assay")["feature1",])
feature1 <- merge(feature1, data.frame(condition=obj$condition), by='row.names')
avg_cond1 <- mean(feature1[which(feature1$condition == "condition1"), "V1"])
avg_cond2 <- mean(feature1[which(feature1$condition == "condition2"), "V1"])
log2(avg_cond1 / avg_cond2)
This does not equal 8.5
Any idea on what I'm doing wrong?
I'm running R 4.0.2 with Seurat 4.0.0
UPDATE: This gitHub issue solves my second question. However, I have a follow up question:
In their avg_logFC calculation they calculate the mean of e^x - 1 (where x is the expression values for a gene) and then log it (line 176/7).
My question is: Why do they calculate the mean of e^x - 1 instead of just the mean of x?
