From af84462f3ea40a8e55a032b57dd034681fabbdfb Mon Sep 17 00:00:00 2001 From: rpoplin Date: Tue, 30 Nov 2010 17:19:54 +0000 Subject: [PATCH] The dev team has decided to change the filter that is added to records that are set to monomorphic by Beagle. It no longer lists the reference allele. Added those filters to the header of the output VCF file. Finally, we no longer use R2=NaN values coming from Beagle. git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@4757 348d0f76-0448-11de-a6fe-93d51630548a --- .../gatk/walkers/beagle/BeagleOutputToVCFWalker.java | 12 +++++++++--- .../gatk/walkers/beagle/BeagleIntegrationTest.java | 4 ++-- 2 files changed, 11 insertions(+), 5 deletions(-) diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/beagle/BeagleOutputToVCFWalker.java b/java/src/org/broadinstitute/sting/gatk/walkers/beagle/BeagleOutputToVCFWalker.java index dfd7077fe..f713e341e 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/beagle/BeagleOutputToVCFWalker.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/beagle/BeagleOutputToVCFWalker.java @@ -86,7 +86,11 @@ public class BeagleOutputToVCFWalker extends RodWalker { hInfo.addAll(VCFUtils.getHeaderFields(getToolkit())); hInfo.add(new VCFFormatHeaderLine("OG",1, VCFHeaderLineType.String, "Original Genotype input to Beagle")); hInfo.add(new VCFInfoHeaderLine("R2", 1, VCFHeaderLineType.Float, "r2 Value reported by Beagle on each site")); - hInfo.add(new VCFInfoHeaderLine("NumGenotypesChanged", 1, VCFHeaderLineType.Integer, "r2 Value reported by Beagle on each site")); + hInfo.add(new VCFInfoHeaderLine("NumGenotypesChanged", 1, VCFHeaderLineType.Integer, "The number of genotypes changed by Beagle")); + hInfo.add(new VCFFilterHeaderLine("BGL_RM_WAS_A", "This 'A' site was set to monomorphic by Beagle")); + hInfo.add(new VCFFilterHeaderLine("BGL_RM_WAS_C", "This 'C' site was set to monomorphic by Beagle")); + hInfo.add(new VCFFilterHeaderLine("BGL_RM_WAS_G", "This 'G' site was set to monomorphic by Beagle")); + hInfo.add(new VCFFilterHeaderLine("BGL_RM_WAS_T", "This 'T' site was set to monomorphic by Beagle")); // Open output file specified by output VCF ROD final List dataSources = this.getToolkit().getRodDataSources(); @@ -307,7 +311,7 @@ public class BeagleOutputToVCFWalker extends RodWalker { filteredVC = new VariantContext("outputvcf", vc_input.getChr(), vc_input.getStart(), vc_input.getEnd(), vc_input.getAlleles(), genotypes, vc_input.getNegLog10PError(), vc_input.filtersWereApplied() ? vc_input.getFilters() : null, vc_input.getAttributes()); else { Set removedFilters = vc_input.filtersWereApplied() ? new HashSet(vc_input.getFilters()) : new HashSet(1); - removedFilters.add(String.format("BGL_RM_WAS_%s/%s",vc_input.getReference().getBaseString(),vc_input.getAlternateAllele(0))); + removedFilters.add(String.format("BGL_RM_WAS_%s",vc_input.getAlternateAllele(0))); filteredVC = new VariantContext("outputvcf", vc_input.getChr(), vc_input.getStart(), vc_input.getEnd(), new HashSet(Arrays.asList(vc_input.getReference())), genotypes, vc_input.getNegLog10PError(), removedFilters, vc_input.getAttributes()); } @@ -324,7 +328,9 @@ public class BeagleOutputToVCFWalker extends RodWalker { } attributes.put("NumGenotypesChanged", numGenotypesChangedByBeagle ); - attributes.put("R2", beagleR2Feature.getR2value().toString() ); + if( !beagleR2Feature.getR2value().equals(Double.NaN) ) { + attributes.put("R2", beagleR2Feature.getR2value().toString() ); + } vcfWriter.add(VariantContext.modifyAttributes(filteredVC,attributes), ref.getBase()); diff --git a/java/test/org/broadinstitute/sting/gatk/walkers/beagle/BeagleIntegrationTest.java b/java/test/org/broadinstitute/sting/gatk/walkers/beagle/BeagleIntegrationTest.java index a4e2b6c91..e7c2c9deb 100755 --- a/java/test/org/broadinstitute/sting/gatk/walkers/beagle/BeagleIntegrationTest.java +++ b/java/test/org/broadinstitute/sting/gatk/walkers/beagle/BeagleIntegrationTest.java @@ -41,7 +41,7 @@ public class BeagleIntegrationTest extends WalkerTest { "-B:beagleR2,BEAGLE " + beagleValidationDataLocation + "inttestbgl.r2 " + "-B:beagleProbs,BEAGLE " + beagleValidationDataLocation + "inttestbgl.gprobs " + "-B:beaglePhased,BEAGLE " + beagleValidationDataLocation + "inttestbgl.phased " + - "-o %s -NO_HEADER", 1, Arrays.asList("93546f4f6a7f5fe24f01a357dcad7c5f")); + "-o %s -NO_HEADER", 1, Arrays.asList("e746763ade40edea56a051dc1dfd6165")); executeTest("test BeagleOutputToVCF", spec); } @@ -72,7 +72,7 @@ public class BeagleIntegrationTest extends WalkerTest { "-B:beagleR2,beagle /humgen/gsa-hpprojects/GATK/data/Validation_Data/EUR_beagle_in_test.r2 "+ "-B:beagleProbs,beagle /humgen/gsa-hpprojects/GATK/data/Validation_Data/EUR_beagle_in_test.gprobs.bgl "+ "-B:beaglePhased,beagle /humgen/gsa-hpprojects/GATK/data/Validation_Data/EUR_beagle_in_test.phased.bgl "+ - "-L 20:1-70000 -o %s -NO_HEADER ",1,Arrays.asList("f8071eda9f762947c6a3404d6d126df2")); + "-L 20:1-70000 -o %s -NO_HEADER ",1,Arrays.asList("71e7ac1118ecef82a7de0b6543973a81")); executeTest("testBeagleChangesSitesToRef",spec); }