From 4a8a6468beb0ce17b378193bcf818e7d30969ace Mon Sep 17 00:00:00 2001
From: chartl <chartl@348d0f76-0448-11de-a6fe-93d51630548a>
Date: Tue, 20 Oct 2009 19:39:32 +0000
Subject: [PATCH] Use read group as a condition for confusion tables. With an
 integration test.

Changed BaseTransitionTable to comparable objects for consistent ordering of output
( e.g. so the integration test doesn't yell so much )




git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@1889 348d0f76-0448-11de-a6fe-93d51630548a
---
 ...seTransitionTableCalculatorJavaWalker.java | 41 +++++++++++++++++--
 ...ionTableCalculatorJavaIntegrationTest.java | 10 ++++-
 2 files changed, 47 insertions(+), 4 deletions(-)

diff --git a/java/src/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaWalker.java b/java/src/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaWalker.java
index cd11176ae..410338a61 100644
--- a/java/src/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaWalker.java
+++ b/java/src/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaWalker.java
@@ -42,6 +42,8 @@ public class BaseTransitionTableCalculatorJavaWalker extends LocusWalker<Referen
     boolean usePileupMismatches = false;
 @Argument(fullName="usePreviousReadBases", doc="Use previous bases of the read as part of the calculation. Will ignore reads if there aren't this many previous bases. Uses the specified number. Defaults to 0", required=false)
     int nPreviousReadBases = 0;
+@Argument(fullName="useReadGroup", doc="Use the group number of the read as a condition of the table.", required = false)
+    boolean useReadGroup = false;
 
     private UnifiedGenotyper ug;
     private ReferenceContextWindow refWindow;
@@ -51,7 +53,7 @@ public class BaseTransitionTableCalculatorJavaWalker extends LocusWalker<Referen
         ug = new UnifiedGenotyper();
         ug.initialize();
         refWindow = new ReferenceContextWindow(nPreviousBases);
-        conditionalTables = new HashSet<BaseTransitionTable>();
+        conditionalTables = new TreeSet<BaseTransitionTable>();
     }
 
     public Integer reduceInit() {
@@ -146,7 +148,9 @@ public class BaseTransitionTableCalculatorJavaWalker extends LocusWalker<Referen
             return false;
         } else if ( read.getBaseQualities()[offset] <= minQualityScore ) {
             return false;
-        } else {
+        } else if ( useSecondaryBase && read.getAttribute("SQ") == null )
+            return false;
+        else {
             return true;
         }
     }
@@ -173,6 +177,10 @@ public class BaseTransitionTableCalculatorJavaWalker extends LocusWalker<Referen
             conditions.add(countMismatches(map.getPileup()));
         }
 
+        if ( useReadGroup ) {
+            conditions.add(read.getReadGroup().getReadGroupId());
+        }
+
         return conditions;
     }
 
@@ -195,6 +203,10 @@ public class BaseTransitionTableCalculatorJavaWalker extends LocusWalker<Referen
             header = header + "\tNumber_of_pileup_mismatches";
         }
 
+        if ( useReadGroup ) {
+            header = header + "\tRead_group";
+        }
+
         return String.format("%s\t%s%n",header,"Counts");
     }
 
@@ -250,7 +262,7 @@ public class BaseTransitionTableCalculatorJavaWalker extends LocusWalker<Referen
 }
 
 
-class BaseTransitionTable {
+class BaseTransitionTable implements Comparable {
 
     private int[][] table;
     private List<Comparable> conditions;
@@ -286,6 +298,29 @@ class BaseTransitionTable {
         }
     }
 
+    public int compareTo(Object obj) {
+        if ( ! ( obj instanceof BaseTransitionTable ) ) {
+            return -1;
+        } else {
+            BaseTransitionTable t = (BaseTransitionTable) obj;
+            if ( this.conditionsMatch(t.conditions) ) {
+                return 0;
+            } else {
+                if ( this.numConditions() == t.numConditions() ) {
+                    ListIterator<Comparable> thisIter = this.conditions.listIterator();
+                    ListIterator<Comparable> thatIter = t.conditions.listIterator();
+                    while ( thisIter.next() == thatIter.next() ) {
+                        // do nothing
+                    }
+                    return thisIter.previous().compareTo(thatIter.previous());
+                } else {
+                    return (this.numConditions() > t.numConditions() ) ? 1 : -1;
+                }
+            }
+        }
+
+    }
+
     public void print( PrintStream out ) {
 
         for ( char observedBase : BaseUtils.BASES ) {
diff --git a/java/test/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaIntegrationTest.java b/java/test/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaIntegrationTest.java
index 25a52554f..ea345d20a 100644
--- a/java/test/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaIntegrationTest.java
+++ b/java/test/org/broadinstitute/sting/playground/gatk/walkers/BaseTransitionTableCalculatorJavaIntegrationTest.java
@@ -21,6 +21,7 @@ public class BaseTransitionTableCalculatorJavaIntegrationTest extends WalkerTest
     public static final String OUTPUT_MD5_LOWCONFIDENTREFTHRESHOLD = "521a9fa7716ed22550c2ba3fe3409070";
     public static final String OUTPUT_MD5_HIGHCONFIDENTREFTHRESHOLD = "8ab6d389fc494881736e9a58126c2f1b";
     public static final String OUTPUT_MD5_ALLARGUMENTS = "f45481946d7a5c70078d432b0baff083";
+    public static final String OUTPUT_MD5_USEREADGROUP = "6dbf18e96f3367c8764b4d1e15eb956b";
     public static final String LOCUS = "1:10,000,000-10,200,000";
     public static final String BAM_FILE = "/humgen/gsa-scr1/GATK_Data/Validation_Data/NA12878.1kg.p2.chr1_10mb_11_mb.SLX.bam";
     public static final String REFERENCE = "/broad/1KG/reference/human_b36_both.fasta";
@@ -64,7 +65,7 @@ public class BaseTransitionTableCalculatorJavaIntegrationTest extends WalkerTest
     public void testBaseTransitionCalculatorJavaHighConfidentRefThreshold() {
         String args = "-T BaseTransitionTableCalculatorJava -o %s -I "+BAM_FILE+" -L "+LOCUS+" -R "+REFERENCE+ " --confidentRefThreshold 8";
         WalkerTest.WalkerTestSpec spec =  new WalkerTest.WalkerTestSpec(args,1,Arrays.asList(OUTPUT_MD5_HIGHCONFIDENTREFTHRESHOLD));
-        executeTest("BaseTransitionTableCalculatorJava: Low Ref Threshold",spec);
+        executeTest("BaseTransitionTableCalculatorJava: High Ref Threshold",spec);
     }
 
     @Test
@@ -73,4 +74,11 @@ public class BaseTransitionTableCalculatorJavaIntegrationTest extends WalkerTest
         WalkerTest.WalkerTestSpec spec =  new WalkerTest.WalkerTestSpec(args,1,Arrays.asList(OUTPUT_MD5_ALLARGUMENTS));
         executeTest("BaseTransitionTableCalculatorJava: Low Ref Threshold, Low Q Score, Low Mapping Quality, Additional Mismatches",spec);
     }
+
+    @Test
+    public void testBaseTransitionTableCalculatorJavaUseReadGroup() {
+        String args = "-T BaseTransitionTableCalculatorJava -o %s -I "+BAM_FILE+" -L "+LOCUS+" -R "+REFERENCE+" --useReadGroup";
+        WalkerTest.WalkerTestSpec spec =  new WalkerTest.WalkerTestSpec(args,1,Arrays.asList(OUTPUT_MD5_USEREADGROUP));
+        executeTest("testBaseTransitionTableCalculatorJava: Use read group", spec);
+    }
 }