uv run python with the skill's Python package installed; no API key or account required. Written for PyDESeq2 0.5.x (current 0.5.4 as of 2026-09), which requires Python >= 3.11.\"\nmetadata:\nskill-author: AlterLab\nversion: \"1.1.0\"\nlast_updated: \"2026-09-23\"PyDESeq2 is a Python implementation of DESeq2 for differential expression analysis with bulk RNA-seq data. It supports complete workflows from data loading through result interpretation, including single-factor and multi-factor designs, Wald tests with multiple-testing correction, optional apeGLM shrinkage, and integration with pandas and AnnData.
\nUse this skill when:
\n| Scenario | \nUse Instead | \n
|---|---|
| Single-cell differential expression or cluster marker genes | \nalterlab-scanpy (markers) or alterlab-scvi-tools (model-based DE) | \n
| Turning FASTQ into the count matrix (salmon/kallisto/STAR quantification, tximport) | \nalterlab-rnaseq-quant | \n
| Generic regression / GLM / mixed models on non-count data | \nalterlab-statsmodels | \n
| Somatic or germline variant calling from sequencing reads | \nalterlab-nf-core-sarek | \n
| Microbiome feature-table differential abundance | \nalterlab-qiime2-amplicon | \n
uv pip install \"pydeseq2>=0.5,<0.6\"\nSystem requirements (pydeseq2 0.5.x, current 0.5.4): Python ≥3.11; numpy ≥2.0, pandas ≥2.2,\nscipy ≥1.12, scikit-learn ≥1.4, anndata ≥0.11, formulaic ≥1.0.2 and formulaic-contrasts ≥0.2\n(parse the ~ design formula and build contrast vectors), matplotlib ≥3.9. These are pulled in\nautomatically as dependencies.
API note (0.4+): parallelism is configured through an inference object, not a bare n_cpus= kwarg:
from pydeseq2.default_inference import DefaultInference\ninference = DefaultInference(n_cpus=8)\ndds = DeseqDataSet(counts=counts_df, metadata=metadata, design=\"~condition\", inference=inference)\nds = DeseqStats(dds, contrast=[\"condition\", \"treated\", \"control\"], inference=inference)\n.T if loaded genes × samples); filter low-count genes (e.g., total reads < 10); drop samples with missing metadata.\"~condition\", \"~batch + condition\"); put adjustment variables before the variable of interest.DeseqDataSet(...).deseq2() runs the full pipeline (size factors → dispersions → LFCs → Cook's outliers).DeseqStats(dds, contrast=[var, test, ref]).summary(); read results_df.ds.lfc_shrink() for visualization/ranking only; p-values stay unshrunken.padj < 0.05, plot volcano/MA, save CSV/pickle.Minimal skeleton:
\nfrom pydeseq2.dds import DeseqDataSet\nfrom pydeseq2.ds import DeseqStats\n\ndds = DeseqDataSet(counts=counts_df, metadata=metadata, design=\"~condition\")\ndds.deseq2()\nds = DeseqStats(dds, contrast=[\"condition\", \"treated\", \"control\"])\nds.summary()\nsignificant = ds.results_df[ds.results_df.padj < 0.05]\nThis skill includes a complete standalone script for standard analyses:
\npython scripts/run_deseq2_analysis.py \\\n --counts counts.csv \\\n --metadata metadata.csv \\\n --design \"~batch + condition\" \\\n --contrast condition treated control \\\n --output results/ \\\n --min-counts 10 --alpha 0.05 --n-cpus 4 --plots\nIt handles data loading/validation, gene+sample filtering, the full DESeq2 pipeline, statistical testing with customizable parameters, result export (CSV, pickle), and optional volcano/MA plots. Refer users to scripts/run_deseq2_analysis.py for batch-processing multiple datasets.
references/pipeline_steps.md.references/interpretation_and_plots.md.references/interpretation_and_plots.md.references/api_reference.md.references/workflow_guide.md..T if needed.\"~batch + condition\").padj < 0.05 (Benjamini-Hochberg FDR), not raw p-values.[variable, test_level, reference_level]. contrast also accepts a raw\nnumpy contrast vector over the design matrix columns for comparisons a three-element list\ncannot express (e.g. interaction terms, averaging several levels).DeseqStats(..., lfc_null=1.0, alt_hypothesis=\"greaterAbs\") asks \"is |LFC| > 1?\" inside the model. That is the statistically\ncorrect way to demand an effect size — filtering a lfc_null=0 result on `abs(log2FoldChange)\n\n\n1` afterwards does not control the FDR for that claim.
\n
DeseqDataSet(..., size_factors_fit_type=\"poscounts\") uses\nthe positive-counts estimator instead of the median-of-ratios default, which fails when no gene\nis detected in every sample. control_genes= restricts size-factor estimation to spike-ins or\nhousekeeping genes.references/pipeline_steps.md — Quick-start, the six pipeline steps with full code (data prep, design, fitting, testing, shrinkage, export), and four common experimental designs.references/interpretation_and_plots.md — Filtering/ranking significant genes, quality metrics, volcano and MA plots, and a troubleshooting guide.references/api_reference.md — Complete PyDESeq2 class/method/parameter and data-structure documentation.references/workflow_guide.md — In-depth complete workflows, data-loading patterns, multi-factor designs, and best practices.