Repository navigation
vdjtools 4.7.0 — the V/J boundary from the germline (#182)
model.germline_boundary reports where the V and J germlines stop, in CDR3 nucleotide space.
infer_nt_batch carries it as v_end_germline / j_start_germline beside v_end / j_start,
which are unchanged.
Scenario.v_end is a property of the argmax recombination history and the right answer to that
question. It is the wrong answer to "where does the V germline stop": maximising P(sequence)
explains N-region nucleotides as templated whenever it can, so the boundary walks outward. That is
the defect #182 reports.
The replacement needs no DP and no recombination model — only the germline and the genetic code.
An exonuclease does not stop on a codon boundary, so the last residue a germline touches is usually
part germline and part N region, and an amino-acid alignment can only round that to a whole residue.
What the amino acid still fixes is which nucleotides are admissible: 17 of the 20 residues open
every one of their codons with the same base, so against a germline GGA (Gly) an observed Glu can
only be GAA or GAG — both begin with the germline's G, so the germline demonstrably reaches one
nucleotide further. Only Ser, Leu and Arg are ambiguous in the first base, and there the codons are
weighed by how far each lets the germline reach. The decision lives in arda.cdr3fix.boundary_nt
(arda 2.31.0, a hard floor); this passes it the model's own cdr3_segment, never arda's anchor
table and never cut_segment.
Measured against isalgo/airr_control's human.trb.ntvj, on the 8,132 VDJdb human TRB junctions
whose boundary every control observation agrees on and whose V/J the bundled olga:human_T_beta
model carries. The control's VEnd/JStart are a maximal nucleotide alignment — verified, not
assumed: recomputing it from the observed cdr3nt reproduces them on 93.5 % of V and 98.8 % of J
observations, the residual being allele polymorphism.
| exact | Scenario |
germline_boundary |
amino-acid alignment alone |
|---|---|---|---|
v.end, VDJdb residues (nt + 1) // 3 |
7,244 (89.08 %) | 7,554 (92.89 %) | 5,839 (71.79 %) |
j.start, VDJdb residues ceil(nt / 3) |
7,483 (92.02 %) | 7,966 (97.96 %) | 7,968 (97.97 %) |
| V boundary, nucleotides | 5,293 (65.09 %) | 6,536 (80.37 %) | 2,355 (28.96 %) |
| J boundary, nucleotides | 5,301 (65.19 %) | 6,056 (74.47 %) | 2,462 (30.27 %) |
The last column is VDJdb's legacy k-mer scanner and arda.cdr3fix 2.30.1 — the same number, because
both answer in whole residues. So the codon decision buys 21 points over rounding and 4 over the
argmax history, and this is the first thing here to beat the k-mer scanner rather than tie it.
j.start cannot move in VDJdb's residue convention and does not: it is the first fully
J-templated residue, the same residue whether the germline reaches one or two nucleotides into the
one before it. The nucleotide answer moves a great deal, and is the one to read.
The two sides are bounded by different things. The V residue count is right on 83.99 % of records
and the codon extension on 95.71 % of those, so V is limited by the protein alignment. The J
count is right on 97.97 % and the extension on 76.02 %, because the nucleotide a J boundary turns
on is the codon's third — the one position the genetic code leaves free.
germline_boundary is computed for every row of infer_nt_batch, including one the DP declined:
the germline explains what it explains whether or not a recombination was found. It declines only
when there is no germline to align to. Cost is 2.67 us/row against infer_nt_batch's 82.81.