Bulk RNA-seq DE with PyDESeq2: load counts, normalize, fit negative binomial models, Wald test (BH-FDR), LFC shrinkage, volcano/MA plots.
SkillScore95Maintained
Automatic
◆Data Analyst pack
Last verified 9/3/2026
What it changes
Passed 21 of 22 test cases with the skill, and 15 of 22 without it.
Without skill68%
With skill95%
Checks it passed that the model alone missed
Calls to_picklable_anndata() on the DeseqDataSet object before pickling.
Uses Wilkinson formula syntax design="~condition" when instantiating DeseqDataSet.
Replaces design_factors with the Wilkinson formula string parameter design="~condition".
and 4 more in the full report
7 fixed14 passed either way0 still failing1 broke
Pass-rate lift
+27 pp
Tokens
+111%
Agent turns
0%
Response time
-6%
Model gemini-3.6-flashMeasured 9/3/2026v2 · latest is v321 of 22 cases passed with skill
Tokens, turns and time describe resource use, not a dollar price.
This benchmark measured an earlier version of the skill.
See it in action
One of 7 test cases the skill fixed. Same model, same prompt, with and without the skill.
1 of 7
Prompt
I am upgrading a legacy script that used PyDESeq2 0.2 with DeseqDataSet(counts, metadata, design_factors='condition'). When running in PyDESeq2 0.4+, how should the design formula parameter be updated?
What was checked, without and with the skill
Checked
Without skill
With skill
Replaces design_factors with the Wilkinson formula string parameter design="~condition".
Missed
Passed
case-18 · gemini-3.6-flash · 9/3/2026 · answers are capped and redacted
What's inside
The instructions your model reads when the skill loads.
SKILL.md393 linessingle file
PyDESeq2 Differential Expression Analysis
Overview
PyDESeq2 is a Python reimplementation of the R DESeq2 package for differential gene expression analysis from bulk RNA-seq count data. It fits negative binomial generalized linear models per gene, estimates dispersion with empirical Bayes shrinkage, and performs Wald tests with Benjamini-Hochberg FDR correction. This skill covers the full pipeline from raw counts to publication-ready result tables and visualizations.
When to Use
Identifying differentially expressed genes between two or more experimental conditions from bulk RNA-seq
Performing two-group comparisons (e.g., treated vs control) with proper statistical testing
Running multi-factor designs that account for batch effects or covariates (e.g., ~batch + condition)
Applying log2 fold change shrinkage (apeGLM) for ranking and visualization
Use omics-plotting SKILL after DE for publication-quality plots of differential expression results
Converting R-based DESeq2 workflows to a pure Python environment
Integrating DE analysis into larger Python bioinformatics pipelines (e.g., with scanpy, pandas)
Use DESeq2 (R/Bioconductor) or edgeR instead for the reference R implementations with the broadest method support and community validation
Data requirements: Raw (unnormalized) integer count matrix (samples x genes) + sample metadata DataFrame
Environment: Python 3.10+; optional matplotlib, seaborn for visualization
pip install pydeseq2 matplotlib seaborn
Pre-flight Interview
Settle these with the user before writing any analysis code.
decisions:
- id: D1
param: design
kind: required
source: data
ask: "Which column of the sample sheet separates the groups you want to compare?"
default: null