Search in sources :

Example 46 with UserException

use of org.broadinstitute.hellbender.exceptions.UserException in project gatk by broadinstitute.

the class MergeVcfs method doWork.

@Override
protected Object doWork() {
    final ProgressLogger progress = new ProgressLogger(logger, 10000);
    final List<String> sampleList = new ArrayList<>();
    final Collection<CloseableIterator<VariantContext>> iteratorCollection = new ArrayList<>(INPUT.size());
    final Collection<VCFHeader> headers = new HashSet<>(INPUT.size());
    VariantContextComparator variantContextComparator = null;
    SAMSequenceDictionary sequenceDictionary = null;
    if (SEQUENCE_DICTIONARY != null) {
        sequenceDictionary = SamReaderFactory.makeDefault().referenceSequence(REFERENCE_SEQUENCE).open(SEQUENCE_DICTIONARY).getFileHeader().getSequenceDictionary();
    }
    for (final File file : INPUT) {
        IOUtil.assertFileIsReadable(file);
        final VCFFileReader fileReader = new VCFFileReader(file, false);
        final VCFHeader fileHeader = fileReader.getFileHeader();
        if (variantContextComparator == null) {
            variantContextComparator = fileHeader.getVCFRecordComparator();
        } else {
            if (!variantContextComparator.isCompatible(fileHeader.getContigLines())) {
                throw new IllegalArgumentException("The contig entries in input file " + file.getAbsolutePath() + " are not compatible with the others.");
            }
        }
        if (sequenceDictionary == null)
            sequenceDictionary = fileHeader.getSequenceDictionary();
        if (sampleList.isEmpty()) {
            sampleList.addAll(fileHeader.getSampleNamesInOrder());
        } else {
            if (!sampleList.equals(fileHeader.getSampleNamesInOrder())) {
                throw new IllegalArgumentException("Input file " + file.getAbsolutePath() + " has sample entries that don't match the other files.");
            }
        }
        headers.add(fileHeader);
        iteratorCollection.add(fileReader.iterator());
    }
    if (CREATE_INDEX && sequenceDictionary == null) {
        throw new UserException("A sequence dictionary must be available (either through the input file or by setting it explicitly) when creating indexed output.");
    }
    final VariantContextWriterBuilder builder = new VariantContextWriterBuilder().setOutputFile(OUTPUT).setReferenceDictionary(sequenceDictionary).clearOptions();
    if (CREATE_INDEX) {
        builder.setOption(Options.INDEX_ON_THE_FLY);
    }
    try (final VariantContextWriter writer = builder.build()) {
        writer.writeHeader(new VCFHeader(VCFUtils.smartMergeHeaders(headers, false), sampleList));
        final MergingIterator<VariantContext> mergingIterator = new MergingIterator<>(variantContextComparator, iteratorCollection);
        while (mergingIterator.hasNext()) {
            final VariantContext context = mergingIterator.next();
            writer.add(context);
            progress.record(context.getContig(), context.getStart());
        }
        CloserUtil.close(mergingIterator);
    }
    return null;
}
Also used : CloseableIterator(htsjdk.samtools.util.CloseableIterator) ArrayList(java.util.ArrayList) VCFFileReader(htsjdk.variant.vcf.VCFFileReader) VariantContext(htsjdk.variant.variantcontext.VariantContext) ProgressLogger(org.broadinstitute.hellbender.utils.runtime.ProgressLogger) VariantContextComparator(htsjdk.variant.variantcontext.VariantContextComparator) SAMSequenceDictionary(htsjdk.samtools.SAMSequenceDictionary) MergingIterator(htsjdk.samtools.util.MergingIterator) VariantContextWriterBuilder(htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder) UserException(org.broadinstitute.hellbender.exceptions.UserException) VariantContextWriter(htsjdk.variant.variantcontext.writer.VariantContextWriter) VCFHeader(htsjdk.variant.vcf.VCFHeader) File(java.io.File) HashSet(java.util.HashSet)

Example 47 with UserException

use of org.broadinstitute.hellbender.exceptions.UserException in project gatk by broadinstitute.

the class ApplyVQSR method trancheIntervalIsValid.

private boolean trancheIntervalIsValid(final String sensitivityLimits) {
    final String[] vals = sensitivityLimits.split("to");
    if (vals.length != 2)
        return false;
    try {
        double lowerLimit = Double.parseDouble(vals[0]);
        //why does our last tranche end with 100+? Is there anything greater than 100 percent?  Really???
        double upperLimit = Double.parseDouble(vals[1].replace("+", ""));
    } catch (NumberFormatException e) {
        throw new UserException("Poorly formatted tranche filter name does not contain two sensitivity interval end points.");
    }
    return true;
}
Also used : UserException(org.broadinstitute.hellbender.exceptions.UserException)

Example 48 with UserException

use of org.broadinstitute.hellbender.exceptions.UserException in project gatk by broadinstitute.

the class ApplyVQSR method doSiteSpecificFiltering.

/**
     * Calculate the filter status for a given VariantContext using the combined data from all alleles at a site
     * @param vc
     * @param recals
     * @param builder   is modified by adding attributes
     * @return a String with the filter status for this site
     */
private String doSiteSpecificFiltering(final VariantContext vc, final List<VariantContext> recals, final VariantContextBuilder builder) {
    VariantContext recalDatum = getMatchingRecalVC(vc, recals, null);
    if (recalDatum == null) {
        throw new UserException("Encountered input variant which isn't found in the input recal file. Please make sure VariantRecalibrator and ApplyRecalibration were run on the same set of input variants. First seen at: " + vc);
    }
    final String lodString = recalDatum.getAttributeAsString(GATKVCFConstants.VQS_LOD_KEY, null);
    if (lodString == null) {
        throw new UserException("Encountered a malformed record in the input recal file. There is no lod for the record at: " + vc);
    }
    final double lod;
    try {
        lod = Double.valueOf(lodString);
    } catch (NumberFormatException e) {
        throw new UserException("Encountered a malformed record in the input recal file. The lod is unreadable for the record at: " + vc);
    }
    builder.attribute(GATKVCFConstants.VQS_LOD_KEY, lod);
    builder.attribute(GATKVCFConstants.CULPRIT_KEY, recalDatum.getAttribute(GATKVCFConstants.CULPRIT_KEY));
    if (recalDatum != null) {
        if (recalDatum.hasAttribute(GATKVCFConstants.POSITIVE_LABEL_KEY))
            builder.attribute(GATKVCFConstants.POSITIVE_LABEL_KEY, true);
        if (recalDatum.hasAttribute(GATKVCFConstants.NEGATIVE_LABEL_KEY))
            builder.attribute(GATKVCFConstants.NEGATIVE_LABEL_KEY, true);
    }
    return generateFilterString(lod);
}
Also used : VariantContext(htsjdk.variant.variantcontext.VariantContext) UserException(org.broadinstitute.hellbender.exceptions.UserException)

Example 49 with UserException

use of org.broadinstitute.hellbender.exceptions.UserException in project gatk by broadinstitute.

the class Concordance method onTraversalSuccess.

@Override
public Object onTraversalSuccess() {
    try (ConcordanceSummaryRecord.Writer concordanceSummaryWriter = ConcordanceSummaryRecord.getWriter(summary)) {
        concordanceSummaryWriter.writeRecord(new ConcordanceSummaryRecord(VariantContext.Type.SNP, snpCounts.get(ConcordanceState.TRUE_POSITIVE).longValue(), snpCounts.get(ConcordanceState.FALSE_POSITIVE).longValue(), snpCounts.get(ConcordanceState.FALSE_NEGATIVE).longValue() + snpCounts.get(ConcordanceState.FILTERED_FALSE_NEGATIVE).longValue()));
        concordanceSummaryWriter.writeRecord(new ConcordanceSummaryRecord(VariantContext.Type.INDEL, indelCounts.get(ConcordanceState.TRUE_POSITIVE).longValue(), indelCounts.get(ConcordanceState.FALSE_POSITIVE).longValue(), indelCounts.get(ConcordanceState.FALSE_NEGATIVE).longValue() + indelCounts.get(ConcordanceState.FILTERED_FALSE_NEGATIVE).longValue()));
    } catch (IOException e) {
        throw new UserException("Encountered an IO exception writing the concordance summary table", e);
    }
    if (truePositivesAndFalsePositivesVcfWriter != null) {
        truePositivesAndFalsePositivesVcfWriter.close();
    }
    if (truePositivesAndFalseNegativesVcfWriter != null) {
        truePositivesAndFalseNegativesVcfWriter.close();
    }
    if (filteredTrueNegativesAndFalseNegativesVcfWriter != null) {
        filteredTrueNegativesAndFalseNegativesVcfWriter.close();
    }
    return "SUCCESS";
}
Also used : IOException(java.io.IOException) UserException(org.broadinstitute.hellbender.exceptions.UserException)

Example 50 with UserException

use of org.broadinstitute.hellbender.exceptions.UserException in project gatk by broadinstitute.

the class CalculateGenotypePosteriors method onTraversalStart.

@Override
public void onTraversalStart() {
    final VariantContextWriterBuilder builder = new VariantContextWriterBuilder().setOutputFile(out).setOutputFileType(VariantContextWriterBuilder.OutputType.VCF);
    if (hasReference()) {
        vcfWriter = builder.setReferenceDictionary(getBestAvailableSequenceDictionary()).setOption(Options.INDEX_ON_THE_FLY).build();
    } else {
        vcfWriter = builder.unsetOption(Options.INDEX_ON_THE_FLY).build();
        logger.info("Can't make an index for output file " + out + " because a reference dictionary is required for creating Tribble indices on the fly");
    }
    sampleDB = initializeSampleDB();
    // Get list of samples to include in the output
    final Map<String, VCFHeader> vcfHeaders = Collections.singletonMap(getDrivingVariantsFeatureInput().getName(), getHeaderForVariants());
    final Set<String> vcfSamples = VcfUtils.getSortedSampleSet(vcfHeaders, GATKVariantContextUtils.GenotypeMergeType.REQUIRE_UNIQUE);
    //Get the trios from the families passed as ped
    if (!skipFamilyPriors) {
        final Set<Trio> trios = sampleDB.getTrios();
        if (trios.isEmpty()) {
            logger.info("No PED file passed or no *non-skipped* trios found in PED file. Skipping family priors.");
            skipFamilyPriors = true;
        }
    }
    final VCFHeader header = vcfHeaders.values().iterator().next();
    if (!header.hasGenotypingData()) {
        throw new UserException("VCF has no genotypes");
    }
    if (header.hasInfoLine(GATKVCFConstants.MLE_ALLELE_COUNT_KEY)) {
        final VCFInfoHeaderLine mleLine = header.getInfoHeaderLine(GATKVCFConstants.MLE_ALLELE_COUNT_KEY);
        if (mleLine.getCountType() != VCFHeaderLineCount.A) {
            throw new UserException("VCF does not have a properly formatted MLEAC field: the count type should be \"A\"");
        }
        if (mleLine.getType() != VCFHeaderLineType.Integer) {
            throw new UserException("VCF does not have a properly formatted MLEAC field: the field type should be \"Integer\"");
        }
    }
    // Initialize VCF header
    final Set<VCFHeaderLine> headerLines = VCFUtils.smartMergeHeaders(vcfHeaders.values(), true);
    headerLines.add(GATKVCFHeaderLines.getFormatLine(GATKVCFConstants.PHRED_SCALED_POSTERIORS_KEY));
    headerLines.add(GATKVCFHeaderLines.getInfoLine(GATKVCFConstants.GENOTYPE_PRIOR_KEY));
    if (!skipFamilyPriors) {
        headerLines.add(GATKVCFHeaderLines.getFormatLine(GATKVCFConstants.JOINT_LIKELIHOOD_TAG_NAME));
        headerLines.add(GATKVCFHeaderLines.getFormatLine(GATKVCFConstants.JOINT_POSTERIOR_TAG_NAME));
    }
    headerLines.addAll(getDefaultToolVCFHeaderLines());
    vcfWriter.writeHeader(new VCFHeader(headerLines, vcfSamples));
    final Map<String, Set<Sample>> families = sampleDB.getFamilies(vcfSamples);
    famUtils = new FamilyLikelihoods(sampleDB, deNovoPrior, vcfSamples, families);
}
Also used : VariantContextWriterBuilder(htsjdk.variant.variantcontext.writer.VariantContextWriterBuilder) UserException(org.broadinstitute.hellbender.exceptions.UserException)

Aggregations

UserException (org.broadinstitute.hellbender.exceptions.UserException)100 File (java.io.File)30 IOException (java.io.IOException)30 SAMSequenceDictionary (htsjdk.samtools.SAMSequenceDictionary)14 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)11 SamReader (htsjdk.samtools.SamReader)10 VariantContext (htsjdk.variant.variantcontext.VariantContext)10 ProgressLogger (org.broadinstitute.hellbender.utils.runtime.ProgressLogger)10 SAMRecord (htsjdk.samtools.SAMRecord)9 IntervalList (htsjdk.samtools.util.IntervalList)9 List (java.util.List)9 SAMFileHeader (htsjdk.samtools.SAMFileHeader)8 ReferenceSequenceFileWalker (htsjdk.samtools.reference.ReferenceSequenceFileWalker)8 SamLocusIterator (htsjdk.samtools.util.SamLocusIterator)8 LogManager (org.apache.logging.log4j.LogManager)8 Logger (org.apache.logging.log4j.Logger)8 GATKException (org.broadinstitute.hellbender.exceptions.GATKException)7 MetricsFile (htsjdk.samtools.metrics.MetricsFile)6 SamRecordFilter (htsjdk.samtools.filter.SamRecordFilter)5 FileNotFoundException (java.io.FileNotFoundException)5