What is pseudoreplication in single-cell RNA-seq?

It's treating individual cells as independent statistical replicates when the true unit of biological replication is the sample they came from, such as a donor, mouse, or patient. Cells from the same sample are correlated, so scoring thousands of cells as thousands of independent observations produces p-values far smaller than the data actually support.

Why does Seurat's FindMarkers give p-values like 1e-200?

FindMarkers runs its default tests, Wilcoxon or t-test, cell by cell, so the effective sample size is the cell count rather than the number of donors. With thousands of cells per group, standard errors shrink enormously and even small, biologically uninteresting differences reach extreme significance.

How do I make pseudobulk data from a Seurat object?

Use `AggregateExpression()` or `PseudobulkExpression()` with `group.by` set to the combination of condition, sample_id, and cell type annotation, summing raw counts within each group. Feed the resulting matrix into DESeq2's `DESeqDataSetFromMatrix()` or edgeR's `DGEList` and `glmQLFit` workflow with sample, not cell, as the unit of observation.

Is the Wilcoxon test in Seurat always wrong to use?

No. It's appropriate for within-dataset marker discovery, like finding genes that distinguish one cluster from another in the samples you already profiled. It becomes pseudoreplication when you use it to claim a condition effect, treatment versus control or disease versus healthy, that you intend to generalize beyond the specific samples in hand.

How many cells do I need per sample for pseudobulk differential expression?

There's no universal number, but the muscat vignette flags sample-cluster combinations under about 50 cells as candidates for exclusion because the aggregated pseudobulk profile gets noisy below that. What matters more than the exact cutoff is having enough distinct biological samples per condition, ideally three or more, not more cells per sample.

Sanity check · Single-Cell RNA-seq

How to Avoid Pseudoreplication in Single-Cell RNA-seq

A p-value of 1e-200 from FindMarkers usually means you tested cells instead of donors, and pseudobulk is the fix a reviewer already knows to ask for.

By Ming "Tommy" Tang, Director of Bioinformatics in Big Pharma · Reviewed September 2026 · 4 min read

You clustered your 10x dataset, annotated the cell types, then ran FindMarkers() between stimulated and control cells inside your favorite cluster. The p-value came back 1e-180. You've never seen a p-value that small in a bulk RNA-seq experiment, and that should make you suspicious rather than pleased.

Here's what's at stake: that gene list isn't a validated set of hits, it's a list of genes correlated with which donor happened to contribute which cells. Send it to the wet lab for qPCR and half of it won't replicate, because the test never asked whether the effect holds across donors, only whether individual cells differ from each other.

In the next hour you can settle this: count how many biological samples actually sit behind each condition, cross-tab cluster membership against sample ID, and rebuild the comparison as pseudobulk per sample so DESeq2 or edgeR treats the donor, not the cell, as the unit of replication.

What it looks like when it's happening

FindMarkers between condition groups returns p-values in the 1e-50 to 1e-300 range for hundreds or thousands of genes
The number of genes passing padj < 0.05 is far larger than you'd ever see in a comparable bulk RNA-seq experiment with the same design
Dropping or subsetting a single donor changes the top DE gene list substantially
A volcano plot shows a wall of points hugging the smallest representable p-value with little gradient of significance
p-values are extreme while log2 fold changes are modest, so significance and effect size are decoupled
One cluster's condition effect is really just one donor: a cross-tab of cluster by sample shows that cluster is 80-90% cells from a single sample
The comparison quietly used cell count as the sample size (thousands) instead of the number of donors or biological samples (often 2-6)

Why it happens

FindMarkers by default runs a Wilcoxon rank-sum test (or a t-test) cell by cell. Every cell in your stimulated cluster counts as one independent observation, so a comparison built from 3 donors and 3,000 cells per condition gets scored as though n were 3,000, not 3. Standard errors shrink with the square root of that inflated n, and p-values collapse toward zero for even small, biologically uninteresting differences.

The biology behind this: cells from the same donor are not independent draws from the same distribution. They share genotype, the physiological state of the donor at the time of sampling, freeze-thaw history, and an ambient RNA background that alone can make up 3-35% of counts per cell and varies by sample. That shared structure means the variance you should be testing against is the variance between donors, not the variance between cells within a donor. A cell-level test silently substitutes the wrong denominator and reports a confidence you don't have.

This compounds with double dipping: you use the same expression matrix to define clusters and then to test for condition differences inside those clusters, which biases p-values downward before pseudoreplication even enters the picture. If you clustered on integrated coordinates (Harmony, scVI), that step further blends donor-specific structure into a shared graph, so clusters look consistent across donors even when donor identity, not treatment, was driving part of the split.

None of this is rescued by sequencing more cells. One single-cell sample profiled across four conditions gives you weaker evidence than three true biological replicates per condition run as bulk RNA-seq, because no amount of cells substitutes for replication at the level that actually varies from experiment to experiment: the donor.

The checks

Run them in order. Each one tells you what healthy looks like and what the problem looks like.

0/7 checked · saved in this browser

01Count your real n, not your cell count
Tabulate cells against sample ID and condition before you interpret any DE result. Run table(seu_obj$sample_id, seu_obj$condition) and count the number of unique sample_id values per condition, not the cell totals in that table.
r
```
table(seu_obj$sample_id, seu_obj$condition)
```
Healthy
Three or more unique biological samples (donors, mice, patients) represented in each condition group you're comparing.
Red flag
Only one or two unique samples per condition, meaning any apparent difference could just be donor identity rather than treatment effect.
02Plot the raw p-value distribution
Pull the raw p_val column from FindMarkers output (not padj) and plot a histogram across all tested genes.
r
```
hist(markers$p_val, breaks = 50)
```
Healthy
A modest bump of small p-values near zero against a broad background near 1, roughly consistent with the number of genes you'd expect to genuinely differ.
Red flag
A spike of hundreds or thousands of genes piled at the smallest representable p-value, often 1e-100 or smaller.
03Cross-tab cluster membership against sample ID
Build a cluster-by-sample contingency table to see how evenly cells are distributed across your biological samples within each cluster you plan to test.
r
```
table(seu_obj$seurat_clusters, seu_obj$sample_id)
```
Healthy
Every cluster has meaningful representation, ideally at least the ~50-cell muscat guideline, from every sample in every condition you're comparing.
Red flag
A cluster where the large majority of cells trace back to a single donor: any condition effect you find there is confounded with that donor.
04Leave-one-donor-out and rerun
Drop one biological sample entirely and rerun the same DE test on the remaining samples, then compare the top genes and their ranking.
r
```
seu_sub <- subset(seu_obj, subset = sample_id != "donor1")
markers_loo <- FindMarkers(seu_sub, ident.1 = "stim", ident.2 = "ctrl", group.by = "condition")
```
Healthy
The top genes and their approximate ranking stay largely stable when one sample is removed.
Red flag
The top-20 gene list reshuffles substantially, or genes that were highly significant disappear once one donor is dropped.
05Check cell counts per sample-cluster combination
Cross-tab clusters by samples and flag any combination below roughly 50 cells, the guideline muscat uses before excluding a combination from pseudobulk aggregation.
r
```
counts <- table(seu_obj$seurat_clusters, seu_obj$sample_id)
which(counts < 50, arr.ind = TRUE)
```
Healthy
The sample-cluster combinations you intend to compare all sit comfortably above the ~50-cell threshold.
Red flag
Several sample-cluster combinations under 50 cells, especially clustered in one condition, which will make that pseudobulk sample noisy or effectively missing.
06Rebuild as pseudobulk and compare results
Aggregate raw UMI counts per sample x cell type x condition, run the standard bulk DE pipeline, and compare the resulting gene list and p-value magnitudes to the original cell-level test.
r
```
pb <- AggregateExpression(
  seu_obj, assays = "RNA", return.seurat = TRUE,
  group.by = c("condition", "sample_id", "seurat_annotations")
)
```
Healthy
A shorter, more modest gene list with p-values in a plausible bulk RNA-seq range, and log2 fold changes that track the direction of the original cell-level result.
Red flag
The pseudobulk analysis returns zero or near-zero significant genes while the cell-level test returned thousands, confirming the original signal was mostly inflated by pseudoreplication.
07Check whether condition is confounded with batch or donor
Cross-tab condition against any known technical batch variable (10x run, sequencing lane, processing day) to see whether the two are entangled.
r
```
table(seu_obj$condition, seu_obj$batch)
```
Healthy
Multiple samples and multiple batches represented within each condition, so condition and batch are crossed, not perfectly aligned.
Red flag
Each condition maps to an entirely separate set of donors or batches, a 1:1 confound that no statistical method can untangle after the fact.

What to do about it

Aggregate to pseudobulk, test with DESeq2 or edgeR

When: You have three or more biological samples per condition and want a defensible cross-condition DE list per cell type.

Sum raw UMI counts per sample x cell type with AggregateExpression(seu_obj, assays = "RNA", return.seurat = TRUE, group.by = c("condition","sample_id","seurat_annotations")), extract the counts matrix, build sample-level metadata, then run DESeqDataSetFromMatrix(countData = pb_counts, colData = sample_metadata, design = ~ condition) followed by DESeq(), or the edgeR DGEList -> calcNormFactors -> estimateDisp -> glmQLFit -> glmQLFTest pipeline. Sum counts, don't average; sum-based aggregation consistently controls false positives better than mean-based aggregation.

Caveat: Your usable sample size drops to the number of biological samples, not cells, so with only 1-2 samples per condition you don't have the degrees of freedom for a real test, and pseudobulk averages away within-sample heterogeneity like rare subclusters.

Use mixed models to keep cell-level resolution

When: You need to preserve within-sample state transitions or subcluster resolution that pseudobulk would average away, and you have enough samples to fit a random effect.

Use muscat's mixed-model differential state methods (MAST_RE, NEBULA-LN) or the mixed-model modes of pbDS(), which model sample as a random effect instead of collapsing counts.

Caveat: Slower to fit and more sensitive to model misspecification than pseudobulk; benchmarking across 18 methods still ranks pseudobulk ahead of mixed models on precision and specificity, so reach for this only when aggregation genuinely loses signal you need.

Fix the study design before you generate data

When: You're still planning the experiment, not analyzing existing data.

Budget for three or more true biological replicates per condition rather than spending the budget on many cells from one or two samples across several conditions. One scRNA-seq sample profiled across four conditions gives less reliable insight than three biological replicates per condition analyzed as bulk or pseudobulk.

Caveat: More samples costs more money and reagent, and if the biology you care about is a rare subtype only visible at single-cell resolution, you may still need scRNA-seq depth even at a lower replicate count.

Drop underpowered sample-cluster combinations

When: Some clusters have too few cells from one or more samples to aggregate into a reliable pseudobulk profile.

Apply the muscat guideline of roughly 50 cells per sample-cluster combination: filter out or explicitly flag any combination below that threshold before aggregating, rather than letting a near-empty pseudobulk sample dilute the test.

Caveat: Aggressive filtering shrinks your effective sample size further and biases the analysis toward abundant cell types, so report which clusters or samples you excluded and why.

Remove donor or batch confounding per cell type with RUV

When: PCA of your pseudobulk samples shows batch- or donor-driven variation stacked on top of the condition effect.

Apply RUV2, RUVIII, or RUV4 separately within each cell type's pseudobulk matrix, using known negative-control genes, or RUVIII_PBPS synthetic pseudoreplicates when you have no true negative controls, then feed the estimated unwanted-variation factors into the DESeq2 or edgeR design.

Caveat: Running RUV across all cell types pooled together instead of per cell type misestimates the unwanted variation and can strip out real condition signal along with the batch effect.

When not to "fix" it

If your comparison is within-dataset marker discovery for cell type annotation, asking which genes distinguish this cluster from that cluster in the samples you actually profiled, cell-level FindMarkers is the right tool as is. You aren't generalizing to a third donor, you're describing the dataset in front of you. Likewise, if you have genuine technical replicates (the same donor's library split across two lanes) and your question is about technical reproducibility of capture or library prep rather than a biological treatment effect, treating those cells as repeated measures is exactly the intended use, not pseudoreplication.

Five things experienced analysts do here

Before opening Seurat, write down n as the number of independent biological units for your comparison; if that number is 1 or 2 per condition, no downstream test recovers it.
Always pseudobulk from raw or RNA-assay counts, never from integrated (Harmony, scVI) or log-normalized coordinates, since those steps were built to blend the variation you need to keep for testing.
Run the cluster-by-sample cross-tab as a standing QC step for every cross-condition DE comparison, the same way you'd check a rarefaction curve, before you look at a single p-value.
Treat any FindMarkers p-value below roughly 1e-50 in a cross-condition test as a prompt to check for pseudoreplication first, then rank real candidates by log2FC and cross-donor consistency instead of p-value alone.
When a paper or vignette shows a cell-level DE workflow, check whether the comparison is within one sample (fine) or across biological conditions (needs pseudobulk or mixed models) before copying it into your own analysis.

Questions people ask

What is pseudoreplication in single-cell RNA-seq?: It's treating individual cells as independent statistical replicates when the true unit of biological replication is the sample they came from, such as a donor, mouse, or patient. Cells from the same sample are correlated, so scoring thousands of cells as thousands of independent observations produces p-values far smaller than the data actually support.
Why does Seurat's FindMarkers give p-values like 1e-200?: FindMarkers runs its default tests, Wilcoxon or t-test, cell by cell, so the effective sample size is the cell count rather than the number of donors. With thousands of cells per group, standard errors shrink enormously and even small, biologically uninteresting differences reach extreme significance.
How do I make pseudobulk data from a Seurat object?: Use AggregateExpression() or PseudobulkExpression() with group.by set to the combination of condition, sample_id, and cell type annotation, summing raw counts within each group. Feed the resulting matrix into DESeq2's DESeqDataSetFromMatrix() or edgeR's DGEList and glmQLFit workflow with sample, not cell, as the unit of observation.
Is the Wilcoxon test in Seurat always wrong to use?: No. It's appropriate for within-dataset marker discovery, like finding genes that distinguish one cluster from another in the samples you already profiled. It becomes pseudoreplication when you use it to claim a condition effect, treatment versus control or disease versus healthy, that you intend to generalize beyond the specific samples in hand.
How many cells do I need per sample for pseudobulk differential expression?: There's no universal number, but the muscat vignette flags sample-cluster combinations under about 50 cells as candidates for exclusion because the aggregated pseudobulk profile gets noisy below that. What matters more than the exact cutoff is having enough distinct biological samples per condition, ideally three or more, not more cells per sample.

Guide · How to Avoid Pseudoreplication in Spatial Transcriptomics
Compare · DESeq2 vs edgeR: Which One Should You Use?
Guide · How to Choose Cell QC Thresholds in Single-Cell RNA-seq
Guide · How to Detect Batch Effects in Bulk RNA-seq
Guide · How to Find and Remove Doublets in Single-Cell RNA-seq
Glossary · False discovery rate (FDR)
Glossary · Log fold change (log2FC)
Glossary · Pseudobulk

Sources

Differential expression testing • Seurat — FindMarkers default cell-level tests, AggregateExpression pseudobulk workflow, and DESeq2 integration
How to create pseudobulk from single-cell RNAseq data — Practical pseudobulk aggregation workflow using Seurat and presto
Benchmarking methods for detecting differential states between conditions from multi-subject single-cell RNA-seq data — Pseudobulk vs mixed model vs naive cell-level method false positive rate comparison
Differential state analysis with muscat — aggregateData() workflow and the ~50 cell per sample-cluster combination guideline
Removal of unwanted variation in pseudobulk analysis of single-cell RNA sequencing data — RUV2/RUVIII applied per cell type to control batch effects confounded with condition
Pseudobulk Expression, PseudobulkExpression • Seurat — Aggregate vs average method options and scale.factor parameter for pseudobulk construction

Part of the Pseudoreplication series.