Repository navigation
Release 2026.3.0 - #111
Open
bencap wants to merge 29 commits into
Open
Release 2026.3.0#111bencap wants to merge 29 commits into
bencap wants to merge 29 commits into
Conversation
Add `vrs_utils.identify_allele` and `normalize_and_identify`, and switch `vrs_map._construct_vrs_allele` to use them in place of direct `ga4gh_identify` calls. The GA4GH Merkle-tree caches sub-object digests on the object after first identification, so any subsequent mutation (notably the pre-map `refgetAccession` swap, normalization, or state coercion) leaves a stale id unless the cached digests are cleared first. Clearing both the location and allele digests before identification ensures the id always reflects current content. - vrs_utils: new module centralizing the digest-correctness invariant - vrs_map: route both the ref-identical and SNV/delins branches through the new helpers; the ref-identical branch also gains a reassigned `normalize` return value as a side-effect - tests/test_vrs_utils: offline coverage for content-addressing, the stale-digest clearing invariant, the normalize+identify pairing, and the malformed-input error path - align: incidental formatter-driven whitespace tweak
Add MappingOutcome so every (variant, level) record carries a typed outcome, distinguishing a benign absence (intronic, no protein consequence) from a genuine failure that error_message alone cannot convey. Derive the preferred (authoritative) layer from the target's assay level rather than always preferring genomic.
…c projection Select a coding transcript for NC_ protein-coding targets via Ensembl locus overlap plus MANE, so they are no longer silently skipped by reverse translation. Project each measured variant onto its deterministically reachable layers (g<->c, nucleotide->p), emitted as typed-outcome records and routed by preferred_layer_only.
…ferred-layer records - Track represented variant IDs at the preferred layer; re-attribute only variants that have no preferred-layer record, avoiding duplicate mapped_scores for variants with both a dead genomic attempt and a measured protein record (e.g. codon-optimised targets). - Synthesize a preferred-layer failure for variants that mapped only at a non-preferred layer (e.g. wild-type p.= on a genomic-preferred target) so every input variant always has exactly one output record. - Extract _map_protein_layer in vrs_map to return (mapping, reason) instead of an ad-hoc error MappedScore; a row that maps at no layer is failed once, layer-agnostically, carrying the detailed reason. - Add TestNullFailureDedup and TestMapProteinLayerReason test coverage.
…ignments
For protein-vs-DNA BLAT alignments, qcoords always increase (protein
reads N→C regardless of genome strand), so they cannot be used to
detect strand. Switch to tcoords direction for protein queries.
Also normalise hit_subranges and hit_range entries with min/max so
they are always in ascending order, which they are not when the target
gene sits on the minus strand.
- Use tcoords direction (not qcoords) for strand detection when
-q=prot is in blat_params
- Wrap hit_subrange and hit_range endpoints in min/max in both
_get_best_match and align_target_to_protein
_get_mapped_reference_sequence had no CDNA branch, so it fell through to the genomic chromosome lookup and returned the NC_ accession as the post_mapped reference for the cdna layer. This caused target_genes post_mapped_metadata to carry NC_000017.11 under the "cdna" key instead of the NM transcript. Add a CDNA branch that resolves the NM/ENST accession from tx_output.nm (preferred, covers NC_/sequence-based targets) or from the target's own accession when it is already an NM_/ENST (cdna-source targets). Returns None rather than a chromosome when no NM is resolvable. Also adds unit tests for all three layer paths.
…port detailed mapping results
There is no underscore in Ensembl IDs so removing it here.
…rotein-id-check fix: Remove underscore from check for Ensembl Protein ID
…ness feat: VRS Correctness
…ant-projection feat: Reverse Translation Mapping Support
- Reindent the all([...]) call in _offset_target_sequence to hug the brackets instead of nesting the list on its own lines - Bind target = scoreset_metadata.target_genes[target_gene] in select_transcripts to shorten downstream accesses
…unterparts Accession-based targets declared with Ensembl protein (ENSP) or transcript (ENST) accessions previously passed through as bare, unresolved accessions or were dropped, leaving downstream annotation inconsistent with the RefSeq accessions used elsewhere in the pipeline. Resolving them to their RefSeq MANE counterpart keeps output consistent regardless of how the target was declared. - add get_gene_symbol_from_ensembl_protein / _transcript lookups as the bridge from a non-RefSeq accession to get_mane_transcripts_for_gene - add _select_refseq_protein_counterpart / _select_refseq_cdna_counterpart, falling back to the bare accession (protein) or None (cDNA, since np is required on TxSelectResult) when no counterpart resolves - skip standalone VRS mapping for a target with no record group of its own (e.g. base-editor score sets where every row groups under the cDNA target) instead of KeyError'ing on records[target_gene] - re-raise the original alignment error instead of losing its traceback via `raise AlignmentError from e`
identify_allele cleared the cached digests but called ga4gh_identify with the default in_place, which only fills an empty id — and every allele reaching it already has one, minted by AlleleTranslator over un-normalized content. Force recomputation with in_place="always", and recompute location.id separately, since ga4gh_identify writes only the id of the object it is handed. Only deletions and duplications in repeat regions were affected; substitutions normalize to themselves, so their pre-normalization id stayed correct by coincidence and hid the bug. Also look up vrs_ref_allele_seq by extension name rather than extensions[0], which was correct only while it was the sole extension, and compute it only for VRS 1.3 VariationDescriptors — VRS 2.x recovers the reference from refgetAccession and coordinates, so computing it for 2.0 alleles cost a seqrepo fetch per allele for a field nothing reads.
|
coverage: 0.0%. remained the same — mavedb-dev into mavedb-main |
…red writes The mapper built accession alleles on whatever sequence cdot assembled from the genome, and stored it under the accession alias inside a transaction it never committed. For NM_007294.3 that sequence (7207 nt) is not NCBI's record (7224 nt), so the mapper minted a different refget, and so a different allele digest, than reverse translation for the same variant, which left duplicate alleles and broke the convergent/projection labelling. - seqfetcher() now puts hgvs's SeqRepo-backed SeqFetcher ahead of the FASTA fallback, matching the API's chain. - RefSeq accessions are read from SeqRepo and never written by a mapping run; a missing one raises ReferenceSequenceNotFoundError. resolve_refget() requires exactly one refget for an accession instead of taking the first of a set. - A missing Ensembl transcript is assembled from each genome build separately, refused unless the builds agree, and stored in SeqRepo without reassigning an existing alias. - Every mapper write to SeqRepo commits immediately and the SQLite busy timeout is raised, so parallel jobs no longer queue behind one process's open transaction. - Golden tests pin the accession -> refget contract against NCBI records held in a fixture that the API repo also carries.
On the pinned ga4gh.vrs 2.0.0a6, Extension lives in ga4gh.core, so test_annotate.py failed at collection. CI only runs test_mavedb_data.py, which hid the break.
The FastAPI server in src/api is network-facing; the floor keeps a resolve from landing on starlette releases with known advisories.
hgdownload.cse.ucsc.edu is a legacy alias whose certificate only covers hgdownload.soe.ucsc.edu. Download blat and hg38.2bit from the soe host over HTTPS, fail on HTTP errors, and check hg38.2bit against the MD5 published in UCSC's md5sum.txt.
The workflow publishes to a PyPI project this repository does not own and has failed on every release.
…rgets Accession targets declared on an NM_ (or an ENST with no RefSeq MANE counterpart) have no transcript selection, so the genomic-layer cdna entry in reference_sequences was never emitted. Downstream consumers reading post-mapped metadata, such as MaveDB's UniProt ID mapping, then found no coding transcript and skipped the target. Resolve the entry's accession with _resolve_cdna_accession, the same resolver used for the cdna identity target_mapping, so the two cannot disagree. Do not overwrite an existing cdna layer.
…he declared protein A declared NM_/ENST accession was swapped for its gene's MANE Select transcript without checking the protein, so a non-MANE isoform could be mapped against a different protein. The swap now requires the MANE transcript's protein to match the declared one (looked up via cdot for cDNA accessions) and otherwise keeps the declared accession with a warning. Accession mapping skips the protein alignment adjustment.
…inutes A protein query against the translated genome took longer than the fixed 600 s when several mapping jobs shared the host. BLAT_TIMEOUT_SECONDS sets the limit, defaulting to 1200.
…king by proximity The Bio.Align port compared the gene start against an alignment's first target coordinate, which on the minus strand is its upper bound. For SCN5A that chose an 84-residue fragment over the full-length hit, and transcript selection then failed. Also log when another hit overlapping the gene aligns more of the query than the proximity pick.
SeqRepo hashes a sequence after stripping '*' and whitespace and uppercasing, but store_sequence hashed the raw string. A protein target ending in a stop '*' was stored under the stop-free digest with a second ga4gh alias naming the raw digest, and refget lookups for that sequence then returned either.
_get_uniprot_ref always read targetGenes[0], so every target in a multi-target score set inherited the first target's UniProt ID. Gene normalization tries that ID first, so in urn:mavedb:00001275-0-2 BRAF resolved to the KRAS locus, its alignment landed on a chr12 fragment with no transcripts, and transcript selection fell through to the UniProt fallback with KRAS's accession.
protein_sequence was only bound on the transcript branch of _select_protein_reference, so a target with no compatible transcripts raised UnboundLocalError instead of falling back to its UniProt reference. Neither select_transcripts nor the mapping router catches that, so the whole score set failed.
This branch has not been deployed
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
3 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Features
MappingOutcomeschema definition and anoutcomefield for per-variant mapping resultsBug Fixes
mapped_scoresrecords at the preferred layerMaintenance
.claudeto.gitignore