Skip to content

Commit 033d92e

Browse files
committed
CDS length calculation bug fix
1 parent 663c62f commit 033d92e

2 files changed

Lines changed: 75 additions & 75 deletions

File tree

biojava-genome/src/main/java/org/biojava/nbio/genome/util/ChromosomeMappingTools.java

Lines changed: 17 additions & 75 deletions
Original file line numberDiff line numberDiff line change
@@ -504,7 +504,6 @@ public static ChromPos getChromPosForward(int cdsPos, List<Integer> exonStarts,
504504
*/
505505
public static int getCDSLengthReverse(List<Integer> exonStarts, List<Integer> exonEnds, int cdsStart, int cdsEnd) {
506506

507-
boolean inCoding = false;
508507
int codingLength = 0;
509508

510509
if (cdsEnd < cdsStart) {
@@ -513,8 +512,6 @@ public static int getCDSLengthReverse(List<Integer> exonStarts, List<Integer> ex
513512
cdsStart = tmp;
514513
}
515514

516-
int lengthExons = 0;
517-
518515
// map reverse
519516
for (int i = exonStarts.size() - 1; i >= 0; i--) {
520517

@@ -526,53 +523,19 @@ public static int getCDSLengthReverse(List<Integer> exonStarts, List<Integer> ex
526523
end = start;
527524
start = tmp;
528525
}
529-
lengthExons += end - start;
530-
531-
if (start <= cdsEnd && end >= cdsEnd) {
532-
inCoding = true;
533-
534-
int tmpstart = start;
535-
if (start < cdsStart) {
536-
tmpstart = cdsStart;
537-
}
538-
codingLength += (cdsEnd - tmpstart);
539-
540-
boolean debug = logger.isDebugEnabled();
541526

542-
if ( debug) {
527+
if ((start < cdsStart && end < cdsStart) || (start > cdsEnd && end > cdsEnd))
528+
continue;
543529

544-
StringBuffer b = new StringBuffer();
530+
if (start < cdsStart)
531+
start = cdsStart;
545532

546-
b.append(" UTR :" + (cdsEnd + 1) + " - " + (end) + newline);
547-
if (tmpstart == start)
548-
b.append(" -> ");
549-
else
550-
b.append(" <-> ");
551-
b.append("Exon :" + tmpstart + " - " + cdsEnd + " | " + (cdsEnd - tmpstart) + " | " + codingLength + " | " + (codingLength % 3) + newline);
552-
// single exon with UTR on both ends
553-
if (tmpstart != start)
554-
b.append(" UTR :" + (cdsStart - 1) + " - " + start + newline);
555-
logger.debug(b.toString());
533+
if (end > cdsEnd)
534+
end = cdsEnd;
556535

557-
}
558-
} else if (start <= cdsStart && end >= cdsStart) {
559-
inCoding = false;
560-
codingLength += (end - cdsStart);
561-
562-
logger.debug(" <- Exon : " + (cdsStart+1) + " - " + end + " | " + (end - cdsStart) + " | " + (codingLength-3) + " | " + (codingLength % 3));
563-
logger.debug(" UTR : " + start + " - " + (cdsStart ));
564-
565-
} else if (inCoding) {
566-
// full exon is coding
567-
codingLength += (end - start);
568-
logger.debug(" Exon : " + start + " - " + end + " | " + (end - start) + " | " + codingLength + " | " + (codingLength % 3));
569-
} else {
570-
// e.g. see UBQLN3
571-
logger.debug(" no translation!");
572-
}
536+
codingLength += (end - start + 1);
573537
}
574-
logger.debug("length exons: " + lengthExons + " codin length: " + (codingLength - 3));
575-
return codingLength - 3;
538+
return codingLength;
576539
}
577540

578541
/**
@@ -585,47 +548,26 @@ public static int getCDSLengthReverse(List<Integer> exonStarts, List<Integer> ex
585548
* @return
586549
*/
587550
public static int getCDSLengthForward(List<Integer> exonStarts, List<Integer> exonEnds, int cdsStart, int cdsEnd) {
588-
boolean inCoding = false;
551+
589552
int codingLength = 0;
590553

591-
int lengthExons = 0;
592-
// map forward
593554
for (int i = 0; i < exonStarts.size(); i++) {
594555

595556
int start = exonStarts.get(i);
596557
int end = exonEnds.get(i);
597-
lengthExons += end - start;
598558

599-
logger.debug("forward exon: " + (start+1) + " - " + end + " | " + (end - start));
559+
if ( (start < cdsStart && end < cdsStart) || (start > cdsEnd && end > cdsEnd) )
560+
continue;
600561

601-
if (start+1 <= cdsStart +1 && end >= cdsStart+1) {
602-
603-
inCoding = true;
604-
codingLength += (end - cdsStart);
562+
if (start < cdsStart)
563+
start = cdsStart;
605564

606-
logger.debug(" UTR : " + start + " - " + (cdsStart ));
607-
logger.debug(" -> Exon : " + (cdsStart+1) + " - " + end + " | " + (end - cdsStart+1) + " | " + codingLength + " | " + (codingLength % 3));
565+
if (end > cdsEnd)
566+
end = cdsEnd;
608567

609-
} else if (start+1 <= cdsEnd && end >= cdsEnd) {
610-
611-
inCoding = false;
612-
codingLength += (cdsEnd - start);
613-
614-
logger.debug(" <- Exon : " + (start +1)+ " - " + cdsEnd + " | " + (cdsEnd - start+1) + " | " + codingLength + " | " + (codingLength % 3));
615-
logger.debug(" UTR : " + cdsEnd + 1 + " - " + end);
616-
617-
} else if (inCoding) {
618-
// full exon is coding
619-
codingLength += (end - start);
620-
621-
logger.debug(" Exon :" + (start+1) + " - " + end + " | " + (end - start+1) + " | " + codingLength + " | " + (codingLength % 3));
622-
}
568+
codingLength += (end - start + 1);
623569
}
624-
625-
logger.debug("length exons: " + Integer.toString(lengthExons));
626-
logger.debug("CDS length:" + Integer.toString((codingLength-3)));
627-
628-
return codingLength-3 ;
570+
return codingLength;
629571
}
630572

631573
/** Extracts the exon boundaries in CDS coordinates. (needs to be divided by 3 to get AA positions)
Lines changed: 58 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,58 @@
1+
package org.biojava.nbio.genome;
2+
3+
import org.biojava.nbio.genome.util.ChromosomeMappingTools;
4+
import org.junit.Test;
5+
6+
import java.util.ArrayList;
7+
import java.util.Arrays;
8+
import java.util.List;
9+
10+
import static junitx.framework.ComparableAssert.assertEquals;
11+
12+
/**
13+
* Created by Yana Valasatava on 8/14/17.
14+
*/
15+
public class TestChromosomeMappingTools {
16+
17+
@Test
18+
public void testGetCDSLengthForward() {
19+
20+
List<Integer> exonStarts = new ArrayList<>(Arrays.asList(10, 30, 50, 70));
21+
List<Integer> exonEnds = new ArrayList<>(Arrays.asList(20, 40, 60, 80));
22+
int cdsStart = 35;
23+
int cdsEnd = 75;
24+
25+
int cdsDesired = 23;
26+
int cdsTest = ChromosomeMappingTools.getCDSLengthForward(exonStarts, exonEnds, cdsStart, cdsEnd);
27+
28+
assertEquals(cdsDesired, cdsTest);
29+
}
30+
31+
@Test
32+
public void testGetCDSLengthReverseAsc() {
33+
34+
List<Integer> exonStarts = new ArrayList<>(Arrays.asList(10, 50, 70));
35+
List<Integer> exonEnds = new ArrayList<>(Arrays.asList(20, 60, 80));
36+
int cdsStart = 55;
37+
int cdsEnd = 75;
38+
39+
int cdsDesired = 12;
40+
int cdsTest = ChromosomeMappingTools.getCDSLengthReverse(exonStarts, exonEnds, cdsStart, cdsEnd);
41+
42+
assertEquals(cdsDesired, cdsTest);
43+
}
44+
45+
@Test
46+
public void testGetCDSLengthReverseDesc() {
47+
48+
List<Integer> exonStarts = new ArrayList<>(Arrays.asList(70, 50, 10));
49+
List<Integer> exonEnds = new ArrayList<>(Arrays.asList(80, 60, 20));
50+
int cdsStart = 75;
51+
int cdsEnd = 55;
52+
53+
int cdsDesired = 12;
54+
int cdsTest = ChromosomeMappingTools.getCDSLengthReverse(exonStarts, exonEnds, cdsStart, cdsEnd);
55+
56+
assertEquals(cdsDesired, cdsTest);
57+
}
58+
}

0 commit comments

Comments
 (0)