Java Code Examples for htsjdk.variant.vcf.VCFConstants

The following examples show how to use htsjdk.variant.vcf.VCFConstants. These examples are extracted from open source projects. You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may check out the related API usage on the sidebar.
Example 1
@Test
public void testLowQualVariantToGQ0HomRef() {
    final ReblockGVCF reblocker = new ReblockGVCF();

    reblocker.dropLowQuals = true;
    final Genotype g = makeG("sample1", LONG_REF, Allele.NON_REF_ALLELE, 200, 100, 200, 11, 0, 37);
    final VariantContext toBeNoCalled = makeDeletionVC("lowQualVar", Arrays.asList(LONG_REF, DELETION, Allele.NON_REF_ALLELE), LONG_REF.length(), g);
    final VariantContext dropped = reblocker.lowQualVariantToGQ0HomRef(toBeNoCalled, toBeNoCalled);
    Assert.assertEquals(dropped, null);

    reblocker.dropLowQuals = false;
    final VariantContext modified = reblocker.lowQualVariantToGQ0HomRef(toBeNoCalled, toBeNoCalled);
    Assert.assertTrue(modified.getAttributes().containsKey(VCFConstants.END_KEY));
    Assert.assertTrue(modified.getAttributes().get(VCFConstants.END_KEY).equals(13));
    Assert.assertTrue(modified.getReference().equals(SHORT_REF));
    Assert.assertTrue(modified.getAlternateAllele(0).equals(Allele.NON_REF_ALLELE));
    Assert.assertTrue(!modified.filtersWereApplied());
    Assert.assertTrue(modified.getLog10PError() == VariantContext.NO_LOG10_PERROR);

    //No-calls were throwing NPEs.  Now they're not.
    final Genotype g2 = makeG("sample1", Allele.NO_CALL,Allele.NO_CALL);
    final VariantContext noData = makeDeletionVC("noData", Arrays.asList(LONG_REF, DELETION, Allele.NON_REF_ALLELE), LONG_REF.length(), g2);
    final VariantContext notCrashing = reblocker.lowQualVariantToGQ0HomRef(noData, noData);
    Assert.assertTrue(notCrashing.getGenotype(0).isNoCall());
}
 
Example 2
Source Project: genomewarp   Source File: VcfToVariant.java    License: Apache License 2.0 6 votes vote down vote up
private static boolean parseInlineGenotypeFields(String field, VariantCall.Builder vcBuilder,
    ListValue.Builder lvBuilder, IntGenotypeFieldAccessors.Accessor accessor, Genotype g) {

  final int[] intValues = accessor.getValues(g);
  if (intValues == null || intValues.length == 0) {
    return false;
  }

  if (field.equals(VCFConstants.GENOTYPE_PL_KEY)) {
    // HTSJDK folds GL's into PL's. We only use PL's to store genotype likelihood.
    for (int i = 0; i < intValues.length; i++) {
      // We add 0.0 to remove the possiblity of getting -0.0.
      vcBuilder.addGenotypeLikelihood(-(double) intValues[i] / 10.0 + 0.0);
    }
    return false;
  }

  for (int i = 0; i < intValues.length; i++) {
    lvBuilder.addValues(Value.newBuilder().setNumberValue(intValues[i]));
  }
  return true;
}
 
Example 3
private void calculateAFs(final Iterable<Genotype> genotypes) {
    final int[] truthAC = new int[CopyNumberTriStateAllele.ALL_ALLELES.size()];
    final int[] callsAC = new int[CopyNumberTriStateAllele.ALL_ALLELES.size()];
    for (final Genotype genotype : genotypes) {
        final List<Allele> alleles = genotype.getAlleles();
        if (alleles.size() > 1) {
            throw new GATKException("unexpected CNV genotype ploidy: " + alleles.size());
        } else if (!alleles.isEmpty() && alleles.get(0).isCalled()) {
            final int index = CopyNumberTriStateAllele.valueOf(alleles.get(0)).index();
            callsAC[index]++;
        }
        final String truthGT = String.valueOf(genotype.getExtendedAttribute(VariantEvaluationContext.TRUTH_GENOTYPE_KEY, VCFConstants.MISSING_VALUE_v4));
        final int truthAlleleIndex = truthGT.equals(VCFConstants.MISSING_VALUE_v4) ? -1 : Integer.parseInt(truthGT);
        if (truthAlleleIndex >= 0) {
            final List<Allele> contextAlleles = getAlleles();
            if (truthAlleleIndex >= contextAlleles.size()) {
                throw new GATKException("unexpected CNV truth genotype makes reference to a non-existent allele: " + truthGT);
            }
            truthAC[CopyNumberTriStateAllele.valueOf(contextAlleles.get(truthAlleleIndex)).index()]++;
        }
        genotype.getAllele(0);
    }
    this.truthAC = truthAC;
    this.callsAC = callsAC;
}
 
Example 4
@Test
public void testAnnotateVariantContextWithPreprocessingValues() {
    final FeatureDataSource<VariantContext> featureDataSource = new FeatureDataSource<>(new File(smallM2Vcf));
    SortedSet<Transition> relevantTransitions = new TreeSet<>();
    relevantTransitions.add(Transition.transitionOf('G', 'T'));

    final Map<Transition, Double> preAdapterQFakeScoreMap = new HashMap<>();
    final double amGTPreAdapterQ = 20.0;
    preAdapterQFakeScoreMap.put(Transition.transitionOf('G', 'T'), amGTPreAdapterQ);  // preAdapterQ suppression will do nothing.

    for (final VariantContext vc : featureDataSource) {
        final VariantContext updatedVariantContext = OrientationBiasFilterer.annotateVariantContextWithPreprocessingValues(vc, relevantTransitions, preAdapterQFakeScoreMap);

        final Genotype genotypeTumor = updatedVariantContext.getGenotype("TUMOR");
        final Genotype genotypeNormal = updatedVariantContext.getGenotype("NORMAL");

        Assert.assertNotEquals(genotypeTumor.getExtendedAttribute(OrientationBiasFilterConstants.FOB, VCFConstants.EMPTY_ALLELE), VCFConstants.EMPTY_ALLELE);
        Assert.assertNotEquals(genotypeTumor.getExtendedAttribute(OrientationBiasFilterConstants.P_ARTIFACT_FIELD_NAME, VCFConstants.EMPTY_ALLELE), VCFConstants.EMPTY_ALLELE);

        assertArtifact(amGTPreAdapterQ, genotypeTumor, Transition.transitionOf('G', 'T'));

        // The NORMAL is always ref/ref in the example file.
        assertNormal(genotypeNormal);
    }
}
 
Example 5
@Test
public void testCalculatePosteriorSamplePlusExternal() {
    final VariantContext testOverlappingBase = makeVC("1", Arrays.asList(Aref,T), makeG("s1",T,T,40,20,0),
            makeG("s2",Aref,T,18,0,24),
            makeG("s3",Aref,T,22,0,12));
    final List<VariantContext> supplTest1 = new ArrayList<>(3);
    supplTest1.add(new VariantContextBuilder(makeVC("2", Arrays.asList(Aref,T))).attribute(GATKVCFConstants.MLE_ALLELE_COUNT_KEY,2).attribute(VCFConstants.ALLELE_NUMBER_KEY,10).make());
    supplTest1.add(new VariantContextBuilder(makeVC("3", Arrays.asList(Aref,T))).attribute(VCFConstants.ALLELE_COUNT_KEY,4).attribute(VCFConstants.ALLELE_NUMBER_KEY,22).make());
    supplTest1.add(makeVC("4", Arrays.asList(Aref,T),
            makeG("s_1",T,T),
            makeG("s_2",Aref,T)));
    final VariantContext test1result = PosteriorProbabilitiesUtils.calculatePosteriorProbs(testOverlappingBase,supplTest1,0,0.001,true,false, false);
    // the counts here are ref=30, alt=14
    final Genotype test1exp1 = makeGwithPLs("t1",T,T,new double[]{-3.370985, -1.415172, -0.01721766});
    final Genotype test1exp2 = makeGwithPLs("t2",Aref,T,new double[]{-1.763792, -0.007978791, -3.010024});
    final Genotype test1exp3 = makeGwithPLs("t3",Aref,T,new double[]{-2.165587, -0.009773643, -1.811819});
    Assert.assertEquals(arraysEq(test1exp1.getPL(),_mleparse((List<Integer>) test1result.getGenotype(0).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
    Assert.assertEquals(arraysEq(test1exp2.getPL(),_mleparse((List<Integer>) test1result.getGenotype(1).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
    Assert.assertEquals(arraysEq(test1exp3.getPL(),_mleparse((List<Integer>) test1result.getGenotype(2).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");

    final VariantContext testNonOverlapping = makeVC("1", Arrays.asList(Aref,T), makeG("s1",T,T,3,1,0));
    final List<VariantContext> other = Arrays.asList(makeVC("2", Arrays.asList(Aref,C),makeG("s2",C,C,10,2,0)));
    final VariantContext test2result = PosteriorProbabilitiesUtils.calculatePosteriorProbs(testNonOverlapping,other,0,0.001,true,false,false);
    final Genotype test2exp1 = makeGwithPLs("SGV",T,T,new double[]{-4.078345, -3.276502, -0.0002661066});
    Assert.assertEquals(arraysEq(test2exp1.getPL(),_mleparse((List<Integer>) test2result.getGenotype(0).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
}
 
Example 6
@Test
public void testGetNumReads(){
    final Allele refAllele = Allele.create("A", true);
    final Allele altAllele = Allele.create("T");
    final Genotype homRefDP = new GenotypeBuilder("sample1", Arrays.asList(refAllele, refAllele)).DP(5).make();
    final Genotype homRefMinDP = new  GenotypeBuilder("sample2", Arrays.asList(refAllele, refAllele))
            .DP(11)
            .attribute(GATKVCFConstants.MIN_DP_FORMAT_KEY, 10).make();
    final Genotype hetAlt = new GenotypeBuilder("sample3", Arrays.asList(refAllele, altAllele)).DP(100).make();
    final VariantContext vc = new VariantContextBuilder().alleles(Arrays.asList(refAllele, altAllele))
            .chr("1")
            .start(15L)
            .stop(15L)
            .genotypes(homRefDP, homRefMinDP, hetAlt)
            .attribute(VCFConstants.DEPTH_KEY, 5 + 10 + 100)
            .make();

    Assert.assertEquals(RMSMappingQuality.getNumOfReads(vc, null), 115 - 5 - 10);
}
 
Example 7
@Test(dataProvider = "updatePLsSACsAndADData")
public void testUpdatePLsAndADData(final VariantContext originalVC,
                                   final VariantContext selectedVC,
                                   final List<Genotype> expectedGenotypes) {
    // initialize cache of allele anyploid indices
    for (final Genotype genotype : originalVC.getGenotypes()) {
        GenotypeLikelihoods.initializeAnyploidPLIndexToAlleleIndices(originalVC.getNAlleles() - 1, genotype.getPloidy());
    }

    final VariantContext selectedVCwithGTs = new VariantContextBuilder(selectedVC).genotypes(originalVC.getGenotypes()).make();

    final GenotypesContext oldGs = selectedVCwithGTs.getGenotypes();
    final GenotypesContext actual = selectedVCwithGTs.getNAlleles() == originalVC.getNAlleles() ? oldGs :
                                    AlleleSubsettingUtils.subsetAlleles(oldGs, 0, originalVC.getAlleles(),
                                                                        selectedVCwithGTs.getAlleles(),
                                                                        GenotypeAssignmentMethod.DO_NOT_ASSIGN_GENOTYPES,
                                                                        originalVC.getAttributeAsInt(VCFConstants.DEPTH_KEY, 0));

    Assert.assertEquals(actual.size(), expectedGenotypes.size());
    for ( final Genotype expected : expectedGenotypes ) {
        final Genotype actualGT = actual.get(expected.getSampleName());
        Assert.assertNotNull(actualGT);
        VariantContextTestUtils.assertGenotypesAreEqual(actualGT, expected);
    }
}
 
Example 8
/**
 * Add physical phase information to the provided variant context
 *
 * @param vc   the variant context
 * @param ID   the ID to use
 * @param phaseGT the phase GT string to use
 * @return phased non-null variant context
 */
private static VariantContext phaseVC(final VariantContext vc, final String ID, final String phaseGT, final int phaseSetID) {
    final List<Genotype> phasedGenotypes = new ArrayList<>();
    for ( final Genotype g : vc.getGenotypes() ) {
        final List<Allele> alleles = g.getAlleles();
        if (phaseGT.equals(phase10) && g.isHet()) Collections.reverse(alleles); // swap the alleles if heterozygous
        final Genotype genotype = new GenotypeBuilder(g)
            .alleles(alleles)
            .phased(true)
            .attribute(GATKVCFConstants.HAPLOTYPE_CALLER_PHASING_ID_KEY, ID)
            .attribute(GATKVCFConstants.HAPLOTYPE_CALLER_PHASING_GT_KEY, phaseGT)
            .attribute(VCFConstants.PHASE_SET_KEY, phaseSetID)
            .make();
        phasedGenotypes.add(genotype);
    }
    return new VariantContextBuilder(vc).genotypes(phasedGenotypes).make();
}
 
Example 9
Source Project: picard   Source File: FingerprintUtils.java    License: MIT License 6 votes vote down vote up
private static VariantContextWriter getVariantContextWriter(final File outputFile,
                                                            final File referenceSequenceFileName,
                                                            final String sample,
                                                            final String source,
                                                            final ReferenceSequenceFile ref) {
    final VariantContextWriter variantContextWriter = new VariantContextWriterBuilder()
            .setReferenceDictionary(ref.getSequenceDictionary())
            .setOutputFile(outputFile).build();

    final Set<VCFHeaderLine> lines = new LinkedHashSet<>();
    lines.add(new VCFHeaderLine("reference", referenceSequenceFileName.getAbsolutePath()));
    lines.add(new VCFHeaderLine("source", source));
    lines.add(new VCFHeaderLine("fileDate", new Date().toString()));

    lines.add(VCFStandardHeaderLines.getFormatLine(VCFConstants.GENOTYPE_PL_KEY));
    lines.add(VCFStandardHeaderLines.getFormatLine(VCFConstants.GENOTYPE_ALLELE_DEPTHS));
    lines.add(VCFStandardHeaderLines.getFormatLine(VCFConstants.DEPTH_KEY));

    final VCFHeader header = new VCFHeader(lines, Collections.singletonList(sample));
    header.setSequenceDictionary(ref.getSequenceDictionary());
    variantContextWriter.writeHeader(header);
    return variantContextWriter;
}
 
Example 10
@Test
public void testLikelihoods(){

    final Allele REF = Allele.create("T", true);
    final Allele ALT = Allele.create("A");

    final int refDepth= 5;
    final int altDepth= 3;
    final int refMQ = 10;
    final int altMQ = 0;
    final List<GATKRead> refReads = IntStream.range(0, refDepth).mapToObj(i -> makeRead(refMQ)).collect(Collectors.toList());
    final List<GATKRead> altReads = IntStream.range(0, altDepth).mapToObj(i -> makeRead(altMQ)).collect(Collectors.toList());
    final AlleleLikelihoods<GATKRead, Allele> likelihoods =
            ArtificialAnnotationUtils.makeLikelihoods("sample1", refReads, altReads, -1.0, -1.0, REF, ALT);

    final VariantContext vc = makeVC();
    final ReferenceContext referenceContext= null;
    final Map<String, Object> annotate = new MappingQualityZero().annotate(referenceContext, vc, likelihoods);
    Assert.assertEquals(annotate.size(), 1, "size");
    Assert.assertEquals(annotate.keySet(), Collections.singleton(VCFConstants.MAPPING_QUALITY_ZERO_KEY), "annots");
    Assert.assertEquals(annotate.get(VCFConstants.MAPPING_QUALITY_ZERO_KEY), String.valueOf(altDepth));

}
 
Example 11
public void writeOutRecalibrationTable(final VariantContextWriter recalWriter, final SAMSequenceDictionary seqDictionary) {
    // we need to sort in coordinate order in order to produce a valid VCF
    Collections.sort( data, VariantDatum.getComparator(seqDictionary) );

    // create dummy alleles to be used
    List<Allele> alleles = Arrays.asList(Allele.create("N", true), Allele.create("<VQSR>", false));

    for( final VariantDatum datum : data ) {
        if (VRAC.useASannotations)
            alleles = Arrays.asList(datum.referenceAllele, datum.alternateAllele); //use the alleles to distinguish between multiallelics in AS mode
        VariantContextBuilder builder = new VariantContextBuilder("VQSR", datum.loc.getContig(), datum.loc.getStart(), datum.loc.getEnd(), alleles);
        builder.attribute(VCFConstants.END_KEY, datum.loc.getEnd());
        builder.attribute(GATKVCFConstants.VQS_LOD_KEY, String.format("%.4f", datum.lod));
        builder.attribute(GATKVCFConstants.CULPRIT_KEY, (datum.worstAnnotation != -1 ? annotationKeys.get(datum.worstAnnotation) : "NULL"));

        if ( datum.atTrainingSite ) builder.attribute(GATKVCFConstants.POSITIVE_LABEL_KEY, true);
        if ( datum.atAntiTrainingSite ) builder.attribute(GATKVCFConstants.NEGATIVE_LABEL_KEY, true);

        recalWriter.add(builder.make());
    }
}
 
Example 12
@Test
public void testCalculatePosteriorSamplePlusExternal() {
    final VariantContext testOverlappingBase = makeVC("1", Arrays.asList(Aref,T), makeG("s1",T,T,40,20,0),
            makeG("s2",Aref,T,18,0,24),
            makeG("s3",Aref,T,22,0,12));
    final List<VariantContext> supplTest1 = new ArrayList<>(3);
    supplTest1.add(new VariantContextBuilder(makeVC("2", Arrays.asList(Aref,T))).attribute(GATKVCFConstants.MLE_ALLELE_COUNT_KEY,2).attribute(VCFConstants.ALLELE_NUMBER_KEY,10).make());
    supplTest1.add(new VariantContextBuilder(makeVC("3", Arrays.asList(Aref,T))).attribute(VCFConstants.ALLELE_COUNT_KEY,4).attribute(VCFConstants.ALLELE_NUMBER_KEY,22).make());
    supplTest1.add(makeVC("4", Arrays.asList(Aref,T),
            makeG("s_1",T,T),
            makeG("s_2",Aref,T)));
    final VariantContext test1result = PosteriorProbabilitiesUtils.calculatePosteriorProbs(testOverlappingBase,supplTest1, 0, defaultOptions);
    // the counts here are ref=30, alt=14
    final Genotype test1exp1 = makeGwithLog10GLs("t1",T,T,new double[]{-3.370985, -1.415172, -0.01721766});
    final Genotype test1exp2 = makeGwithLog10GLs("t2",Aref,T,new double[]{-1.763792, -0.007978791, -3.010024});
    final Genotype test1exp3 = makeGwithLog10GLs("t3",Aref,T,new double[]{-2.165587, -0.009773643, -1.811819});
    Assert.assertEquals(arraysEq(test1exp1.getPL(),_mleparse((List<Integer>) test1result.getGenotype(0).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
    Assert.assertEquals(arraysEq(test1exp2.getPL(),_mleparse((List<Integer>) test1result.getGenotype(1).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
    Assert.assertEquals(arraysEq(test1exp3.getPL(),_mleparse((List<Integer>) test1result.getGenotype(2).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");

    final VariantContext testNonOverlapping = makeVC("1", Arrays.asList(Aref,T), makeG("s1",T,T,3,1,0));
    final List<VariantContext> other = Arrays.asList(makeVC("2", Arrays.asList(Aref,C),makeG("s2",C,C,10,2,0)));
    final VariantContext test2result = PosteriorProbabilitiesUtils.calculatePosteriorProbs(testNonOverlapping,other,0, defaultOptions);
    final Genotype test2exp1 = makeGwithLog10GLs("SGV",T,T,new double[]{-4.078345, -3.276502, -0.0002661066});
    Assert.assertEquals(arraysEq(test2exp1.getPL(),_mleparse((List<Integer>) test2result.getGenotype(0).getAnyAttribute(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY))), "");
}
 
Example 13
public static VariantContext makeHomRef(final String contig, final int start, final int GQ, final int end) {
    final VariantContextBuilder vcb = new VariantContextBuilder("test", contig, start, end, ALLELES);
    final GenotypeBuilder gb = new GenotypeBuilder(SAMPLE_NAME, Arrays.asList(REF, REF));
    gb.GQ(GQ);
    gb.DP(10);
    gb.AD(new int[]{1, 2});
    gb.PL(new int[]{0, 10, 100});
    vcb.attribute(VCFConstants.END_KEY, end);
    return vcb.genotypes(gb.make()).make();
}
 
Example 14
private static VCFHeader getHeader() {
    final Set<VCFHeaderLine> headerlines = new LinkedHashSet<>();
    VCFStandardHeaderLines.addStandardFormatLines(headerlines, true,
                                                  VCFConstants.GENOTYPE_KEY,
                                                  VCFConstants.GENOTYPE_QUALITY_KEY,
                                                  VCFConstants.GENOTYPE_PL_KEY, VCFConstants.DEPTH_KEY);
    final SAMSequenceDictionary dict = new SAMSequenceDictionary(
            Collections.singletonList(new SAMSequenceRecord("1", 100)));
    final VCFHeader header = new VCFHeader(headerlines, Collections.singleton(SAMPLE));
    header.setSequenceDictionary(dict);
    return header;
}
 
Example 15
@Test
public void testDescriptions() throws Exception {
    final RMSMappingQuality cov = new RMSMappingQuality();
    Assert.assertEquals(cov.getDescriptions().size(), 1);
    Assert.assertEquals(cov.getDescriptions().get(0).getID(), VCFConstants.RMS_MAPPING_QUALITY_KEY);
    Assert.assertEquals(cov.getRawDescriptions().size(), 1);
    Assert.assertEquals(cov.getRawDescriptions().get(0).getID(), GATKVCFConstants.RAW_MAPPING_QUALITY_WITH_DEPTH_KEY);
    RMSMappingQuality annotationClass = new RMSMappingQuality();
    Assert.assertEquals(annotationClass.getPrimaryRawKey(), GATKVCFConstants.RAW_MAPPING_QUALITY_WITH_DEPTH_KEY);
    Assert.assertEquals(annotationClass.getKeyNames(), Sets.newHashSet(VCFConstants.RMS_MAPPING_QUALITY_KEY, GATKVCFConstants.RAW_MAPPING_QUALITY_WITH_DEPTH_KEY));
}
 
Example 16
private static double[] attributeValueToDoubleArray(final Object value, final String key, final Supplier<double[]> defaultResult, final double missingValue) {
    Utils.nonNull(key);
    final ToDoubleFunction<Object> doubleConverter = o -> {
        if (o == null) {
            return missingValue;
        } else {
            final String s = String.valueOf(o);
            if (s.equals(VCFConstants.MISSING_VALUE_v4)) {
                return missingValue;
            } else {
                try {
                    return Double.parseDouble(s);
                } catch (final NumberFormatException ex) {
                    throw new GATKException(String.format("INFO annotation '%s' contains a non-double value '%s'", key, s), ex);
                }
            }
        }
    };

    if (value == null) {
        return defaultResult.get();
    } else if (value.getClass().isArray()) {
        final double[] result = new double[Array.getLength(value)];
        for (int i = 0; i < result.length; i++) {
            result[i] = doubleConverter.applyAsDouble(String.valueOf(Array.get(value, i)));
        }
        return result;
    } else if (value.getClass().isAssignableFrom(Iterable.class)) {
        return StreamSupport.stream(((Iterable<?>)value).spliterator(), false)
                .mapToDouble(doubleConverter).toArray();
    } else { // as a last resort with transform it into an String and try to parse an array out of it.
        return Stream.of(String.valueOf(value).trim().replaceAll("\\[|\\]", "")
                .split(VCFConstants.INFO_FIELD_ARRAY_SEPARATOR))
                .mapToDouble(doubleConverter).toArray();
    }
}
 
Example 17
private static int[] attributeValueToIntArray(final Object value, final String key, final Supplier<int[]> defaultResult, final int missingValue) {
    Utils.nonNull(key);
    final ToIntFunction<Object> intConverter = o -> {
        if (o == null) {
            return missingValue;
        } else {
            final String s = String.valueOf(o);
            if (s.equals(VCFConstants.MISSING_VALUE_v4)) {
                return missingValue;
            } else {
                try {
                    return Integer.parseInt(s);
                } catch (final NumberFormatException ex) {
                    throw new GATKException(String.format("INFO annotation '%s' contains a non-int value '%s'", key, s), ex);
                }
            }
        }
    };

    if (value == null) {
        return defaultResult.get();
    } else if (value.getClass().isArray()) {
        final int[] result = new int[Array.getLength(value)];
        for (int i = 0; i < result.length; i++) {
            result[i] = intConverter.applyAsInt(String.valueOf(Array.get(value, i)));
        }
        return result;
    } else if (value.getClass().isAssignableFrom(Iterable.class)) {
        return StreamSupport.stream(((Iterable<?>)value).spliterator(), false)
                .mapToInt(intConverter).toArray();
    } else { // as a last resort with transform it into an String and try to parse an array out of it.
        return Stream.of(String.valueOf(value).trim().replaceAll("\\[|\\]", "")
                .split(VCFConstants.INFO_FIELD_ARRAY_SEPARATOR))
                .mapToInt(intConverter).toArray();
    }
}
 
Example 18
@Test
public void testNullLikelihoods() throws Exception {
    final VariantContext vc= makeVC();
    final ReferenceContext referenceContext= null;
    final InfoFieldAnnotation cov = new MappingQualityZero();
    final Map<String, Object> annotate = cov.annotate(referenceContext, vc, null);
    Assert.assertTrue(annotate.isEmpty());

    Assert.assertEquals(cov.getDescriptions().size(), 1);
    Assert.assertEquals(cov.getDescriptions().get(0).getID(), VCFConstants.MAPPING_QUALITY_ZERO_KEY);
}
 
Example 19
public String getFilterString() {
    if (filters.isEmpty()) {
        return VCFConstants.PASSES_FILTERS_v4;
    } else {
        return filters.stream().sorted().map(Enum<EvaluationFilter>::name)
                .collect(Collectors.joining(VCFConstants.FILTER_CODE_SEPARATOR));
    }
}
 
Example 20
static final VariantContextBuilder makeInsertion(final String chr, final int pos, final int end, final int svLen,
                                                 final Allele refAllele) {

    return new VariantContextBuilder().chr(chr).start(pos).stop(end)
            .alleles(Arrays.asList(refAllele, INS_SYMB_ALLELE))
            .id(makeID(SimpleSVType.SupportedType.INS.name(), chr, pos, chr, end, ""))
            .attribute(VCFConstants.END_KEY, end)
            .attribute(SVLEN, svLen)
            .attribute(SVTYPE, SimpleSVType.SupportedType.INS.name());
}
 
Example 21
private static boolean isMissing(final Object value) {
    if ( value == null ) { return true; }
    else if ( value.equals(VCFConstants.MISSING_VALUE_v4) ) { return true; }
    else if ( value instanceof List ) {
        // handles the case where all elements are null or the list is empty
        for ( final Object elt : (List)value) {
            if (elt != null) {
                return false;
            }
        }
        return true;
    } else {
        return false;
    }
}
 
Example 22
@DataProvider(name = "GoodBandPartitionData")
public Object[][] makeBandPartitionData() {
    return new Object[][]{
            {Collections.singletonList(1), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,MAX_GENOTYPE_QUAL+1))},
            {Arrays.asList(1, 2, 3), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,2), Range.closedOpen(2,3),Range.closedOpen(3,MAX_GENOTYPE_QUAL+1))},
            {Arrays.asList(1, 10), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,10), Range.closedOpen(10, MAX_GENOTYPE_QUAL+1))},
            {Arrays.asList(1, 10, 30), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,10), Range.closedOpen(10, 30),Range.closedOpen(30,MAX_GENOTYPE_QUAL+1))},
            {Arrays.asList(1, 10, MAX_GENOTYPE_QUAL - 1), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,10), Range.closedOpen(10, MAX_GENOTYPE_QUAL - 1),Range.closedOpen(MAX_GENOTYPE_QUAL - 1,MAX_GENOTYPE_QUAL+1))},
            {Arrays.asList(1, 10, MAX_GENOTYPE_QUAL), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,10), Range.closedOpen(10, MAX_GENOTYPE_QUAL), Range.closedOpen(MAX_GENOTYPE_QUAL, MAX_GENOTYPE_QUAL+1))},
            {Arrays.asList(1, 10, MAX_GENOTYPE_QUAL + 1), Arrays.asList(Range.closedOpen(0,1), Range.closedOpen(1,10), Range.closedOpen(10, MAX_GENOTYPE_QUAL+1))},
            {Collections.singletonList(VCFConstants.MAX_GENOTYPE_QUAL + 1), Arrays.asList(Range.closedOpen(0,MAX_GENOTYPE_QUAL+1))}
    };
}
 
Example 23
Source Project: gatk   Source File: SVVCFWriter.java    License: BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
@VisibleForTesting
static VCFHeader getVcfHeader(final SAMSequenceDictionary referenceSequenceDictionary) {
    final Set<VCFHeaderLine> headerLines = new HashSet<>(GATKSVVCFHeaderLines.getSymbAltAlleleLines());
    headerLines.addAll(GATKSVVCFHeaderLines.getInfoLines());
    headerLines.add(VCFStandardHeaderLines.getInfoLine(VCFConstants.END_KEY));
    headerLines.addAll(GATKSVVCFHeaderLines.getFormatLines());
    headerLines.addAll(GATKSVVCFHeaderLines.getFilterLines());
    final VCFHeader header = new VCFHeader(new VCFHeader( headerLines ));
    header.setSequenceDictionary(referenceSequenceDictionary);
    return header;
}
 
Example 24
private boolean alleleFrequencyInRange(final VariantContext vc) {
    if (!vc.hasAttribute(VCFConstants.ALLELE_FREQUENCY_KEY)) {
        return false;
    } else {
        final double alleleFrequency = vc.getAttributeAsDouble(VCFConstants.ALLELE_FREQUENCY_KEY, -1.0);
        return minPopulationAlleleFrequency < alleleFrequency && alleleFrequency < maxPopulationAlleleFrequency;
    }
}
 
Example 25
public PileupSummary(final VariantContext vc, final ReadPileup pileup) {
    contig = vc.getContig();
    position = vc.getStart();
    alleleFrequency = vc.getAttributeAsDouble(VCFConstants.ALLELE_FREQUENCY_KEY, 0);
    final byte altBase = vc.getAlternateAllele(0).getBases()[0];
    final byte refBase = vc.getReference().getBases()[0];
    final int[] baseCounts = pileup.getBaseCounts();
    altCount = baseCounts[BaseUtils.simpleBaseToBaseIndex(altBase)];
    refCount = baseCounts[BaseUtils.simpleBaseToBaseIndex(refBase)];
    totalCount = (int) MathUtils.sum(baseCounts);
    otherAltsCount = totalCount - altCount - refCount;
}
 
Example 26
/**
 * lookup the depth from the VC DP field or calculate by summing the depths of the genotypes
 */
private static int calculateVCDepth(VariantContext vc) {
    if ( vc.hasAttribute(VCFConstants.DEPTH_KEY) ) {
        return vc.getAttributeAsInt(VCFConstants.DEPTH_KEY, 0);
    } else { // handle the gVCF case from the HaplotypeCaller
        return vc.getGenotypes().stream()
                .mapToInt(ReferenceConfidenceVariantContextMerger::getBestDepthValue)
                .sum();
    }
}
 
Example 27
/**
 * Note that {@code delRange} is expected to be pre-process to VCF spec compatible,
 * e.g. if chr1:101-200 is deleted, then {@code delRange} should be chr1:100-200
 * @param delRange
 */
@VisibleForTesting
static VariantContextBuilder makeDeletion(final SimpleInterval delRange, final Allele refAllele) {

    return new VariantContextBuilder()
            .chr(delRange.getContig()).start(delRange.getStart()).stop(delRange.getEnd())
            .alleles(Arrays.asList(refAllele, altSymbAlleleDel))
            .id(makeID(SimpleSVType.SupportedType.DEL.name(), delRange.getContig(), delRange.getStart(), delRange.getEnd()))
            .attribute(VCFConstants.END_KEY, delRange.getEnd())
            .attribute(SVLEN, - delRange.size() + 1)
            .attribute(SVTYPE, SimpleSVType.SupportedType.DEL.name());
}
 
Example 28
/**
 * Remove the stale attributes from the merged set
 *
 * @param attributes the attribute map
 */
private static void removeStaleAttributesAfterMerge(final Map<String, Object> attributes) {
    attributes.remove(VCFConstants.ALLELE_COUNT_KEY);
    attributes.remove(VCFConstants.ALLELE_FREQUENCY_KEY);
    attributes.remove(VCFConstants.ALLELE_NUMBER_KEY);
    attributes.remove(GATKVCFConstants.MLE_ALLELE_COUNT_KEY);
    attributes.remove(GATKVCFConstants.MLE_ALLELE_FREQUENCY_KEY);
    attributes.remove(VCFConstants.END_KEY);
}
 
Example 29
/**
 * Writes an appropriate VCF header, given our arguments, to the provided writer
 *
 * @param vcfWriter writer to which the header should be written
 */
public void writeHeader( final VariantContextWriter vcfWriter, final SAMSequenceDictionary sequenceDictionary,
                         final Set<VCFHeaderLine>  defaultToolHeaderLines) {
    Utils.nonNull(vcfWriter);

    final Set<VCFHeaderLine> headerInfo = new HashSet<>();
    headerInfo.addAll(defaultToolHeaderLines);

    headerInfo.addAll(genotypingEngine.getAppropriateVCFInfoHeaders());
    // all annotation fields from VariantAnnotatorEngine
    headerInfo.addAll(annotationEngine.getVCFAnnotationDescriptions());
    // all callers need to add these standard annotation header lines
    headerInfo.add(GATKVCFHeaderLines.getInfoLine(GATKVCFConstants.DOWNSAMPLED_KEY));
    headerInfo.add(GATKVCFHeaderLines.getInfoLine(GATKVCFConstants.MLE_ALLELE_COUNT_KEY));
    headerInfo.add(GATKVCFHeaderLines.getInfoLine(GATKVCFConstants.MLE_ALLELE_FREQUENCY_KEY));
    // all callers need to add these standard FORMAT field header lines
    VCFStandardHeaderLines.addStandardFormatLines(headerInfo, true,
            VCFConstants.GENOTYPE_KEY,
            VCFConstants.GENOTYPE_QUALITY_KEY,
            VCFConstants.DEPTH_KEY,
            VCFConstants.GENOTYPE_PL_KEY);

    if ( ! hcArgs.doNotRunPhysicalPhasing ) {
        headerInfo.add(GATKVCFHeaderLines.getFormatLine(GATKVCFConstants.HAPLOTYPE_CALLER_PHASING_ID_KEY));
        headerInfo.add(GATKVCFHeaderLines.getFormatLine(GATKVCFConstants.HAPLOTYPE_CALLER_PHASING_GT_KEY));
    }

    // FILTER fields are added unconditionally as it's not always 100% certain the circumstances
    // where the filters are used.  For example, in emitting all sites the lowQual field is used
    headerInfo.add(GATKVCFHeaderLines.getFilterLine(GATKVCFConstants.LOW_QUAL_FILTER_NAME));

    if ( emitReferenceConfidence() ) {
        headerInfo.addAll(referenceConfidenceModel.getVCFHeaderLines());
    }

    final VCFHeader vcfHeader = new VCFHeader(headerInfo, sampleSet);
    vcfHeader.setSequenceDictionary(sequenceDictionary);
    vcfWriter.writeHeader(vcfHeader);
}
 
Example 30
private void applyContaminationFilter(final M2FiltersArgumentCollection MTFAC, final VariantContext vc, final Collection<String> filters) {
    final Genotype tumorGenotype = vc.getGenotype(tumorSample);
    final double[] alleleFractions = GATKProtectedVariantContextUtils.getAttributeAsDoubleArray(tumorGenotype, VCFConstants.ALLELE_FREQUENCY_KEY,
            () -> new double[] {1.0}, 1.0);
    final double maxFraction = MathUtils.arrayMax(alleleFractions);
    if (maxFraction < contamination) {
        filters.add(CONTAMINATION_FILTER_NAME);
    }
}