Skip to content

Add gene-annotation resource, source-qualified gene matrices, and STAR count+TPM aggregation (v0.2.2) - #18

Merged
Fuki-Kudoh merged 3 commits into
mainfrom
codex/investigate-issue-17
Jul 21, 2026
Merged

Add gene-annotation resource, source-qualified gene matrices, and STAR count+TPM aggregation (v0.2.2)#18
Fuki-Kudoh merged 3 commits into
mainfrom
codex/investigate-issue-17

Conversation

@Fuki-Kudoh

Copy link
Copy Markdown
Owner

Motivation

  • Provide gene symbols and exon-union lengths as a validated reference resource so gene-level matrices can include human-readable GeneName and enable TPM calculation from STAR counts.
  • Emit source-qualified gene-level matrices (e.g. salmon_gene_tpm.tsv) while preserving legacy v0.2.x unqualified outputs for compatibility.
  • Aggregate native STAR ReadsPerGene.out.tab across samples and compute WulfRNA-calculated gene-level TPMs using exon-union lengths and proper stranded-column selection.

Description

  • Require and validate a new reference file combined_gene_annotation.tsv and add GeneAnnotation dataclass plus parse_gene_annotation and validate_gene_resource_consistency to enforce presence, uniqueness, and positive exon-union gene_length_bp.
  • Produce source-qualified gene matrices abundance/<quantifier>_gene_expected_counts.tsv and abundance/<quantifier>_gene_tpm.tsv that include gene_id and GeneName, while still writing legacy abundance/gene_expected_counts.tsv and abundance/gene_tpm.tsv for v0.2.x compatibility.
  • Implement STAR-side handling: copy per-sample ReadsPerGene.out.tab to abundance/star_gene_counts/, validate and aggregate integer counts to abundance/star_gene_counts.tsv, and compute WulfRNA TPMs into abundance/star_gene_tpm.tsv using the stranded column selection and exon-union lengths with validation on TPM sums.
  • Update resume/manifest logic, phase output expectations, and final validation to include the new files and to fingerprint combined_gene_annotation.tsv; bump package version to 0.2.2 and update README.md, CHANGELOG.md, CITATION.cff, and docs.

Testing

  • Added unit tests tests/test_aggregation.py covering transcript aggregation, legacy vs source-qualified output, and STAR stranded-column behavior, and extended tests/test_reference_validation.py and tests/test_star_integration.py to validate the new gene-annotation handling and STAR TPM flow.
  • Ran the test suite with pytest -q which executed the new and updated tests (test_aggregation.py, test_reference_validation.py, test_star_integration.py) and all tests passed.

Codex Task

@Fuki-Kudoh

Copy link
Copy Markdown
Owner Author

Thanks for implementing the main scope of Issue #17. The overall design looks good, and the CI passes on Python 3.10 and 3.11. Before merging, I would like to address the following points.

1. Add a reproducible way to generate combined_gene_annotation.tsv

This PR makes combined_gene_annotation.tsv mandatory for all runs, including Salmon-only and kallisto-only workflows. However, there is currently no utility or documented procedure for generating it from the source GTF.

Please add a script such as:

python scripts/build_gene_annotation.py \
    --gtf combined.gtf \
    --output combined_gene_annotation.tsv

The script should:

  • read exon records from the GTF;
  • group intervals by gene_id;
  • merge overlapping exon intervals so overlapping bases are counted only once;
  • calculate the exon-union length as gene_length_bp;
  • extract gene_name as GeneName;
  • write NA when gene_name is unavailable;
  • preserve the gene IDs exactly as used by the STAR index and combined_tx2gene.tsv;
  • fail clearly on structurally inconsistent annotations, such as one gene ID occurring on multiple chromosomes.

Please also document how this file should be generated during reference preparation.

Without this addition, existing v0.2.1 reference directories cannot be upgraded reproducibly and will fail immediately under v0.2.2.

2. Expand tests for resume behavior and STAR error handling

The current tests cover the basic source-qualified matrices, annotation parsing, and STAR stranded-column selection. Please also add tests for the following important acceptance cases:

Resume and migration

  • A v0.2.1-style completed work directory reruns only aggregate because the new aggregate outputs are missing.
  • STAR alignment and Salmon/kallisto quantification are not rerun during that migration.
  • Changing combined_gene_annotation.tsv forces aggregate and downstream phases only.
  • A second unchanged invocation skips the completed aggregate phase.
  • Copied per-sample STAR files are included in aggregate completion checks.

STAR input validation

  • inconsistent gene sets between samples;
  • duplicate STAR gene IDs within a sample;
  • non-integer selected counts;
  • rows with fewer or more than four columns;
  • missing gene annotations;
  • zero total RPK;
  • multiple genes and multiple samples with known expected TPM values.

These tests are particularly important because resume compatibility is a major part of this change.

3. Validate the serialized STAR TPM values

The current TPM-sum validation checks the in-memory values immediately after normalization:

TPM = RPK / total_RPK * 1_000_000

Because of the normalization formula, this sum will almost always be one million by construction. The actual file is subsequently rounded to six decimal places, so the current check does not validate the values that users receive.

Please validate the rounded or serialized TPM values instead. Either approach is acceptable:

  • round each TPM value to six decimal places before validating the column sum; or
  • write the file, read it back, and validate the serialized values.

Please document the absolute tolerance used for this check and add a test with enough genes to exercise cumulative rounding error.

Apart from these points, the main implementation matches the intended design: source-qualified matrices, GeneName, legacy outputs, STAR count collection, strandedness-aware count selection, STAR-derived TPM, annotation fingerprinting, and aggregate output validation all look appropriate.

@Fuki-Kudoh

Copy link
Copy Markdown
Owner Author

Thanks for addressing the previous review points. The reference builder, resume/migration tests, STAR validation cases, and serialized TPM-sum validation now look good.

I found one remaining blocker in scripts/build_gene_annotation.py.

Support unquoted GTF attributes

The current attribute parser requires every GTF attribute to have a quoted value:

ATTRIBUTE_RE = re.compile(r'^\s*([^\\s]+)\\s+"([^"]*)"\s*$')

However, GENCODE GTF files commonly contain attributes with unquoted numeric values, for example:

exon_number 1;
level 2;

An otherwise valid exon record such as the following will therefore fail with Malformed GTF attribute:

chr1	HAVANA	exon	100	200	.	+	.	gene_id "ENSMUSG00000000001.1"; transcript_id "ENSMUST00000000001.1"; gene_name "Gnai3"; exon_number 1; level 2; tag "basic";

Because WulfRNA references are expected to be built from GENCODE-style annotations, this would prevent the new builder from working with the intended real-world input.

Please either:

  1. parse both quoted and unquoted GTF attribute values; or
  2. extract only the quoted attributes needed by this script, namely gene_id and gene_name, while ignoring unrelated attributes.

The second approach may be simpler and more robust because the builder does not need to validate fields such as level or exon_number.

Please also add a test containing a realistic GENCODE-style attribute field with:

exon_number 1;
level 2;

and confirm that the expected gene_id, GeneName, and exon-union length are produced.

After this parser compatibility issue is fixed and CI passes, I think PR #18 will be ready to merge.

@Fuki-Kudoh
Fuki-Kudoh marked this pull request as ready for review July 21, 2026 07:55
@Fuki-Kudoh
Fuki-Kudoh merged commit 3e0155e into main Jul 21, 2026
2 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant