fix: predictCoding uses CDS-mapped width for exon/intron boundary deletions (#83)

A deletion starting in an exon and extending into the intron produced
incorrect REFCODON/VARCODON because .getRefCodons() and the frameshift
calculation used the genomic width of the variant rather than the
transcript-space (CDSLOC) width.

For example, a 51bp genomic deletion with only 37bp overlapping the CDS
would compute cend using 51, causing the reference codon to extend
incorrectly into the next exon's sequence.

Fix: replace width(txlocal) (genomic width) with
width(mcols(txlocal)$CDSLOC) (CDS-mapped width) in both:
- .getRefCodons(): codon boundary calculation
- frameshift detection: refwidth for length-change check

Fixes #83

Author: jmg421

R/methods-predictCoding.R

@@ -113,8 +113,9 @@ setMethod("predictCoding", c("VRanges", "TxDb", "ANY", "missing"),
         mcols(txlocal)$varAllele <- va
     }
 
-    ## frameshift
-    refwidth <- width(txlocal)
+    ## frameshift: use CDS-mapped width, not genomic width, because a
+    ## deletion extending into an intron has genomic width > CDS overlap (#83).
+    refwidth <- width(mcols(txlocal)$CDSLOC)
     altallele <- mcols(txlocal)$varAllele
     fmshift <- abs(width(altallele) - refwidth) %% 3 != 0 
     if (any(fmshift))
@@ -200,10 +201,13 @@ setMethod("predictCoding", c("VRanges", "TxDb", "ANY", "missing"),
 .getRefCodons <- function(txlocal, altpos, seqSource, cdsbytx)
 { 
     ## adjust codon end for 
-    ## - width of the reference sequence
+    ## - width of the reference sequence in transcript space
     ## - position of alt allele substitution in the codon
+    ## Use CDSLOC width (transcript-space) not genomic width, because
+    ## a deletion extending into an intron has genomic width >> CDS width (#83).
+    cds_width <- width(mcols(txlocal)$CDSLOC)
     cstart <- ((start(mcols(txlocal)$CDSLOC) - 1L) %/% 3L) * 3L + 1L
-    cend <- cstart + (((altpos + width(txlocal) - 2L) %/% 3L) * 3L + 2L)
+    cend <- cstart + (((altpos + cds_width - 2L) %/% 3L) * 3L + 2L)
     txord <- match(mcols(txlocal)$TXID, names(cdsbytx))
     txseqs <- extractTranscriptSeqs(seqSource, cdsbytx[txord])
     DNAStringSet(substring(txseqs, cstart, cend))