What is Dr. Zubair Khalid's research focus?

Dr. Zubair Khalid specializes in molecular virology, mRNA vaccine development, and computational biology, with a focus on avian pathogens like IBDV and Avian Reovirus.

Where is Dr. Zubair Khalid currently working?

Dr. Zubair Khalid is a Postdoctoral Research Associate at the University of Maryland (UMD), specifically within the Department of Animal and Avian Sciences.

Deseq2 Single Cell Rna Seq

DESeq2 is a differential expression (DE) analysis tool originally designed for bulk RNA sequencing data. When applied to single cell RNA sequencing (scRNA seq), it is most effective after aggregating individual cell counts into pseudobulk samples. This approach allows researchers to compare expression between groups of samples while leveraging DESeq2’s robust normalization and statistical modeling. You should use this guide if you are a computational biologist or bioinformatician analyzing scRNA seq data with multiple biological replicates, and you need a reliable method for identifying genes that change between conditions. The Bioconductor package provides extensive documentation for DESeq2 Bioconductor.

Pseudobulk aggregation sums raw counts from all cells belonging to the same sample or experimental condition. This step converts the sparse single cell matrix into a bulk like count table suitable for DESeq2. The Galaxy Training Network offers a clear workflow for scRNA seq pseudobulk analysis Galaxy Training Network. By working with aggregated counts, you reduce dropout noise and focus on sample level variation, which is the scale at which biological replicates differ.

At a Glance

Aspect	Description
Tool	DESeq2 (R package from Bioconductor)
Input	Pseudobulk count matrix (sample x gene) and sample metadata
Output	Log2 fold changes, adjusted p values, normalized counts
Best Use Case	Comparing two or more conditions across biological replicates in scRNA seq experiments
Limitations	Loses per cell resolution, requires at least 2 3 samples per condition, assumes independent samples

Decision Criteria

Before using DESeq2 with single cell data, evaluate these criteria:

Biological replicates are essential. DESeq2 models dispersion across replicates. A single sample per condition is not acceptable unless you use a different dispersion estimate method which is still not recommended. Aim for at least three samples per condition.
Number of cells per sample matters. Pseudobulk samples with very few cells (e.g., under 50) may produce unstable estimates. A practical lower bound is 100 200 cells per sample, though guidelines vary by tissue and platform.
Dropout and sparsity are mitigated by aggregation. If your study focuses on rare subpopulations or requires cell level resolution, consider single cell specific DE methods such as MAST or SCTransform instead.
Batch effects must be addressed. If your samples were processed in different batches, include batch as a covariate in the DESeq2 formula. Batch effects at single cell level can be corrected with tools like Harmony before aggregation, but ensure the correction does not remove biological signal.

The NCBI Sequence Read Archive provides extensive metadata to help you design experiments with proper replication NCBI Sequence Read Archive.

Practical Workflow

Step 1: Obtain a Raw Count Matrix per Cell

Start with a standard scRNA seq pipeline (e.g., Cell Ranger from 10x Genomics or STARsolo). The output is a genes x cells count matrix. Keep only the raw unique molecular identifier (UMI) counts, do not apply normalization yet.

Step 2: Aggregate Counts into Pseudobulk Samples

For each biological sample (e.g., each individual or each well in a plate based experiment), sum the UMI counts across all cells belonging to that sample. The result is a gene x sample matrix. Use the R package Seurat, scran, or DropletUtils to perform this aggregation. Ensure your metadata maps each cell to its original sample identifier.

Step 3: Create a DESeqDataSet

Load the pseudobulk count matrix into R and build a DESeqDataSet object. The design formula should include the condition of interest and any known covariates (e.g., batch, sex). Example:

dds <- DESeqDataSetFromMatrix(countData = psbulk_counts,
                              colData = metadata,
                              design = ~ batch + condition)

Step 4: Run the DESeq2 Pipeline

Apply the DESeq() function. This performs normalization, dispersion estimation, and statistical testing. For scRNA seq pseudobulk, use the default settings. If you have many genes with very low counts, using minReplicatesForReplace = Inf can avoid replacement of outliers. The EMBL EBI Training materials explain the underlying model in detail EMBL EBI Training.

Step 5: Filter and Extract Results

After analysis, filter results using a reasonable adjusted p value threshold (e.g., 0.05) and a minimum log2 fold change (e.g., 1). Use results(dds) to extract the comparison of interest. The lfcShrink() function can provide more stable fold change estimates.

Step 6: Perform Quality Checks

Dispersion plot: Check that the dispersion estimates follow the fitted curve. Excessive scatter may indicate unmodeled batch effects or sample outliers.
P value histogram: A uniform distribution for high p values is expected. A spike near zero indicates true positives, a spike near one may suggest model misspecification.
MA plot: Visualize the relationship between mean expression and log2 fold change. Genes with low counts should be filtered out.

A study on vascular smooth muscle cells in thromboembolic hypertension used similar aggregation steps before differential expression analysis, demonstrating the applicability of this workflow PubMed article 10.

Common Mistakes

Treating individual cells as replicates. This inflates sample size and underestimates variance, leading to false positives. Always aggregate by biological sample.
Using normalized or imputed counts as input for DESeq2. DESeq2 expects raw integer counts. Normalized values break its internal modeling. Work from the original UMI matrix.
Insufficient filtering of lowly expressed genes. Remove genes with very low total counts across samples (e.g., fewer than 10 counts in total). Alternatively, let DESeq2’s independent filtering handle this, but pre filtering can stabilize dispersion estimates.
Ignoring batch effects in the design formula. If batches are confounded with conditions, you may not be able to separate biological and technical variation. Plan experiments to avoid confounding, or use a paired design.
Interpreting fold changes from pseudobulk as per cell expression changes. Pseudobulk reflects the total expression in a sample. A gene may appear upregulated due to an increase in cell proportion rather than per cell expression. Consider complementary analyses like gene set enrichment or cell type proportion testing.

One colorectal cancer study used DESeq2 to identify DEGs from bulk data, but the same principles of careful filtering apply to pseudobulk scRNA seq PubMed article 8.

Limits and Uncertainty

Pseudobulk analysis with DESeq2 is a powerful technique but has inherent limitations:

Loss of single cell resolution. You cannot study heterogeneity within a sample. If your question concerns rare subpopulations or cell state transitions, use a method that works on individual cells.
Assumption of independence between samples. DESeq2 assumes each sample is an independent observation. This is valid if samples come from distinct individuals. For paired designs (e.g., pre/post treatment from the same subject), you must account for pairing in the design formula.
Dropout events are not explicitly modeled. Aggregation reduces dropout impact, but very lowly expressed genes may still be missed. Interpreting a non significant gene does not mean it is not expressed in a small subset of cells.
Dispersion estimates can be unstable with few samples. The method works best with at least three samples per group. With two samples per condition, you may accept higher false discovery rates.

A preprint on vascular cognitive impairment used pseudobulk analysis to identify receptor dysregulation, highlighting the need for caution when drawing conclusions from modest sample sizes PubMed article 11.

Frequently Asked Questions

Q: Can I run DESeq2 directly on single cell raw counts without aggregating?
No. Doing so treats each cell as an independent replicate, which violates the assumption of independent samples and underestimates biological variance. The result is many false positive DEGs.

Q: How many cells per sample are needed for a reliable pseudobulk analysis?
There is no fixed number, but we recommend at least 100 cells per sample after quality filtering. Samples with fewer than 50 cells often produce noisy count profiles and unstable dispersion estimates.

Q: Should I use sizeFactors from DESeq2 or compute normalization separately?
Use DESeq2’s internal median of ratios normalization. It accounts for library size differences across pseudobulk samples. Avoid applying CPM or TPM normalization before DESeq2, as it interferes with the model.

Q: Can I use DESeq2 for differential expression between cell types within a single sample?
No. DESeq2 compares groups of biological replicates. For within sample cell type comparisons, use methods designed for single cell data such as MAST or Wilcoxon rank sum tests applied to the per cell expression matrix.

Zubair Khalid