Fix for Adam K's reported bug: we weren't handling reads that were entirely insertions properly in LIBS. Specifically, the event bases were off-by-one (which was disasterous in Adam's case with a 1bp read). Added a unit test to cover this case.

This commit is contained in:
Eric Banks 2012-08-20 23:12:41 -04:00
parent 5b1781fdac
commit 40d5efc804
2 changed files with 43 additions and 11 deletions

View File

@ -301,6 +301,7 @@ public class LocusIteratorByState extends LocusIterator {
final CigarOperator op = state.getCurrentCigarOperator(); // current cigar operator final CigarOperator op = state.getCurrentCigarOperator(); // current cigar operator
final CigarElement nextElement = state.peekForwardOnGenome(); // next cigar element final CigarElement nextElement = state.peekForwardOnGenome(); // next cigar element
final CigarElement lastElement = state.peekBackwardOnGenome(); // last cigar element final CigarElement lastElement = state.peekBackwardOnGenome(); // last cigar element
final boolean isSingleElementCigar = nextElement == lastElement;
final CigarOperator nextOp = nextElement.getOperator(); // next cigar operator final CigarOperator nextOp = nextElement.getOperator(); // next cigar operator
final CigarOperator lastOp = lastElement.getOperator(); // last cigar operator final CigarOperator lastOp = lastElement.getOperator(); // last cigar operator
final int readOffset = state.getReadOffset(); // the base offset on this read final int readOffset = state.getReadOffset(); // the base offset on this read
@ -308,7 +309,7 @@ public class LocusIteratorByState extends LocusIterator {
final boolean isBeforeDeletion = nextOp == CigarOperator.DELETION; final boolean isBeforeDeletion = nextOp == CigarOperator.DELETION;
final boolean isAfterDeletion = lastOp == CigarOperator.DELETION; final boolean isAfterDeletion = lastOp == CigarOperator.DELETION;
final boolean isBeforeInsertion = nextOp == CigarOperator.INSERTION; final boolean isBeforeInsertion = nextOp == CigarOperator.INSERTION;
final boolean isAfterInsertion = lastOp == CigarOperator.INSERTION; final boolean isAfterInsertion = lastOp == CigarOperator.INSERTION && !isSingleElementCigar;
final boolean isNextToSoftClip = nextOp == CigarOperator.S || (state.getGenomeOffset() == 0 && read.getSoftStart() != read.getAlignmentStart()); final boolean isNextToSoftClip = nextOp == CigarOperator.S || (state.getGenomeOffset() == 0 && read.getSoftStart() != read.getAlignmentStart());
int nextElementLength = nextElement.getLength(); int nextElementLength = nextElement.getLength();
@ -328,8 +329,10 @@ public class LocusIteratorByState extends LocusIterator {
else { else {
if (!filterBaseInRead(read, location.getStart())) { if (!filterBaseInRead(read, location.getStart())) {
String insertedBaseString = null; String insertedBaseString = null;
if (nextOp == CigarOperator.I) if (nextOp == CigarOperator.I) {
insertedBaseString = new String(Arrays.copyOfRange(read.getReadBases(), readOffset + 1, readOffset + 1 + nextElement.getLength())); final int insertionOffset = isSingleElementCigar ? 0 : 1;
insertedBaseString = new String(Arrays.copyOfRange(read.getReadBases(), readOffset + insertionOffset, readOffset + insertionOffset + nextElement.getLength()));
}
pile.add(new PileupElement(read, readOffset, false, isBeforeDeletion, isAfterDeletion, isBeforeInsertion, isAfterInsertion, isNextToSoftClip, insertedBaseString, nextElementLength)); pile.add(new PileupElement(read, readOffset, false, isBeforeDeletion, isAfterDeletion, isBeforeInsertion, isAfterInsertion, isNextToSoftClip, insertedBaseString, nextElementLength));
size++; size++;

View File

@ -104,7 +104,7 @@ public class LocusIteratorByStateUnitTest extends BaseTest {
after.setBaseQualities(new byte[] {20,20,20,20,20,20,20,20,20,20}); after.setBaseQualities(new byte[] {20,20,20,20,20,20,20,20,20,20});
after.setCigarString("10M"); after.setCigarString("10M");
List<SAMRecord> reads = Arrays.asList(before,during,after); List<SAMRecord> reads = Arrays.asList(before, during, after);
// create the iterator by state with the fake reads and fake records // create the iterator by state with the fake reads and fake records
li = makeLTBS(reads,readAttributes); li = makeLTBS(reads,readAttributes);
@ -134,9 +134,9 @@ public class LocusIteratorByStateUnitTest extends BaseTest {
// create a test version of the Reads object // create a test version of the Reads object
ReadProperties readAttributes = createTestReadProperties(); ReadProperties readAttributes = createTestReadProperties();
SAMRecord indelOnlyRead = ArtificialSAMUtils.createArtificialRead(header,"indelOnly",0,firstLocus,76); SAMRecord indelOnlyRead = ArtificialSAMUtils.createArtificialRead(header, "indelOnly", 0, firstLocus, 76);
indelOnlyRead.setReadBases(Utils.dupBytes((byte)'A',76)); indelOnlyRead.setReadBases(Utils.dupBytes((byte)'A',76));
indelOnlyRead.setBaseQualities(Utils.dupBytes((byte)'@',76)); indelOnlyRead.setBaseQualities(Utils.dupBytes((byte) '@', 76));
indelOnlyRead.setCigarString("76I"); indelOnlyRead.setCigarString("76I");
List<SAMRecord> reads = Arrays.asList(indelOnlyRead); List<SAMRecord> reads = Arrays.asList(indelOnlyRead);
@ -148,10 +148,10 @@ public class LocusIteratorByStateUnitTest extends BaseTest {
// and considers it to be an indel-containing read. // and considers it to be an indel-containing read.
Assert.assertTrue(li.hasNext(),"Should have found a whole-indel read in the normal base pileup without extended events enabled"); Assert.assertTrue(li.hasNext(),"Should have found a whole-indel read in the normal base pileup without extended events enabled");
AlignmentContext alignmentContext = li.next(); AlignmentContext alignmentContext = li.next();
Assert.assertEquals(alignmentContext.getLocation().getStart(),firstLocus,"Base pileup is at incorrect location."); Assert.assertEquals(alignmentContext.getLocation().getStart(), firstLocus, "Base pileup is at incorrect location.");
ReadBackedPileup basePileup = alignmentContext.getBasePileup(); ReadBackedPileup basePileup = alignmentContext.getBasePileup();
Assert.assertEquals(basePileup.getReads().size(),1,"Pileup is of incorrect size"); Assert.assertEquals(basePileup.getReads().size(),1,"Pileup is of incorrect size");
Assert.assertSame(basePileup.getReads().get(0),indelOnlyRead,"Read in pileup is incorrect"); Assert.assertSame(basePileup.getReads().get(0), indelOnlyRead, "Read in pileup is incorrect");
} }
/** /**
@ -168,7 +168,7 @@ public class LocusIteratorByStateUnitTest extends BaseTest {
leadingRead.setCigarString("1M75I"); leadingRead.setCigarString("1M75I");
SAMRecord indelOnlyRead = ArtificialSAMUtils.createArtificialRead(header,"indelOnly",0,secondLocus,76); SAMRecord indelOnlyRead = ArtificialSAMUtils.createArtificialRead(header,"indelOnly",0,secondLocus,76);
indelOnlyRead.setReadBases(Utils.dupBytes((byte)'A',76)); indelOnlyRead.setReadBases(Utils.dupBytes((byte) 'A', 76));
indelOnlyRead.setBaseQualities(Utils.dupBytes((byte)'@',76)); indelOnlyRead.setBaseQualities(Utils.dupBytes((byte)'@',76));
indelOnlyRead.setCigarString("76I"); indelOnlyRead.setCigarString("76I");
@ -177,7 +177,7 @@ public class LocusIteratorByStateUnitTest extends BaseTest {
fullMatchAfterIndel.setBaseQualities(Utils.dupBytes((byte)'@',76)); fullMatchAfterIndel.setBaseQualities(Utils.dupBytes((byte)'@',76));
fullMatchAfterIndel.setCigarString("75I1M"); fullMatchAfterIndel.setCigarString("75I1M");
List<SAMRecord> reads = Arrays.asList(leadingRead,indelOnlyRead,fullMatchAfterIndel); List<SAMRecord> reads = Arrays.asList(leadingRead, indelOnlyRead, fullMatchAfterIndel);
// create the iterator by state with the fake reads and fake records // create the iterator by state with the fake reads and fake records
li = makeLTBS(reads, createTestReadProperties()); li = makeLTBS(reads, createTestReadProperties());
@ -204,7 +204,36 @@ public class LocusIteratorByStateUnitTest extends BaseTest {
numAlignmentContextsFound++; numAlignmentContextsFound++;
} }
Assert.assertEquals(numAlignmentContextsFound,2,"Found incorrect number of alignment contexts"); Assert.assertEquals(numAlignmentContextsFound, 2, "Found incorrect number of alignment contexts");
}
/**
* Test to make sure that reads supporting only an indel (example cigar string: 76I) do
* not negatively influence the ordering of the pileup.
*/
@Test
public void testWholeIndelReadTest2() {
final int firstLocus = 44367788, secondLocus = firstLocus + 1;
SAMRecord read = ArtificialSAMUtils.createArtificialRead(header,"read",0,secondLocus,1);
read.setReadBases(Utils.dupBytes((byte) 'A', 1));
read.setBaseQualities(Utils.dupBytes((byte) '@', 1));
read.setCigarString("1I");
List<SAMRecord> reads = Arrays.asList(read);
// create the iterator by state with the fake reads and fake records
li = makeLTBS(reads, createTestReadProperties());
while(li.hasNext()) {
AlignmentContext alignmentContext = li.next();
ReadBackedPileup p = alignmentContext.getBasePileup();
Assert.assertTrue(p.getNumberOfElements() == 1);
PileupElement pe = p.iterator().next();
Assert.assertTrue(pe.isBeforeInsertion());
Assert.assertFalse(pe.isAfterInsertion());
Assert.assertEquals(pe.getEventBases(), "A");
}
} }
private static ReadProperties createTestReadProperties() { private static ReadProperties createTestReadProperties() {