Add tximport-style transcript-length normalization and pytximport support

.github/workflows/pr_validation.yml

@@ -46,8 +46,23 @@ jobs:
         run: |
           uv run ruff format --check pydeseq2
       - name: Type check with mypy
+        if: matrix.python == '3.11'
         run: |
           uv run mypy -p pydeseq2
+      - name: Test SciPy 1.13 sparse compatibility
+        if: matrix.python == '3.11'
+        run: |
+          uv run --no-project --isolated --python 3.11 \
+            --with-editable . \
+            --with 'numpy==2.0.0' \
+            --with 'scipy==1.13.0' \
+            --with 'anndata==0.11.0' \
+            --with 'pandas==2.2.2' \
+            --with 'zarr<3' \
+            --with 'pytest==8.4.0' \
+            python -m pytest -q \
+            tests/test_transcript_length_normalization.py \
+            tests/test_utils.py
       - name: Test with pytest
         run: |
           uv run coverage run -m pytest

README.md

@@ -24,6 +24,7 @@ In case there is a feature you would particularly like to be implemented, feel f
   - [Installation](#installation)
     - [Requirements](#requirements)
   - [Getting started](#getting-started)
+    - [Transcript-length normalization](#transcript-length-normalization)
     - [Documentation](#documentation)
     - [Data](#data)
   - [Contributing](#contributing)
@@ -93,6 +94,63 @@ Please don't hesitate to open an issue in case you encounter any issue due to po
 The [Getting Started](https://pydeseq2.readthedocs.io/en/latest/auto_examples/index.html) section of the documentation
 contains downloadable examples on how to use PyDESeq2.
 
+### Transcript-length normalization
+
+For unscaled estimated gene counts imported from transcript-level quantifiers
+(`countsFromAbundance="no"` in tximport), pass the matching average transcript lengths
+alongside the counts. Both matrices must use samples as rows and genes as columns. When
+using data frames, their labels and ordering must be identical; arrays must already use
+the same ordering:
+
+```python
+dds = DeseqDataSet(
+    counts=estimated_counts,
+    metadata=metadata,
+    transcript_lengths=average_transcript_lengths,
+    design="~condition",
+)
+dds.deseq2()
+```
+
+[`pytximport`](https://github.com/complextissue/pytximport), described in [4], can
+produce a compatible AnnData object directly:
+
+```python
+from pydeseq2.dds import DeseqDataSet
+from pytximport import tximport
+
+sample_names = ["sample_1", "sample_2"]
+txi = tximport(
+    [f"{sample}/quant.sf" for sample in sample_names],
+    data_type="salmon",
+    transcript_gene_map=transcript_gene_map,
+    counts_from_abundance=None,
+    output_type="anndata",
+)
+txi.obs = metadata.loc[sample_names].copy()
+dds = DeseqDataSet(adata=txi, design="~condition")
+dds.deseq2()
+```
+
+Use only unscaled estimated counts (`counts_from_abundance=None`); the file list and
+metadata rows must use the same sample order. For backed AnnData, call
+`adata.to_memory()` first and ensure sufficient memory. Length-source precedence is
+deterministic: explicit `transcript_lengths`, `adata.layers["avg_tx_length"]`, then
+compatible pytximport fields; PyDESeq2 copies `adata.obsm["length"]` into the canonical
+layer. Because that matrix has no gene labels and AnnData does not align its second
+`obsm` dimension with `adata.var`, pass it before gene-axis subsetting/reordering or
+apply the same selection and order to `adata.obsm["length"]`.
+
+Following the
+[`tximport`](https://bioconductor.org/packages/release/bioc/html/tximport.html) and
+[`DESeq2`](https://bioconductor.org/packages/release/bioc/html/DESeq2.html)
+workflow [1, 3], PyDESeq2 rounds estimated counts and combines transcript-length
+offsets with median-of-ratios library-size normalization. Average lengths and the
+resulting sample-by-gene factors are available in `dds.layers["avg_tx_length"]`
+and `dds.layers["normalization_factors"]`. The `ratio` and `poscounts`
+size-factor methods are supported; `iterative` is not currently compatible with
+transcript-length offsets.
+
 
 ### Documentation
 
@@ -159,6 +217,16 @@ In Dec 2025, the maintenance of PyDESeq2 was taken over by the scverse community
         Bioinformatics, 35(12), 2084-2092.
         <https://academic.oup.com/bioinformatics/article/35/12/2084/5159452>
 
+[3] Soneson, C., Love, M. I., & Robinson, M. D. (2015). "Differential analyses
+        for RNA-seq: transcript-level estimates improve gene-level inferences."
+        F1000Research, 4:1521.
+        <https://doi.org/10.12688/f1000research.7563.1>
+
+[4] Kuehl, M., Wong, M. N., Wanner, N., Bonn, S., & Puelles, V. G. (2024).
+        "Gene count estimation with pytximport enables reproducible analysis of
+        bulk RNA sequencing data in Python." Bioinformatics, 40(12), btae700.
+        <https://doi.org/10.1093/bioinformatics/btae700>
+
 ## License
 
 PyDESeq2 is released under an [MIT license](https://github.com/owkin/PyDESeq2/blob/main/LICENSE).

docs/source/conf.py

@@ -56,6 +56,7 @@
     "python": ("https://docs.python.org/3", None),
     "numpy": ("https://numpy.org/doc/stable/", None),
     "pandas": ("https://pandas.pydata.org/docs/", None),
+    "scipy": ("https://docs.scipy.org/doc/scipy/", None),
     "anndata": ("https://anndata.readthedocs.io/en/latest/", None),
 }