diff --git a/java/src/org/broadinstitute/sting/playground/gatk/walkers/BeagleOutputToVCFWalker.java b/java/src/org/broadinstitute/sting/playground/gatk/walkers/BeagleOutputToVCFWalker.java index 9056d917c..61991deba 100755 --- a/java/src/org/broadinstitute/sting/playground/gatk/walkers/BeagleOutputToVCFWalker.java +++ b/java/src/org/broadinstitute/sting/playground/gatk/walkers/BeagleOutputToVCFWalker.java @@ -204,7 +204,7 @@ public class BeagleOutputToVCFWalker extends RodWalker { ArrayList beagleGenotypePairs = beaglePhasedFeature.getGenotypes().get(sample); Allele originalAlleleA = g.getAllele(0); - Allele originalAlleleB = g.getAllele(1); + Allele originalAlleleB = (g.getAlleles().size() == 2) ? g.getAllele(1) : g.getAllele(0); // hack to deal with no-call genotypes // We have phased genotype in hp. Need to set the isRef field in the allele. List alleles = new ArrayList(); @@ -213,6 +213,9 @@ public class BeagleOutputToVCFWalker extends RodWalker { String alleleB = beagleGenotypePairs.get(1); byte[] r = alleleA.getBytes(); + + //System.out.println(context.getLocation() + " : " + alleleA + " " + alleleB); + byte rA = r[0]; Boolean isRefA = (refByte == rA); diff --git a/java/src/org/broadinstitute/sting/playground/gatk/walkers/ProduceBeagleInputWalker.java b/java/src/org/broadinstitute/sting/playground/gatk/walkers/ProduceBeagleInputWalker.java index 8277469b4..8a7876cfa 100755 --- a/java/src/org/broadinstitute/sting/playground/gatk/walkers/ProduceBeagleInputWalker.java +++ b/java/src/org/broadinstitute/sting/playground/gatk/walkers/ProduceBeagleInputWalker.java @@ -148,6 +148,26 @@ public class ProduceBeagleInputWalker extends RodWalker { beagleGenotypesWriter.format("%c %c ", a, b); } } + else if (genotype.isCalled() && !genotype.hasAttribute(VCFConstants.GENOTYPE_LIKELIHOODS_KEY)) { + // hack to deal with input VCFs with no genotype likelihoods. Just assume the called genotype + // is confident. This is useful for Hapmap and 1KG release VCFs. + double AA = 0.02; + double AB = 0.02; + double BB = 0.02; + + if (genotype.isHomRef()) { AA = 0.96; } + else if (genotype.isHet()) { AB = 0.96; } + else if (genotype.isHomVar()) { BB = 0.96; } + + beagleWriter.printf("%.2f %.2f %.2f ", AA, AB, BB); + + if (beagleGenotypesWriter != null) { + char a = genotype.getAllele(0).toString().charAt(0); + char b = genotype.getAllele(0).toString().charAt(0); + + beagleGenotypesWriter.format("%c %c ", a, b); + } + } else { beagleWriter.print("0.33 0.33 0.33 "); // write 1/3 likelihoods for uncalled genotypes. if (beagleGenotypesWriter != null)