zzh
/
gatk-3.8
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

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
This commit is contained in:
rpoplin 2010-11-30 17:19:54 +00:00
parent 21256909bb
commit af84462f3e
2 changed files with 11 additions and 5 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
java/src/org/broadinstitute/sting/gatk/walkers/beagle/BeagleOutputToVCFWalker.java
View File

@ -86,7 +86,11 @@ public class BeagleOutputToVCFWalker extends RodWalker<Integer, Integer> {
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<ReferenceOrderedDataSource> dataSources = this.getToolkit().getRodDataSources();
@ -307,7 +311,7 @@ public class BeagleOutputToVCFWalker extends RodWalker<Integer, Integer> {
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<String> removedFilters = vc_input.filtersWereApplied() ? new HashSet<String>(vc_input.getFilters()) : new HashSet<String>(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<Allele>(Arrays.asList(vc_input.getReference())), genotypes, vc_input.getNegLog10PError(), removedFilters, vc_input.getAttributes());
}
@ -324,7 +328,9 @@ public class BeagleOutputToVCFWalker extends RodWalker<Integer, Integer> {
}
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());

4
java/test/org/broadinstitute/sting/gatk/walkers/beagle/BeagleIntegrationTest.java
View File

@ -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);
}