edge-rs 0.0.1

Negative binomial differential expression for bulk and single-cell RNA-seq: the edgeR numerical stack and NEBULA mixed models.
Documentation
[![CI](https://github.com/GregorLueg/edge-rs/actions/workflows/test.yml/badge.svg)](https://github.com/GregorLueg/edge-rs/actions/workflows/test.yml)
[![Crates.io Version](https://img.shields.io/crates/v/edge-rs.svg)](https://crates.io/crates/edge-rs)
[![docs.rs](https://img.shields.io/docsrs/edge-rs)](https://docs.rs/edge-rs)
[![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://opensource.org/licenses/MIT)

# edge-rs

## Description

Negative binomial differential expression for bulk and single-cell RNA-seq. This
is the edgeR numerical stack in Rust: normalisation, Cox-Reid dispersion
estimation, the Levenberg-damped NB GLM, quasi-likelihood weights, the exact
test and `diffSpliceDGE`. On top of that sit the limma routines edgeR leans on
(`squeezeVar`, the F-distribution fits, lowess and locfit smoothing,
`voomLmFit`) and NEBULA, a negative binomial gamma mixed model for single cell.

No R, no Python, no BLAS to hunt down. CPU only.

## Install

```sh
cargo add edge-rs
```

SIMD dispatch happens at runtime, so a stock build already picks up NEON on
aarch64 and AVX2 on x86. The 512-bit kernels are the exception: they are gated
at compile time and need

```sh
RUSTFLAGS="-C target-cpu=x86-64-v4" cargo build --release
```

A bare `+avx512f` is not enough, since `wide` only ships its real 512-bit types
under the full feature set.

## Quick start

The classic GLM chain: filter, normalise, estimate dispersions, fit, test, rank.
Counts are dense and gene-major throughout, `n_genes` rows of `n_samples`
values, row-major, so one gene is a contiguous slice.

```rust
use edge_rs::core::dgelist::DgeList;
use edge_rs::core::expression::ave_log_cpm;
use edge_rs::core::filtering::filter_by_expr;
use edge_rs::core::normalisation::{NormMethod, calc_norm_factors};
use edge_rs::dispersion::estimate::estimate_disp;
use edge_rs::glm::fit::glm_fit;
use edge_rs::glm::test::{GlmTestInput, Tested, glm_lrt};
use edge_rs::prelude::*;
use edge_rs::results::{SortBy, top_tags};

let dge = DgeList::new(counts, n_genes, n_samples, Some(group))?;

let keep = filter_by_expr(
    &dge.counts, dge.n_genes, dge.n_samples,
    None, dge.group.as_deref(), None, None,
)?;
let mut dge = dge.subset_genes(&keep)?;

dge.norm_factors = calc_norm_factors(
    &dge.counts, dge.n_genes, dge.n_samples,
    None, NormMethod::Tmm, None, None,
)?;

// Log of lib.size * norm.factors, one value per sample.
let offset = dge.offset()?;
let abundance = ave_log_cpm(
    &dge.counts, dge.n_genes, dge.n_samples,
    None, Some(&offset), 2.0, None,
)?;

let disp = estimate_disp(
    &dge.counts, dge.n_genes, dge.n_samples,
    &design, n_coef, &offset, None, Some(&abundance), None,
)?;
let dispersion = Recycled::by_gene(disp.tagwise.expect("tagwise is on by default"));

let fit = glm_fit(
    &dge.counts, dge.n_genes, dge.n_samples,
    &design, n_coef, &dispersion, &offset, None, 0.125,
)?;

let input = GlmTestInput {
    counts: &dge.counts,
    n_genes: dge.n_genes,
    n_samples: dge.n_samples,
    design: &design,
    n_coef,
    dispersion: &dispersion,
    offset: &offset,
    weights: None,
    log_cpm: Some(&abundance),
};
let lrt = glm_lrt(&input, &fit, &Tested::Coef(vec![n_coef - 1]), None)?;

let top = top_tags(
    &lrt.log_fc,
    lrt.log_cpm.as_deref().unwrap(),
    Some(&lrt.statistic),
    &lrt.p_value,
    10,
    SortBy::PValue,
    1.0,
)?;
```

Swap `glm_fit` for `glm_ql_fit` and `glm_lrt` for `glm_ql_ftest` to get the
quasi-likelihood pipeline instead. `exact_test` covers the pre-GLM two-group
path, `voom_lmfit` the limma route, and `nebula` the single-cell one.

Fits and containers are generic over `EdgeFloat`, so single-cell counts can be
held as `f32` and halve the memory. Likelihoods, Cox-Reid determinants, the
optimisers and every p-value run in `f64` regardless: an NB deviance is a
difference of large logs and loses edgeR parity in `f32` long before it saves
anything.

## What is not in here

This is not a general limma port. Only the pieces `estimateDisp`, `glmQLFit` and
`voomLmFit` actually reach for are implemented, so the voom chain stops at the
moderated variances. There is no `eBayes` and no limma `topTable` (but might be
implemented in the future...). On the single-cell side, NEBULA is `NBGMM` only:
`PMM` needs Poisson-gamma kernels this crate does not have, and `NBLMM` has no
golden to validate against.

## Parity

The end-to-end suites gate against committed CSV fixtures generated by
`tests/r/generate_fixtures.R` against R `4.5.1` with edgeR `4.8.2`, limma
`3.66.0`, statmod `1.5.2` and nebula `1.5.8`. CI never runs the R; the fixtures
are in the repository and the Rust reads them.

Tolerances are measured rather than guessed. Every one carries the worst error
actually observed, and the test harness prints that table on demand:

```sh
EDGE_RS_TOL_REPORT=1 cargo test --release --test 'e2e_*' -- --nocapture
```

## Deviations

The port's working reference was
[edgePython](https://github.com/pachterlab/edgePython), itself a port of edgeR,
limma and nebula. Where it disagreed with the R and C++ it came from, upstream
won: edgeR and limma are what users compare against, and their output is what
the test suite gates on. NEBULA is ported from the nebula package's own C++
rather than from the Python.

Rarer, and worth stating plainly: a handful of places where edgeR or limma are
themselves wrong, each verified by reproducing the fault in the installed
package. Reproducing a bug faithfully is not parity worth having.

Both sets are written up in [docs/UPSTREAM_DEVIATIONS.md](docs/UPSTREAM_DEVIATIONS.md),
with the upstream file and line, what this crate does instead, and a test that
would fail if the behaviour drifted back. The edgePython entries were read
against version 0.2.6, commit `1e572ae`, and any of them may well have been
fixed upstream since.

## Licence

MIT License

Copyright (c) 2026 Gregor Lueg

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.