Merged bug fix from Stable into Unstable

public/java/src/org/broadinstitute/sting/gatk/walkers/annotator/InbreedingCoeff.java

@@ -4,6 +4,7 @@ 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.annotator.interfaces.ActiveRegionBasedAnnotation;
+import org.broadinstitute.sting.gatk.walkers.Walker;
 import org.broadinstitute.sting.gatk.walkers.annotator.interfaces.AnnotatorCompatibleWalker;
 import org.broadinstitute.sting.gatk.walkers.annotator.interfaces.InfoFieldAnnotation;
 import org.broadinstitute.sting.gatk.walkers.annotator.interfaces.StandardAnnotation;
@@ -16,10 +17,7 @@ import org.broadinstitute.sting.utils.variantcontext.Genotype;
 import org.broadinstitute.sting.utils.variantcontext.GenotypesContext;
 import org.broadinstitute.sting.utils.variantcontext.VariantContext;
 
-import java.util.Arrays;
+import java.util.*;
-import java.util.HashMap;
-import java.util.List;
-import java.util.Map;
 
 
 /**
@@ -33,8 +31,12 @@ import java.util.Map;
 public class InbreedingCoeff extends InfoFieldAnnotation implements StandardAnnotation, ActiveRegionBasedAnnotation {
 
     private static final int MIN_SAMPLES = 10;
+    private Set<String> founderIds;
 
     public Map<String, Object> annotate(RefMetaDataTracker tracker, AnnotatorCompatibleWalker walker, ReferenceContext ref, Map<String, AlignmentContext> stratifiedContexts, VariantContext vc) {
+        //If available, get the founder IDs and cache them. the IC will only be computed on founders then.
+        if(founderIds == null)
+            founderIds =  ((Walker)walker).getSampleDB().getFounderIds();
         return calculateIC(vc);
     }
 
@@ -43,7 +45,7 @@ public class InbreedingCoeff extends InfoFieldAnnotation implements StandardAnno
     }
 
     private Map<String, Object> calculateIC(final VariantContext vc) {
-        final GenotypesContext genotypes = vc.getGenotypes();
+        final GenotypesContext genotypes = (founderIds == null || founderIds.isEmpty()) ? vc.getGenotypes() : vc.getGenotypes(founderIds);
         if ( genotypes == null || genotypes.size() < MIN_SAMPLES )
             return null;
 

public/java/src/org/broadinstitute/sting/gatk/walkers/phasing/PhaseByTransmission.java

@@ -82,6 +82,9 @@ public class PhaseByTransmission extends RodWalker<HashMap<Byte,Integer>, HashMa
     @Argument(shortName = "prior",required = false,fullName = "DeNovoPrior", doc="Prior for de novo mutations. Default: 1e-8")
     private double deNovoPrior=1e-8;
 
+    @Argument(shortName = "fatherAlleleFirst",required = false,fullName = "FatherAlleleFirst", doc="Ouputs the father allele as the first allele in phased child genotype. i.e. father|mother rather than mother|father.")
+    private boolean fatherFAlleleFirst=false;
+
     @Output
     protected VCFWriter vcfWriter = null;
 
@@ -183,12 +186,15 @@ public class PhaseByTransmission extends RodWalker<HashMap<Byte,Integer>, HashMa
             ArrayList<Allele> parentPhasedAlleles = new ArrayList<Allele>(2);
             ArrayList<Allele> childPhasedAlleles = new ArrayList<Allele>(2);
 
-            //If there is a possible phasing between the mother and child => phase
+            //If there is a possible phasing between the parent and child => phase
             int childTransmittedAlleleIndex = childAlleles.indexOf(parentAlleles.get(0));
             if(childTransmittedAlleleIndex > -1){
                 trioPhasedGenotypes.put(parent, new Genotype(DUMMY_NAME, parentAlleles, Genotype.NO_LOG10_PERROR, null, null, true));
                 childPhasedAlleles.add(childAlleles.remove(childTransmittedAlleleIndex));
-                childPhasedAlleles.add(childAlleles.get(0));
+                if(parent.equals(FamilyMember.MOTHER))
+                        childPhasedAlleles.add(childAlleles.get(0));
+                else
+                        childPhasedAlleles.add(0,childAlleles.get(0));
                 trioPhasedGenotypes.put(FamilyMember.CHILD, new Genotype(DUMMY_NAME, childPhasedAlleles, Genotype.NO_LOG10_PERROR, null, null, true));
             }
             else if((childTransmittedAlleleIndex = childAlleles.indexOf(parentAlleles.get(1))) > -1){
@@ -196,7 +202,10 @@ public class PhaseByTransmission extends RodWalker<HashMap<Byte,Integer>, HashMa
                 parentPhasedAlleles.add(parentAlleles.get(0));
                 trioPhasedGenotypes.put(parent, new Genotype(DUMMY_NAME, parentPhasedAlleles, Genotype.NO_LOG10_PERROR, null, null, true));
                 childPhasedAlleles.add(childAlleles.remove(childTransmittedAlleleIndex));
-                childPhasedAlleles.add(childAlleles.get(0));
+                if(parent.equals(FamilyMember.MOTHER))
+                    childPhasedAlleles.add(childAlleles.get(0));
+                else
+                    childPhasedAlleles.add(0,childAlleles.get(0));
                 trioPhasedGenotypes.put(FamilyMember.CHILD, new Genotype(DUMMY_NAME, childPhasedAlleles, Genotype.NO_LOG10_PERROR, null, null, true));
             }
             //This is a Mendelian Violation => Do not phase
@@ -296,6 +305,14 @@ public class PhaseByTransmission extends RodWalker<HashMap<Byte,Integer>, HashMa
             else{
                 phaseFamilyAlleles(mother, father, child);
             }
+
+            //If child should phased genotype should be father first, then swap the alleles
+            if(fatherFAlleleFirst && trioPhasedGenotypes.get(FamilyMember.CHILD).isPhased()){
+                ArrayList<Allele> childAlleles = new ArrayList<Allele>(trioPhasedGenotypes.get(FamilyMember.CHILD).getAlleles());
+                childAlleles.add(childAlleles.remove(0));
+                trioPhasedGenotypes.put(FamilyMember.CHILD,new Genotype(DUMMY_NAME,childAlleles,Genotype.NO_LOG10_PERROR,null,null,true));
+            }
+
         }
 
         /**

public/java/test/org/broadinstitute/sting/gatk/walkers/annotator/VariantAnnotatorIntegrationTest.java

@@ -198,4 +198,14 @@ public class VariantAnnotatorIntegrationTest extends WalkerTest {
         executeTest("Testing ChromosomeCounts annotation with PED file", spec);
     }
 
+    @Test
+    public void testInbreedingCoeffPed() {
+        final String MD5 = "7f1314fada5cb1f35ba1996f8a7a686b";
+        WalkerTestSpec spec = new WalkerTestSpec(
+                "-T VariantAnnotator -R " + b37KGReference + " -A InbreedingCoeff --variant:vcf " + validationDataLocation + "ug.random50000.subset300bp.chr1.family.vcf" +
+                        " -L " + validationDataLocation + "ug.random50000.subset300bp.chr1.family.vcf -NO_HEADER -ped " + validationDataLocation + "ug.random50000.family.ped -o %s", 1,
+                Arrays.asList(MD5));
+        executeTest("Testing InbreedingCoeff annotation with PED file", spec);
+    }
+
 }

public/java/test/org/broadinstitute/sting/gatk/walkers/phasing/PhaseByTransmissionIntegrationTest.java

@@ -29,7 +29,7 @@ public class PhaseByTransmissionIntegrationTest extends WalkerTest {
                         "-o %s"
                 ),
                 2,
-                Arrays.asList("16fefda693156eadf1481fd9de23facb","9418a7a6405b78179ca13a67b8bfcc14")
+                Arrays.asList("d54a142d68dca54e478c13f9a0e4c95c","1a37fcc93a73429f9065b942ab771233")
         );
         executeTest("testTrueNegativeMV", spec);
     }
@@ -48,7 +48,7 @@ public class PhaseByTransmissionIntegrationTest extends WalkerTest {
                         "-o %s"
                 ),
                 2,
-                Arrays.asList("14cf1d21a54d8b9fb506df178b634c56","efc66ae3d036715b721f9bd35b65d556")
+                Arrays.asList("883ea7fd2b200c4b7fa95a4f7aa15931","7b1f5309c3d4f4aa7e9061f288dceb68")
         );
         executeTest("testTruePositiveMV", spec);
     }
@@ -67,7 +67,7 @@ public class PhaseByTransmissionIntegrationTest extends WalkerTest {
                         "-o %s"
                 ),
                 2,
-                Arrays.asList("f9b0fae9fe1e0f09b883a292b0e70a12","398724bc1e65314cc5ee92706e05a3ee")
+                Arrays.asList("e812d62a3449b74b6948ee7deb8a0790","d00922496759e84c66a4b5e222e36997")
         );
         executeTest("testFalsePositiveMV", spec);
     }
@@ -86,7 +86,7 @@ public class PhaseByTransmissionIntegrationTest extends WalkerTest {
                         "-o %s"
                 ),
                 2,
-                Arrays.asList("b8d1aa3789ce77b45430c62d13ee3006","a1a333e08fafb288cda0e7711909e1c3")
+                Arrays.asList("e3c572f933a40e1878a2cfa52049517a","60e4f0be344fb944ab3378f9ab27da64")
         );
         executeTest("testSpecialCases", spec);
     }
@@ -108,7 +108,7 @@ public class PhaseByTransmissionIntegrationTest extends WalkerTest {
                         "-o %s"
                 ),
                 2,
-                Arrays.asList("7201ce7cc47db5840ac6b647709f7c33","c11b5e7cd7459d90d0160f917eff3b1e")
+                Arrays.asList("b42af3b73a2cb38cfc92f8047dd686b3","a69c3f9c005e852b44c29ab25e87ba0d")
         );
         executeTest("testPriorOption", spec);
     }
@@ -128,9 +128,30 @@ public class PhaseByTransmissionIntegrationTest extends WalkerTest {
                         "-o %s"
                 ),
                 1,
-                Arrays.asList("398724bc1e65314cc5ee92706e05a3ee")
+                Arrays.asList("d00922496759e84c66a4b5e222e36997")
         );
         executeTest("testMVFileOption", spec);
     }
 
+    //Test when running with the fatherAlleleFirst option
+    @Test
+    public void testFatherAlleleFirst() {
+        WalkerTestSpec spec = new WalkerTestSpec(
+                buildCommandLine(
+                        "-T PhaseByTransmission",
+                        "-NO_HEADER",
+                        "-R " + b37KGReference,
+                        "--variant " + TPTest,
+                        "-ped "+ goodFamilyFile,
+                        "-L 1:10109-10315",
+                        "-mvf %s",
+                        "-o %s",
+                        "-fatherAlleleFirst"
+                ),
+                2,
+                Arrays.asList("c158a3816357597543ef85c4478c41e8","4f8daca19c8f31bd87850c124f91e330")
+        );
+        executeTest("testFatherAlleleFirst", spec);
+    }
+
 }