PrintCoverageWalker functionality moved to DepthOfCoverageWalker. Added integration tests.

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@2247 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
ebanks 2009-12-03 17:23:59 +00:00
parent 608fa4cc3a
commit c2017cc91b
3 changed files with 36 additions and 89 deletions

View File

@ -39,14 +39,9 @@ import java.util.*;
*/
@By(DataSource.REFERENCE)
public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoCInfo, DepthOfCoverageWalker.DoCInfo> {
enum printType {
NONE,
COMPACT,
DETAILED
}
@Argument(fullName="printStyle", shortName = "s", doc="Printing style: NONE, COMPACT, or DETAILED", required=false)
protected printType printStyle = printType.COMPACT;
@Argument(fullName="printBaseCounts", shortName = "bases", doc="Print individual base counts (A,C,G,T only)", required=false)
protected boolean printBaseCounts = false;
@Argument(fullName="minMAPQ", shortName = "minMAPQ", doc="If provided, we will also list read counts with MAPQ >= this value at a locus in coverage",required=false)
protected int excludeMAPQBelowThis = -1;
@ -107,6 +102,9 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
header.append("\tcoverage_atleast_MQ");
header.append(excludeMAPQBelowThis);
}
if ( printBaseCounts ) {
header.append("\tA_count\tC_count\tG_count\tT_count");
}
if ( byReadGroup ) {
for ( String rg : readGroupNames ) {
header.append("\tcoverage_for_");
@ -139,6 +137,12 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
else if ( p.isDeletion() )
nDeletionReads++;
if ( printBaseCounts ) {
int baseIndex = BaseUtils.simpleBaseToBaseIndex(p.getBase());
if ( baseIndex != -1 )
info.baseCounts[baseIndex]++;
}
SAMReadGroupRecord readGroup = p.getRead().getReadGroup();
if ( readGroup == null )
continue;
@ -187,6 +191,10 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
if ( value.totalCoverage >= minDepthForPercentage ) {
sum.minDepthCoveredLoci++;
}
if ( printBaseCounts ) {
for (int baseIndex = 0; baseIndex < BaseUtils.BASES.length; baseIndex++ )
sum.baseCounts[baseIndex] += value.baseCounts[baseIndex];
}
if ( byReadGroup ) {
for ( String rg : readGroupNames ) {
int oldDepth = sum.depthByReadGroup.get(rg);
@ -220,6 +228,9 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
header.append("\tpercent_loci_covered_atleast_depth");
header.append(minDepthForPercentage);
}
if ( printBaseCounts ) {
header.append("\tA_count\tC_count\tG_count\tT_count");
}
if ( byReadGroup ) {
for ( String rg : readGroupNames ) {
header.append("\tcoverage_for_");
@ -321,10 +332,12 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
sb.append("\t");
}
sb.append((info.totalCoverage - info.numDeletions));
if ( printAverageCoverage ) {
sb.append("\t");
sb.append(String.format("%.2f", ((double)(info.totalCoverage - info.numDeletions)) / totalBases));
}
if ( excludeMAPQBelowThis > 0 ) {
sb.append("\t");
sb.append((info.totalCoverage - info.numBadMQReads));
@ -333,10 +346,19 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
sb.append(String.format("%.2f", ((double)(info.totalCoverage - info.numBadMQReads)) / totalBases));
}
}
if ( minDepthForPercentage >= 0 ) {
sb.append("\t");
sb.append(String.format("%.2f", ((double)info.minDepthCoveredLoci) / totalBases));
}
if ( printBaseCounts ) {
for (int baseIndex = 0; baseIndex < BaseUtils.BASES.length; baseIndex++ ) {
sb.append("\t");
sb.append(String.format("%8d", info.baseCounts[baseIndex]));
}
}
if ( byReadGroup ) {
for ( String rg : readGroupNames ) {
sb.append("\t");
@ -360,10 +382,15 @@ public class DepthOfCoverageWalker extends LocusWalker<DepthOfCoverageWalker.DoC
public int numBadMQReads = 0;
public int minDepthCoveredLoci = 0;
public int[] baseCounts = null;
public HashMap<String, Integer> depthByReadGroup = null;
public HashMap<String, Integer> depthBySample = null;
public DoCInfo() {
if ( printBaseCounts ) {
baseCounts = new int[4];
}
if ( byReadGroup ) {
depthByReadGroup = new HashMap<String, Integer>();
for ( String readGroupName : readGroupNames )

View File

@ -1,81 +0,0 @@
package org.broadinstitute.sting.playground.gatk.walkers;
import net.sf.samtools.SAMRecord;
import org.broadinstitute.sting.gatk.contexts.AlignmentContext;
import org.broadinstitute.sting.gatk.contexts.ReferenceContext;
import org.broadinstitute.sting.gatk.refdata.RefMetaDataTracker;
import org.broadinstitute.sting.gatk.walkers.By;
import org.broadinstitute.sting.gatk.walkers.DataSource;
import org.broadinstitute.sting.gatk.walkers.LocusWalker;
import java.util.List;
@By(DataSource.REFERENCE)
public class PrintCoverageWalker extends LocusWalker<Integer, Integer> {
public void initialize() {
out.println("chrom position reference a c g t");
}
public String walkerType() { return "ByLocus"; }
// Do we actually want to operate on the context?
public boolean filter(RefMetaDataTracker tracker, ReferenceContext ref, AlignmentContext context) {
return true; // We are keeping all the reads
}
public Integer map(RefMetaDataTracker tracker, ReferenceContext ref, AlignmentContext context) {
int aCount = 0;
int cCount = 0;
int gCount = 0;
int tCount = 0;
char upRef = Character.toUpperCase(ref.getBase());
List<SAMRecord> reads = context.getReads();
for ( int i = 0; i < reads.size(); i++ )
{
SAMRecord read = reads.get(i);
if (read.getNotPrimaryAlignmentFlag() ||
read.getDuplicateReadFlag() ||
read.getReadUnmappedFlag() ||
read.getMappingQuality() <= 0
) {
continue;
}
int offset = context.getOffsets().get(i);
char base = (char)read.getReadBases()[offset];
if (base == 'a' || base == 'A') { aCount++; }
if (base == 'c' || base == 'C') { cCount++; }
if (base == 'g' || base == 'G') { gCount++; }
if (base == 't' || base == 'T') { tCount++; }
}
out.printf("%s %d %s %d %d %d %d\n",
context.getContig(),
context.getPosition(),
upRef,
aCount,
cCount,
gCount,
tCount);
return -1;
}
// Given result of map function
public Integer reduceInit() {
return 0;
}
public Integer reduce(Integer value, Integer sum) {
return value + sum;
}
@Override
public void onTraversalDone(Integer result) {
}
}

View File

@ -20,7 +20,8 @@ public class DepthOfCoverageIntegrationTest extends WalkerTest {
expectations.put("-bySample", "93358437153b4d65bdff747e33de1d63");
expectations.put("-byRG", "777e8427eb4bdad300b23800cb7b0592");
expectations.put("-histogram", "96f15e1d9d598d48191e20ee84715d46");
expectations.put("-minMAPQ 1 -bySample -byRG -minDepth 8 -histogram", "783b0bc83c54d883efa8383a379ff17b");
expectations.put("-bases", "baafcb2b90098cad1c5950da9e9932a6");
expectations.put("-minMAPQ 1 -bySample -byRG -minDepth 8 -histogram -bases", "eb27297504e003ce7bb42300dfd436fe");
}
@Test