Added --fdr_filter_level to ApplyVariantCuts so that you can create beautiful tranche plots and also decide which tranche level to filter at. The previous version always filtered at the smallest tranche. The tranche filter names are appropriately added to the VCF header.

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@3747 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
rpoplin 2010-07-08 22:44:10 +00:00
parent 760aaeda88
commit 5d39cd5db8
2 changed files with 18 additions and 16 deletions

View File

@ -65,6 +65,9 @@ public class ApplyVariantCuts extends RodWalker<Integer, Integer> {
private String TRANCHE_FILENAME = "optimizer.dat.tranches";
@Argument(fullName="outputVCFFile", shortName="outputVCF", doc="The output filtered VCF file", required=true)
private String OUTPUT_FILENAME = "optimizer.vcf";
@Argument(fullName="fdr_filter_level", shortName="fdr_filter_level", doc="The FDR level at which to start filtering.", required=true)
private double FDR_FILTER_LEVEL = 0.0;
/////////////////////////////
// Private Member Variables
@ -87,8 +90,10 @@ public class ApplyVariantCuts extends RodWalker<Integer, Integer> {
for( final String line : new XReadLines(new File( TRANCHE_FILENAME )) ) {
if( !firstLine ) {
final String[] vals = line.split(",");
qCuts.add(Double.parseDouble(vals[2]));
filterName.add(vals[4]);
if(Double.parseDouble(vals[0]) >= FDR_FILTER_LEVEL) {
qCuts.add(Double.parseDouble(vals[2]));
filterName.add(vals[4]);
}
}
firstLine = false;
}
@ -119,11 +124,10 @@ public class ApplyVariantCuts extends RodWalker<Integer, Integer> {
}
}
// BUGBUG : why doesn't this work? VCFFilterHeaderLine is protected?
//for( int iii = 1; iii < filterName.size(); iii++ ) {
// hInfo.add(new VCFFilterHeaderLine(filterName.get(iii)), String.format("FDR tranche level at qual " + qCuts.get(iii)));
//}
//hInfo.add(new VCFFilterHeaderLine(filterName.get(0)+"+"), String.format("FDR tranche level at qual > " + qCuts.get(filterName.size()-1)));
for( int iii = 1; iii < filterName.size(); iii++ ) {
hInfo.add(new VCFFilterHeaderLine(filterName.get(iii), String.format("FDR tranche level at qual " + qCuts.get(iii))));
}
hInfo.add(new VCFFilterHeaderLine(filterName.get(0)+"+", String.format("FDR tranche level at qual > " + qCuts.get(filterName.size()-1))));
final VCFHeader vcfHeader = new VCFHeader(hInfo, samples);
vcfWriter.writeHeader(vcfHeader);

View File

@ -155,15 +155,13 @@ public class VariantRecalibrator extends RodWalker<ExpandingArrayList<VariantDat
// Set up default values for the FDR tranches if necessary
if( FDR_TRANCHES == null ) {
FDR_TRANCHES = new Double[3];
FDR_TRANCHES[0] = 12.5;
FDR_TRANCHES[1] = 20.0;
FDR_TRANCHES[2] = 30.0;
//FDR_TRANCHES[0] = 1.0;
//FDR_TRANCHES[1] = 5.0;
//FDR_TRANCHES[2] = 10.0;
//FDR_TRANCHES[3] = 20.0;
//FDR_TRANCHES[4] = 30.0;
FDR_TRANCHES = new Double[6];
FDR_TRANCHES[0] = 0.1;
FDR_TRANCHES[1] = 1.0;
FDR_TRANCHES[2] = 5.0;
FDR_TRANCHES[3] = 12.5;
FDR_TRANCHES[4] = 25.0;
FDR_TRANCHES[5] = 35.0;
}
}