+ *
+ * @author Kiran Garimella, Mark DePristo
+ */
+public class ErrorRatePerCycle extends LocusWalker {
+ @Output PrintStream out;
+ @Argument(fullName="min_base_quality_score", shortName="mbq", doc="Minimum base quality required to consider a base for calling", required=false)
+ public Integer MIN_BASE_QUAL = 0;
+ @Argument(fullName="min_mapping_quality_score", shortName="mmq", doc="Minimum read mapping quality required to consider a read for calling", required=false)
+ public Integer MIN_MAPPING_QUAL = 20;
+
+ private GATKReport report;
+ private GATKReportTable table;
+ private final static String reportName = "ErrorRatePerCycle";
+ private final static String reportDescription = "The error rate per sequenced position in the reads";
+
+ /**
+ * Allows us to use multiple records for the key (read group x cycle)
+ */
+ private static class TableKey implements Comparable {
+ final String readGroup;
+ final int cycle;
+
+ private TableKey(final String readGroup, final int cycle) {
+ this.readGroup = readGroup;
+ this.cycle = cycle;
+ }
+
+ @Override
+ public int compareTo(final TableKey tableKey) {
+ final int scmp = readGroup.compareTo(tableKey.readGroup);
+ if ( scmp == 0 )
+ return Integer.valueOf(cycle).compareTo(tableKey.cycle);
+ else
+ return scmp;
+ }
+ }
+
+ public void initialize() {
+ report = new GATKReport();
+ report.addTable(reportName, reportDescription);
+ table = report.getTable(reportName);
+ table.addPrimaryKey("key", false);
+ table.addColumn("readgroup", 0);
+ table.addColumn("cycle", 0);
+ table.addColumn("mismatches", 0);
+ table.addColumn("counts", 0);
+ table.addColumn("qual", 0);
+ table.addColumn("errorrate", 0.0f, "%.2e");
+ }
+
+ public Integer map(RefMetaDataTracker tracker, ReferenceContext ref, AlignmentContext context) {
+ for ( final PileupElement p : context.getBasePileup() ) {
+ final GATKSAMRecord read = p.getRead();
+ final int offset = p.getOffset();
+ final boolean firstOfPair = ! read.getReadPairedFlag() || read.getFirstOfPairFlag();
+
+ if ( firstOfPair && read.getMappingQuality() >= MIN_MAPPING_QUAL && p.getQual() >= MIN_BASE_QUAL ) {
+ final byte readBase = p.getBase();
+ final byte refBase = ref.getBase();
+ final int cycle = offset;
+
+ if ( BaseUtils.isRegularBase(readBase) && BaseUtils.isRegularBase(refBase) ) {
+ final TableKey key = new TableKey(read.getReadGroup().getReadGroupId(), cycle);
+
+ if ( ! table.containsKey(key) ) {
+ table.set(key, "cycle", cycle);
+ table.set(key, "readgroup", read.getReadGroup().getReadGroupId());
+ }
+
+ table.increment(key, "counts");
+ if (readBase != refBase)
+ table.increment(key, "mismatches");
+ }
+ }
+ }
+
+ return null;
+ }
+
+ public Integer reduceInit() { return null; }
+
+ public Integer reduce(Integer value, Integer sum) { return null; }
+
+ public void onTraversalDone(Integer sum) {
+ for ( final Object key : table.getPrimaryKeys() ) {
+ final int mismatches = (Integer)table.get(key, "mismatches");
+ final int count = (Integer)table.get(key, "counts");
+ final double errorRate = (mismatches + 1) / (1.0*(count + 1));
+ final int qual = QualityUtils.probToQual(1-errorRate, 0.0);
+ table.set(key, "qual", qual);
+ table.set(key, "errorrate", errorRate);
+ }
+
+ report.print(out);
+ }
+}
diff --git a/public/java/src/org/broadinstitute/sting/gatk/walkers/diagnostics/ReadGroupProperties.java b/public/java/src/org/broadinstitute/sting/gatk/walkers/diagnostics/ReadGroupProperties.java
index d7a48d321..14985907d 100644
--- a/public/java/src/org/broadinstitute/sting/gatk/walkers/diagnostics/ReadGroupProperties.java
+++ b/public/java/src/org/broadinstitute/sting/gatk/walkers/diagnostics/ReadGroupProperties.java
@@ -94,9 +94,6 @@ import java.util.Map;
*
* @author Mark DePristo
*/
-
-
-
public class ReadGroupProperties extends ReadWalker {
@Output
public PrintStream out;
diff --git a/public/java/test/org/broadinstitute/sting/gatk/walkers/diagnostics/ErrorRatePerCycleIntegrationTest.java b/public/java/test/org/broadinstitute/sting/gatk/walkers/diagnostics/ErrorRatePerCycleIntegrationTest.java
new file mode 100644
index 000000000..accb9c0cf
--- /dev/null
+++ b/public/java/test/org/broadinstitute/sting/gatk/walkers/diagnostics/ErrorRatePerCycleIntegrationTest.java
@@ -0,0 +1,41 @@
+/*
+ * Copyright (c) 2012, The Broad Institute
+ *
+ * Permission is hereby granted, free of charge, to any person
+ * obtaining a copy of this software and associated documentation
+ * files (the "Software"), to deal in the Software without
+ * restriction, including without limitation the rights to use,
+ * copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following
+ * conditions:
+ *
+ * The above copyright notice and this permission notice shall be
+ * included in all copies or substantial portions of the Software.
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+ * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
+ * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+ * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
+ * HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
+ * WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
+ * OTHER DEALINGS IN THE SOFTWARE.
+ */
+
+package org.broadinstitute.sting.gatk.walkers.diagnostics;
+
+import org.broadinstitute.sting.WalkerTest;
+import org.testng.annotations.Test;
+
+import java.util.Arrays;
+
+public class ErrorRatePerCycleIntegrationTest extends WalkerTest {
+ @Test
+ public void basicTest() {
+ WalkerTestSpec spec = new WalkerTestSpec(
+ "-T ErrorRatePerCycle -R " + b37KGReference + " -I " + b37GoodBAM + " -L 20:10,000,000-10,100,000 -o %s",
+ 1,
+ Arrays.asList("0cc212ecb6df300e321784039ff29f13"));
+ executeTest("ErrorRatePerCycle:", spec);
+ }
+}
\ No newline at end of file