This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Statistical test of difference between two conditions and account for the batch

Hi all,

I work on a scRNA-seq dataset. I computed a module score (AddModuleScore of Seurat) and I want to test if the difference in this score is statistically significant between two conditions.

I have three batches and each contain the two conditions (KO_batch1, WT_batch1, KO_batch2, WT_batch2, KO_batch3, WT_batch3).

I would normally use a Wilcoxon or t-test between condition 1 and condition 2, but doing so I wouldn't account for the batches. It would be testing all cells independently while they're not independent and therefore I would get really low p-values while I shouldn't. I could aggregate the scores, average them per batch and then test if the difference in batches is significant. But I would have a very low power as I would test only three differences somehow (3 batches: module score in condition2 vs. condition1).

So I wonder what is the best to do here? I've heard about linear mixed effect model but I am confused about what they are and how to use them and if they're suited here.

A big thank you for your help:)

significance-testing test seurat scrna-seq batch
• 10,385 views • link
updated by Kevin Blighe Kevin Blighe • written by npont • npont

Are the batches biological replicates?

• 0 views • link
written by ATpoint ATpoint

Great, then I would make a pseudobulk analysis. Test contrasts with a DE tool of choice, for example edgeR, and then use geneset enrichment analysis, for example camera from limma, to test whether there is difference for these genesets (aka modules) across conditions. This is a lot more statistically robust than going with these module scores directly.

• 0 views • link
written by ATpoint ATpoint

Thanks a lot very helpful :) Do you know anything about linear mixed effect models applied to pseudobulk DGE also?

• 0 views • link
written by npont • npont

limma and dream, read about these packages if you want to do that.

• 0 views • link
written by ATpoint ATpoint

1 answer

I agree with the pseudobulk approach as a robust way to handle this, avoiding pseudoreplication from treating cells as independent.

First, aggregate your counts to pseudobulk level per sample (i.e., per batch-condition pair: KO_batch1, WT_batch1, etc.). Use functions like AggregateExpression in Seurat or the muscat package for this.

Then, perform differential expression with a tool that accounts for the paired design (since batches are biological replicates with paired conditions). edgeR or DESeq2 work well; include batch as a pairing factor in the model.

For the module, instead of testing the score directly, run gene set enrichment on the DE results using camera from limma or fgsea. This tests if your gene set is enriched in the differential signal between KO and WT.

If you prefer linear mixed effects on pseudobulk, check out dream from the variancePartition package—it's designed for RNA-seq and handles random effects like batch elegantly. Read the vignettes for dream and limma on Bioconductor.

This keeps things statistically sound with your n=3 replicates. Avoid cell-level tests without correction.

Kevin

• 0 views • link
written by Kevin Blighe Kevin Blighe

Log in to answer this question.

More posts like this