pyedger

A pure-Python port of Bioconductor edgeR (Robinson, McCarthy & Smyth, Bioinformatics 2010) — negative-binomial models for differential expression of count data.

• No rpy2, no R install — the edgeR negative-binomial GLM workflow reimplemented in NumPy / SciPy, vectorised over genes
• Matches Bioconductor edgeR 4.10 numerically (see Numerical parity)
• The canonical pipelines: DGEList → calcNormFactors → estimateDisp → glmQLFit → glmQLFTest and the classic exactTest
• TMM normalization, common / trended / tagwise dispersion, GLM and quasi-likelihood F-tests, filterByExpr, cpm / aveLogCPM
• Both Python-style (glm_fit, estimate_disp, top_tags) and R-style (glmFit, estimateDisp, topTags) names exported

This is a standalone mirror of the implementation developed in omicverse, where it powers the edgeR differential-expression backend of ov.bulk / pyDEG.

Install

pip install pyedger

Quick start

import numpy as np
import pyedger

# counts: genes x samples raw count matrix; group: per-sample condition labels
dge = pyedger.DGEList(counts=counts, group=group)
dge = pyedger.calcNormFactors(dge)                       # TMM normalization
keep = pyedger.filterByExpr(dge, group=group)
dge = dge[keep]

# Quasi-likelihood F-test workflow (the recommended edgeR pipeline)
dge = pyedger.estimateDisp(dge, design)
fit = pyedger.glmQLFit(dge, design)
qlf = pyedger.glmQLFTest(fit, coef=1)
res = pyedger.topTags(qlf, n=np.inf)
res.head()

Classic exact test

dge = pyedger.estimateDisp(dge, design)
et = pyedger.exactTest(dge, pair=("control", "treated"))
pyedger.topTags(et)

API

Numerical parity with edgeR

tests/test_r_parity.py checks pyedger against edgeR 4.10.5 (R 4.5.3, limma 3.68.5, locfit 1.5-9.12) on a paired design with zero-count library pairs, a Milo-style design and limma hyperparameter cases; the reference values are regenerated by tests/data/make_reference.R. The GLM kernels follow edgeR's compiled code step by step and run vectorised over all genes:

glmQLFit(legacy=False) without a dispersion and estimateGLMCommonDisp estimate the NB dispersion with edgeR's optimize(tol=1e-5) over dispersion^(1/4); R's own value is only defined to ~4e-5 relative, and agreement there is at that level (p-values ~1e-4 relative on the test data). On real data the agreement is usually much closer: on four GSE274113 pseudobulk ATAC comparisons (~170k peaks x 26 samples, ~ library + condition, default new QL) p-values agreed to < 1e-7 relative with identical FDR < 0.1 calls.

Defaults follow edgeR 4.10: glmQLFit(legacy=False, abundance_trend=True), estimateDisp(prior_df=None) with the chooseLowessSpan span.

Citation

Robinson, M.D., McCarthy, D.J., Smyth, G.K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26(1), 139–140 (2010).

…and acknowledge omicverse / this repo for the Python port.

License

LGPL-3.0-or-later — matches the upstream Bioconductor package.