Add tximport-style transcript-length normalization and pytximport support by gitbenlewis · Pull Request #451 · scverse/PyDESeq2

gitbenlewis · 2026-07-16T21:36:04Z

Reference Issue or PRs

Closes #305.
Closes #359.
Related to #412.

What does your PR implement? Be specific.

Adds opt-in DESeq2-style transcript-length normalization for dense estimated gene counts.

Pass average transcript lengths explicitly:

dds = DeseqDataSet(
    counts=estimated_counts,
    metadata=metadata,
    transcript_lengths=average_transcript_lengths,
    design="~condition",
)

For AnnData from pytximport:

dds = DeseqDataSet(
    adata=txi,
    transcript_lengths=txi.obsm["length"],
    design="~condition",
)

The implementation rounds estimated counts and combines transcript-length offsets with library-size normalization. Gene-specific factors are used throughout dispersion and LFC fitting, Wald tests, shrinkage, Cook’s outlier replacement/refitting, and VST. Both ratio and poscounts estimation are supported.

Correction requires explicit transcript_lengths; stored metadata does not activate it. Explicit correction of abundance-scaled counts is rejected, and unscaled pytximport input without opt-in produces a warning.

Usage examples are included in the documentation.

Validation

43 transcript-length normalization tests pass locally.
Full local suite: 108 tests pass, including two download tests rerun with network access.
Committed DESeq2-generated fixtures check ratio and poscounts factors and normalized counts using the existing synthetic dataset. A small R regeneration script is included.
Additional manual DESeq2 comparisons checked zero counts, control genes, and downstream fitting on the synthetic dataset.
Lint, formatting, package build, and documentation checks pass.

Local tests used Python 3.12 and pandas 3.0.5. Current GitHub test CI with pandas 3.1.1 fails on deprecation warnings involving AnnData’s future.infer_string access and an existing rename(inplace=True) test call. These compatibility fixes are not included.

Current limitations

Iterative size-factor estimation is not supported with length offsets.
VST transformation of external counts with separate lengths is not supported.
pytximport’s obsm["length"] gene dimension must remain synchronized when genes are subset or reordered.
Committed R reference fixtures cover normalization; downstream tests do not yet provide full end-to-end R reference coverage.

Related work and acknowledgements

This work began with the generic transcript_lengths API and includes explicit pytximport input support following discussion in #412.

Thank you to @maltekuehl for the earlier implementation and design work in #412, and to @Zethson for encouraging this work. The tximport, DESeq2, and pytximport work is credited in the documentation.

gitbenlewis · 2026-07-29T21:12:10Z

Hi @BorisMuzellec — just checking in on PR #451 when you have a chance. All current checks are passing. Could you share an update on the review, or let me know if there’s anything else you’d like me to address? Thanks!
@grst

sriramdevadas · 2026-09-20T22:46:35Z

@Zethson, would it be useful if I verified this against R's DESeqDataSetFromTximport on a public salmon dataset and contributed the R reference outputs as fixtures in the repo's existing test format? Glad to do it, or to stay out of the way if you have it covered. (cc @gitbenlewis)

gitbenlewis · 2026-09-21T20:32:35Z

@Zethson, would it be useful if I verified this against R's DESeqDataSetFromTximport on a public salmon dataset and contributed the R reference outputs as fixtures in the repo's existing test format? Glad to do it, or to stay out of the way if you have it covered. (cc @gitbenlewis)

Thanks! I have a separate PyDESeq2–R DESeq2 parity-testing repository, where I’ve added real-Salmon coverage using six public GEUVADIS samples imported independently through R tximport and pytximport. Import and normalization checks pass, and both PyDESeq2 input interfaces agree.
There are some remaining numerical parity differences that appear unrelated to this PR’s import and normalization changes, with details in the diagnosis document. I’ve started addressing these on a separate dispersion-parity branch. It contains partial improvements; full parity is not yet achieved. I suggest keeping that work in a separate PR.
The R reference outputs are now included as fixtures in the external repository. An independent check and help integrating them into PyDESeq2’s existing test format would still be very welcome.

@Zethson merge conflicts resolved : )

…layout

codecov-commenter · 2026-09-21T21:48:15Z

Codecov Report

❌ Patch coverage is 99.19028% with 2 lines in your changes missing coverage. Please review.
✅ Project coverage is 89.77%. Comparing base (2ca3b74) to head (c742803).
⚠️ Report is 1 commits behind head on main.

Files with missing lines	Patch %	Lines
src/pydeseq2/dds.py	99.01%	2 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #451      +/-   ##
==========================================
+ Coverage   86.05%   89.77%   +3.72%     
==========================================
  Files          15       15              
  Lines        1276     1486     +210     
==========================================
+ Hits         1098     1334     +236     
+ Misses        178      152      -26

Files with missing lines	Coverage Δ
src/pydeseq2/default_inference.py	`97.82% <100.00%> (+1.72%)`	⬆️
src/pydeseq2/dispersions.py	`96.93% <100.00%> (+0.09%)`	⬆️
src/pydeseq2/ds.py	`90.14% <100.00%> (+0.14%)`	⬆️
src/pydeseq2/inference.py	`100.00% <100.00%> (ø)`
src/pydeseq2/misc.py	`90.90% <100.00%> (+4.24%)`	⬆️
src/pydeseq2/dds.py	`92.71% <99.01%> (+5.09%)`	⬆️

... and 1 file with indirect coverage changes

🚀 New features to boost your workflow:

❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

…y coverage

gitbenlewis · 2026-11-03T20:05:48Z

@BorisMuzellec, could you review this when you have a chance? I’ve resolved the conflicts with main and fixed the CI failures. @Zethson, if another reviewer would be better suited, could you help route this? Thanks!

Zethson · 2026-11-03T21:08:24Z

Sorry, @BorisMuzellec is the right person to review this. I can see @maltekuehl potentially help out - I think it's potentially also in your interest?

maltekuehl

First of all, thank you for working on this. I have left comments across the files.

As with many PRs in the LLM-age, this one seems to be mixing many different bug fixes, features and other changes at the same time, making review somewhat difficult. > 1000 LOC tests don't help that, either.

Could you please split it up into different PRs?

a minimal implementation of opt-in transcript-length normalization
sparse data support
unrelated bug fixes

When it comes to tests and the preserve fit part, I believe sizeable chunks can probably also be dropped entirely, unless there are good reasons to keep them. Same goes for docstring or README re-formatting, please keep the diff to what is actually needed at a time.

Happy to review a more focussed PR again!

maltekuehl · 2026-11-04T19:57:27Z

+          uv run --no-project --isolated --python 3.12 \
+            --with-editable . \
+            --with 'numpy==2.2.0' \
+            --with 'scipy==1.15.0' \
+            --with 'anndata==0.11.0' \
+            --with 'pandas==2.3.0' \
+            --with 'zarr<3' \
+            --with 'pytest==8.4.0' \
+            --with 'coverage>=7.10' \


Is there a specific reason for these overrides? We should aim to follow SPEC0 with PyDESeq2 and 3.12 is outside the inclusion window now.

A dependency and Python version update should perhaps be a separate PR to align with SPEC0.

Removed the overrides and restored the upstream test workflow. Python and dependency-version changes are outside the scope of this PR.

A dependency and Python version update should perhaps be a separate PR to align with SPEC0.

Agreed. I’ve removed those changes from this PR so that any version-policy update can be reviewed separately.

maltekuehl · 2026-11-04T19:59:19Z

 warnings.simplefilter("ignore", FutureWarning)


+def _rounded_counts(counts: Any) -> Any:


I would suggest moving this out of dds.py to misc.py since it is not dataset related. Also, is there a good reason to use Any here instead of a union type?

With sparse support removed, the helper was unnecessary. I removed it and inlined the two calls using counts.round() and np.rint().

maltekuehl · 2026-11-04T20:01:10Z

    design
-        Model design.
-        Can be either a pandas DataFrame representing a design matrix, or a formulaic formula in the format ``'x + z'`` or ``'~x+z'``.
-        If a design matrix is provided, DeseqStats built from this DeseqDataSet will only support contrasts in the form of numeric vectors.
+        Model design. Can be either a pandas DataFrame representing a design matrix, or
+        a formulaic formula in the format ``'x + z'`` or ``'~x+z'``.
+        If a design matrix is provided, DeseqStats built from this DeseqDataSet will
+        only support contrasts in the form of numeric vectors.


Could we keep out purely cosmetic docstring changes for the initial PR to get a better understanding of what actually changed? Or is this a mismatch in the linter/formatter?

Reverted the purely cosmetic docstring changes in c742803

maltekuehl · 2026-11-04T20:04:24Z

+        assert not hasattr(dds, attr)
+
+
+def test_transcript_length_normalization_matches_deseq2_formula():


This test feels circular and repetitive, no? Instead of proving a correct implementation of the original formula by repeating the formula in the test, I would suggest comparing against a fixture generated with DESeq2 to prove concordance and to prevent drift. E.g., some of the typical GEUVADIS test data or a subsample thereof. This should also eliminate the need for some of the other functions in this test file and reduce LOCs.

Replaced the formula-based expectations with DESeq2-generated fixtures for ratio and poscounts normalization factors and normalized counts in edd40f7. These reuse the repository’s existing synthetic dataset, with a small R regeneration script included.

maltekuehl · 2026-11-04T20:06:41Z

+1187 on a single test file seems excessive. Are there things you could use @pytest.mark.parametrize for? And instead of testing correctness of every single part, it is often easier to test end-to-end against a fixture to detect drift.

Especially when working with LLMs, please provide strict guardrails on what needs a test and how. Integration tests, correctness tests and end-to-end tests are much preferable over implementation and unit tests, generally.

Reduced the file from 1,187 to 493 lines across bb8b50d, 0e9dfc4, edd40f7, and b162f02 by removing sparse tests, parametrizing repeated cases, and consolidating pipeline checks. The independent R fixtures cover normalization; downstream checks remain Python integration tests rather than full end-to-end R comparisons.

removed redundant assertions and consolidated the integration checks in b162f02, alongside the independent normalization fixtures added in edd40f7. The remaining moments test uses small explicit inputs and hand-calculated expectations rather than repeating the implementation.

maltekuehl · 2026-11-04T20:18:12Z

+from scipy.sparse import sparray
+from scipy.sparse import spmatrix
+
+type CountMatrix = np.ndarray | spmatrix | sparray


Perhaps we should consider moving types to a _types file that can be imported from so that we don't run the risk of circular imports between the consumers of shared types in the future.

Removed the added shared matrix types along with sparse support in dfd7217, so those cross-module type dependencies are no longer part of this PR.

maltekuehl · 2026-11-04T20:26:33Z

+            )
+            has_transcript_lengths = has_transcript_lengths or has_pytximport_fields
+            counts_from_abundance = adata.uns.get("counts_from_abundance")
+            if counts_from_abundance is not None and has_transcript_lengths:


Using pytximport with counts_from_abundance="length_scaled_tpm" is perfectly fine to use with PyDESeq2 out-of-the-box and it carries .obsm["length"], so raising here breaks existing, desired workflows. It should be up to the user which workflow they want to use, we should not force them into the DESeqDataSetFromTximport workflow.

Without an opt-in argument, we would also silently change the output for users where counts_from_abundance is None. In that case, I think a warning would be fair though, since using output from that configuration just like that with PyDESeq2 is not recommended.

Length correction now requires an explicit transcript_lengths= argument following 0e9dfc4. Stored length metadata alone no longer activates correction or causes scaled inputs to be rejected. The error applies only when explicit correction is combined with abundance-scaled counts; ordinary integer-count validation remains unchanged.

Added explicit opt-in in 0e9dfc4 and the suggested warning in b162f02. When counts_from_abundance is explicitly None and no lengths are supplied, the constructor warns without enabling correction. Inputs without that metadata key and abundance-scaled inputs do not receive the warning.

maltekuehl · 2026-11-04T20:32:11Z

+            if adata.isbacked:
+                raise ValueError(
+                    "DeseqDataSet requires an in-memory AnnData object. Call "
+                    "adata.to_memory() before initialization."


Is there a reason to have that here and just here and to include it in this PR? Is this a condition specific to functionality introduced in this PR or just something general?

This was a general restriction rather than a requirement of length normalization. Removed the guard and its dedicated rejection test in b162f02. The guide retains the recommendation to load backed data into memory first.

maltekuehl · 2026-11-04T20:33:43Z

+            tests/test_transcript_length_normalization.py \
+            tests/test_stats.py


Why is it running just these files?

Removed the file selection and restored the upstream workflow in dfd7217, including its normal test-suite invocation.

maltekuehl · 2026-11-04T20:36:56Z

+        if preserve_inherited_fit and adata is not None:
+            source_design = adata.obsm.get("design_matrix")
+            normalization_control_mask = self._make_control_mask(control_genes)
+            fit_settings = (
+                ("fit_type", fit_type),
+                ("min_mu", min_mu),
+                ("min_disp", min_disp),
+                ("max_disp", np.maximum(max_disp, self.n_obs)),
+                ("refit_cooks", refit_cooks),
+                ("min_replicates", min_replicates),
+                ("beta_tol", beta_tol),
+            )
+            preserve_inherited_fit = (
+                isinstance(source_design, pd.DataFrame)


Could you explain this entire block up to L505 here, how it relates to the PR, why it's needed and the reasons behind the implementation? Without diving deep into DESeq2 again, I cannot follow what the purpose is.

Removed the inherited-fit reuse logic in 0e9dfc4. The remaining cleanup prevents copied size factors, dispersions, and LFCs from surviving a new normalization choice. It preserves user metadata while clearing fit-derived results, so subsequent fitting cannot silently reuse stale values. Tests cover switching between ordinary and length normalization, including construction from exported AnnData.

gitbenlewis · 2026-11-05T16:28:57Z

First of all, thank you for working on this. I have left comments across the files.

As with many PRs in the LLM-age, this one seems to be mixing many different bug fixes, features and other changes at the same time, making review somewhat difficult. > 1000 LOC tests don't help that, either.

Could you please split it up into different PRs?

a minimal implementation of opt-in transcript-length normalization

sparse data support

unrelated bug fixes

When it comes to tests and the preserve fit part, I believe sizeable chunks can probably also be dropped entirely, unless there are good reasons to keep them. Same goes for docstring or README re-formatting, please keep the diff to what is actually needed at a time.

Happy to review a more focussed PR again!

@maltekuehl . thank you for taking the time to review this.

I appreciate the detailed feedback and will carefully work through each comment.

I agree that the PR has grown beyond its original scope and become difficult to review. I'll start working on spliting this up into the 3 separate PRs you suggested.

Thanks agian for offering to review a more focused version!

Restore upstream count handling, CI configuration, and compatibility behavior while retaining dense transcript-length normalization. Sparse test cleanup and explicit opt-in behavior remain pending.

Remove sparse-only tests and helpers, and use dense inputs in mixed normalization tests. Drop the assertion tied to the reverted VST fix.

Require transcript_lengths to enable length correction. Remove automatic metadata detection and inherited-fit reuse. Update tests and documentation for explicit input and fresh fitting.

Replace ratio and poscounts formula tests with R-generated references using upstream synthetic data. Include a small regeneration script.

Inline rounding, remove the backed-input guard, and consolodate tests. Warn on unscaled input without opt-in and move the guide into docs

gitbenlewis added 2 commits July 14, 2026 00:56

Add transcript-length normalization support

3e9eb78

Add pytximport AnnData compatibility

5096277

gitbenlewis requested a review from BorisMuzellec as a code owner July 16, 2026 21:36

gitbenlewis mentioned this pull request Jul 16, 2026

Normalization factors and DESeqDataSetFromTximport equivalent #412

Draft

gitbenlewis added 6 commits July 22, 2026 19:48

CI: run mypy only on Python 3.11

384c0d3

Docs: link SciPy sparse matrix type

68453f9

Fix sparse count normalization and fit state

3094471

Fix CountMatrix type alias

366cf54

Fix sparse count handling and inherited fit state

6270295

Reduce sparse support code duplication

0e56de2

Merge upstream main and adapt transcript-length normalization to new …

1281923

…layout

gitbenlewis added 4 commits September 21, 2026 22:07

Cover normalization validation errors and include sparse compatibilit…

9b42a30

…y coverage

Merge branch 'main' into main

1a12926

Fix pandas 3.1 deprecation failures in CI

322993b

Merge upstream main and resolve documentation configuration conflict

e3ce9b0

maltekuehl self-requested a review November 4, 2026 19:55

maltekuehl suggested changes Nov 4, 2026

View reviewed changes

gitbenlewis added 6 commits November 5, 2026 17:58

revert cosmetic docstring changes

c742803

remove sparse support and unrelated fixes from normalization PR

dfd7217

Restore upstream count handling, CI configuration, and compatibility behavior while retaining dense transcript-length normalization. Sparse test cleanup and explicit opt-in behavior remain pending.

remove sparse tests from transcript-length normalization suite

bb8b50d

Remove sparse-only tests and helpers, and use dense inputs in mixed normalization tests. Drop the assertion tied to the reverted VST fix.

make transcript-length normalization explicitly opt-in

0e9dfc4

Require transcript_lengths to enable length correction. Remove automatic metadata detection and inherited-fit reuse. Update tests and documentation for explicit input and fresh fitting.

test transcript-length normalization against DESeq2 references

edd40f7

Replace ratio and poscounts formula tests with R-generated references using upstream synthetic data. Include a small regeneration script.

simplify code and tests

b162f02

Inline rounding, remove the backed-input guard, and consolodate tests. Warn on unscaled input without opt-in and move the guide into docs

gitbenlewis requested a review from maltekuehl November 5, 2026 23:28

This was referenced Nov 7, 2026

[auto-optimize]: Add unit tests for the DESeq2 pipeline and its solvers #492

Open

[auto-optimize]: Speed up deseq2, summary and lfc_shrink #493

Open

		warnings.simplefilter("ignore", FutureWarning)


		def _rounded_counts(counts: Any) -> Any:

		assert not hasattr(dds, attr)


		def test_transcript_length_normalization_matches_deseq2_formula():

		tests/test_transcript_length_normalization.py \
		tests/test_stats.py

Conversation

gitbenlewis commented Jul 16, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Reference Issue or PRs

What does your PR implement? Be specific.

Validation

Current limitations

Related work and acknowledgements

Uh oh!

gitbenlewis commented Jul 29, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

sriramdevadas commented Sep 20, 2026

Uh oh!

gitbenlewis commented Sep 21, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

codecov-commenter commented Sep 21, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Codecov Report

Uh oh!

gitbenlewis commented Nov 3, 2026

Uh oh!

Zethson commented Nov 3, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

maltekuehl left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

gitbenlewis commented Jul 16, 2026 •

edited

Loading

gitbenlewis commented Jul 29, 2026 •

edited

Loading

gitbenlewis commented Sep 21, 2026 •

edited

Loading

codecov-commenter commented Sep 21, 2026 •

edited

Loading

Zethson commented Nov 3, 2026 •

edited

Loading