Skip to content
22 changes: 22 additions & 0 deletions schema.json
Original file line number Diff line number Diff line change
Expand Up @@ -95,6 +95,17 @@
"title": "MappedReferenceSequence",
"type": "object"
},
"MappingOutcome": {
"description": "Per-record outcome for one (variant, annotation level) pair.\n\nThe mapper's output is a complete accounting: for every variant and every\nannotation level in that variant's deterministically-reachable set, there is one\nrecord carrying its outcome -- never a silent omission. This field is uniform across\nmeasured (assay-level) and projected (deterministic non-assay) records so the two can\nbe treated identically by consumers; it distinguishes a benign absence from a genuine\nfailure, which a populated ``error_message`` alone cannot.\n\n- ``MAPPED`` -- a VRS allele was produced (``pre_mapped``/``post_mapped`` populated).\n- ``INTRONIC`` -- the variant's coding projection is intronic: no VRS-representable\n coding form and no protein consequence. Benign (``error_message`` is ``None``).\n- ``NO_PROTEIN_CONSEQUENCE`` -- the protein layer was reachable but yields no\n projectable protein change (e.g. UTR). Benign (``error_message`` is ``None``).\n- ``FAILED`` -- the mapping/projection genuinely failed (mis-selected transcript,\n projection error, unresolvable reference contig). ``error_message`` carries detail.",
"enum": [
"mapped",
"intronic",
"no_protein_consequence",
"failed"
],
"title": "MappingOutcome",
"type": "string"
},
"Number": {
"description": "Define VRS 1.3 Number.",
"properties": {
Expand Down Expand Up @@ -295,6 +306,17 @@
"default": null,
"title": "Near Gap"
},
"outcome": {
"anyOf": [
{
"$ref": "#/$defs/MappingOutcome"
},
{
"type": "null"
}
],
"default": null
},
"post_mapped": {
"anyOf": [
{
Expand Down
32 changes: 26 additions & 6 deletions src/dcd_mapping/align.py
Original file line number Diff line number Diff line change
Expand Up @@ -919,7 +919,22 @@ def _get_best_match(
tcoords = coords[0]
qcoords = coords[1]

strand = Strand.POSITIVE if int(qcoords[0]) <= int(qcoords[-1]) else Strand.NEGATIVE
protein_vs_dna = "-q=prot" in blat_params.get("target_args", "")

# For cDNA queries the strand is read from qcoords direction: cDNA on the
# negative strand is reverse-complemented, so its qcoords decrease.
# For protein queries qcoords always increase (protein reads N→C regardless
# of genome strand), so qcoords direction is uninformative — use tcoords
# instead (they decrease when the gene is on the minus strand).
if protein_vs_dna:
strand = (
Strand.POSITIVE if int(tcoords[0]) <= int(tcoords[-1]) else Strand.NEGATIVE
)
else:
strand = (
Strand.POSITIVE if int(qcoords[0]) <= int(qcoords[-1]) else Strand.NEGATIVE
)

q_start = int(qcoords.min())
q_end = int(qcoords.max())

Expand All @@ -943,10 +958,9 @@ def _get_best_match(
if ts == te or qs == qe:
continue

hit_subranges.append(SequenceRange(start=ts, end=te))
hit_subranges.append(SequenceRange(start=min(ts, te), end=max(ts, te)))
query_subranges.append(SequenceRange(start=min(qs, qe), end=max(qs, qe)))

protein_vs_dna = "-q=prot" in blat_params.get("target_args", "")
alignment_qc = _build_alignment_qc(best_aln, protein_vs_dna=protein_vs_dna)

return AlignmentResult(
Expand All @@ -958,7 +972,10 @@ def _get_best_match(
coverage=coverage,
query_range=SequenceRange(start=q_start, end=q_end),
query_subranges=query_subranges,
hit_range=SequenceRange(start=int(tcoords[0]), end=int(tcoords[-1])),
hit_range=SequenceRange(
start=min(int(tcoords[0]), int(tcoords[-1])),
end=max(int(tcoords[0]), int(tcoords[-1])),
),
hit_subranges=hit_subranges,
score=float(_scores[id(best_aln)]),
next_best_score=next_best,
Expand Down Expand Up @@ -1251,7 +1268,7 @@ def align_target_to_protein(
qs, qe = int(qcoords[i]), int(qcoords[i + 1])
if ts == te or qs == qe:
continue
hit_subranges.append(SequenceRange(start=ts, end=te))
hit_subranges.append(SequenceRange(start=min(ts, te), end=max(ts, te)))
query_subranges.append(SequenceRange(start=min(qs, qe), end=max(qs, qe)))

# Attach full sequences so _build_alignment_qc can do per-base mismatch
Expand All @@ -1265,7 +1282,10 @@ def align_target_to_protein(
result = AlignmentResult(
query_range=SequenceRange(start=int(qcoords.min()), end=int(qcoords.max())),
query_subranges=query_subranges,
hit_range=SequenceRange(start=int(tcoords[0]), end=int(tcoords[-1])),
hit_range=SequenceRange(
start=min(int(tcoords[0]), int(tcoords[-1])),
end=max(int(tcoords[0]), int(tcoords[-1])),
),
hit_subranges=hit_subranges,
percent_identity=_blat_style_identity(
best_counts.identities,
Expand Down
Loading
Loading