Search in sources :

Example 11 with PCACoveragePoN

use of org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN in project gatk-protected by broadinstitute.

the class NormalizeSomaticReadCounts method doWork.

@Override
protected Object doWork() {
    if (!new HDF5Library().load(null)) {
        //Note: passing null means using the default temp dir.
        throw new UserException.HardwareFeatureException("Cannot load the required HDF5 library. " + "HDF5 is currently supported on x86-64 architecture and Linux or OSX systems.");
    }
    IOUtils.canReadFile(ponFile);
    try (final HDF5File hdf5PoNFile = new HDF5File(ponFile)) {
        final PCACoveragePoN pon = new HDF5PCACoveragePoN(hdf5PoNFile, logger);
        final TargetCollection<Target> targetCollection = readTargetCollection(targetFile);
        final ReadCountCollection proportionalCoverageProfile = readInputReadCounts(readCountsFile, targetCollection);
        final PCATangentNormalizationResult tangentNormalizationResult = pon.normalize(proportionalCoverageProfile);
        ;
        tangentNormalizationResult.write(getCommandLine(), tangentNormalizationOutFile, preTangentNormalizationOutFile, betaHatsOutFile, fntOutFile);
        return "SUCCESS";
    }
}
Also used : HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) HDF5Library(org.broadinstitute.hdf5.HDF5Library) HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN) PCATangentNormalizationResult(org.broadinstitute.hellbender.tools.pon.coverage.pca.PCATangentNormalizationResult) HDF5File(org.broadinstitute.hdf5.HDF5File)

Example 12 with PCACoveragePoN

use of org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN in project gatk by broadinstitute.

the class NormalizeSomaticReadCountsIntegrationTest method assertFactorNormalizedValues.

private void assertFactorNormalizedValues(final ReadCountCollection input, final ReadCountCollection factorNormalized) {
    try (final HDF5File ponReader = new HDF5File(TEST_PON)) {
        final PCACoveragePoN pon = new HDF5PCACoveragePoN(ponReader);
        final double[] targetFactors = pon.getTargetFactors();
        final List<String> ponTargets = pon.getTargetNames();
        final Map<String, Integer> ponTargetIndexes = new HashMap<>(ponTargets.size());
        for (int i = 0; i < ponTargets.size(); i++) {
            ponTargetIndexes.put(ponTargets.get(i), i);
        }
        final RealMatrix inputCounts = input.counts();
        final RealMatrix factorNormalizedCounts = factorNormalized.counts();
        for (int i = 0; i < factorNormalizedCounts.getRowDimension(); i++) {
            final double factor = targetFactors[ponTargetIndexes.get(factorNormalized.targets().get(i).getName())];
            final double[] inputValues = inputCounts.getRow(i);
            final double[] outputValues = factorNormalizedCounts.getRow(i);
            for (int j = 0; j < inputValues.length; j++) {
                final double expected = inputValues[j] / factor;
                Assert.assertEquals(outputValues[j], expected, 0.0000001, "" + i + " , " + j);
            }
        }
    }
}
Also used : HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN) HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) HDF5File(org.broadinstitute.hdf5.HDF5File)

Example 13 with PCACoveragePoN

use of org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN in project gatk by broadinstitute.

the class CreatePanelOfNormals method runPipeline.

@Override
protected void runPipeline(final JavaSparkContext ctx) {
    if (!new HDF5Library().load(null)) {
        //Note: passing null means using the default temp dir.
        throw new UserException.HardwareFeatureException("Cannot load the required HDF5 library. " + "HDF5 is currently supported on x86-64 architecture and Linux or OSX systems.");
    }
    if (blacklistOutFile == null) {
        blacklistOutFile = new File(outFile + BLACKLIST_FILE_APPEND);
    }
    if (targetWeightsOutFile == null) {
        targetWeightsOutFile = new File(outFile + TARGET_WEIGHTS_FILE_APPEND);
    }
    // Check parameters and load values to meet the backend PoN creation interface
    validateArguments();
    final TargetCollection<Target> targets = targetArguments.readTargetCollection(true);
    final OptionalInt numberOfEigensamples = parseNumberOfEigensamples(numberOfEigensamplesString);
    // Create the PoN, including QC, if specified.
    if (!isNoQc && !dryRun) {
        logger.info("QC:  Beginning creation of QC PoN...");
        final File outputQCFile = IOUtils.createTempFile("qc-pon-", ".hd5");
        HDF5PCACoveragePoNCreationUtils.create(ctx, outputQCFile, HDF5File.OpenMode.READ_WRITE, inputFile, targets, new ArrayList<>(), targetFactorThreshold, maximumPercentZerosInColumn, maximumPercentZerosInTarget, columnExtremeThresholdPercentile, outlierTruncatePercentileThresh, OptionalInt.of(NUM_QC_EIGENSAMPLES), dryRun);
        logger.info("QC:  QC PoN created...");
        logger.info("QC:  Collecting suspicious samples...");
        try (final HDF5File ponReader = new HDF5File(outputQCFile, HDF5File.OpenMode.READ_ONLY)) {
            final PCACoveragePoN qcPoN = new HDF5PCACoveragePoN(ponReader);
            final List<String> failingSampleNames = CoveragePoNQCUtils.retrieveSamplesWithArmLevelEvents(qcPoN, ctx);
            ParamUtils.writeStringListToFile(failingSampleNames, blacklistOutFile);
            // If no suspicious samples were found, just redo the PoN reduction to save time.
            if (failingSampleNames.size() != 0) {
                logger.info("QC:  Suspicious sample list created...");
                logger.info("Creating final PoN with " + failingSampleNames.size() + " suspicious samples removed...");
                HDF5PCACoveragePoNCreationUtils.create(ctx, outFile, HDF5File.OpenMode.CREATE, inputFile, targets, failingSampleNames, targetFactorThreshold, maximumPercentZerosInColumn, maximumPercentZerosInTarget, columnExtremeThresholdPercentile, outlierTruncatePercentileThresh, numberOfEigensamples, dryRun);
            } else {
                logger.info("QC:  No suspicious samples found ...");
                logger.info("Creating final PoN only redo'ing the reduction step ...");
                HDF5PCACoveragePoNCreationUtils.redoReduction(ctx, numberOfEigensamples, outputQCFile, outFile, HDF5File.OpenMode.CREATE);
            }
        }
    } else {
        logger.info("Creating PoN directly (skipping QC)...");
        HDF5PCACoveragePoNCreationUtils.create(ctx, outFile, HDF5File.OpenMode.CREATE, inputFile, targets, new ArrayList<>(), targetFactorThreshold, maximumPercentZerosInColumn, maximumPercentZerosInTarget, columnExtremeThresholdPercentile, outlierTruncatePercentileThresh, numberOfEigensamples, dryRun);
    }
    if (!dryRun) {
        logger.info("Writing target weights file to " + targetWeightsOutFile + "...");
        writeTargetWeightsFile(outFile, targetWeightsOutFile);
    }
    logger.info("Done...");
}
Also used : HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) HDF5Library(org.broadinstitute.hdf5.HDF5Library) HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN) OptionalInt(java.util.OptionalInt) HDF5File(org.broadinstitute.hdf5.HDF5File) File(java.io.File) HDF5File(org.broadinstitute.hdf5.HDF5File)

Example 14 with PCACoveragePoN

use of org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN in project gatk-protected by broadinstitute.

the class NormalizeSomaticReadCountsIntegrationTest method assertFactorNormalizedValues.

private void assertFactorNormalizedValues(final ReadCountCollection input, final ReadCountCollection factorNormalized) {
    try (final HDF5File ponReader = new HDF5File(TEST_PON)) {
        final PCACoveragePoN pon = new HDF5PCACoveragePoN(ponReader);
        final double[] targetFactors = pon.getTargetFactors();
        final List<String> ponTargets = pon.getTargetNames();
        final Map<String, Integer> ponTargetIndexes = new HashMap<>(ponTargets.size());
        for (int i = 0; i < ponTargets.size(); i++) {
            ponTargetIndexes.put(ponTargets.get(i), i);
        }
        final RealMatrix inputCounts = input.counts();
        final RealMatrix factorNormalizedCounts = factorNormalized.counts();
        for (int i = 0; i < factorNormalizedCounts.getRowDimension(); i++) {
            final double factor = targetFactors[ponTargetIndexes.get(factorNormalized.targets().get(i).getName())];
            final double[] inputValues = inputCounts.getRow(i);
            final double[] outputValues = factorNormalizedCounts.getRow(i);
            for (int j = 0; j < inputValues.length; j++) {
                final double expected = inputValues[j] / factor;
                Assert.assertEquals(outputValues[j], expected, 0.0000001, "" + i + " , " + j);
            }
        }
    }
}
Also used : HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN) HDF5PCACoveragePoN(org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN) HDF5File(org.broadinstitute.hdf5.HDF5File)

Aggregations

HDF5File (org.broadinstitute.hdf5.HDF5File)14 HDF5PCACoveragePoN (org.broadinstitute.hellbender.tools.pon.coverage.pca.HDF5PCACoveragePoN)14 PCACoveragePoN (org.broadinstitute.hellbender.tools.pon.coverage.pca.PCACoveragePoN)14 Array2DRowRealMatrix (org.apache.commons.math3.linear.Array2DRowRealMatrix)8 RealMatrix (org.apache.commons.math3.linear.RealMatrix)8 File (java.io.File)4 HDF5Library (org.broadinstitute.hdf5.HDF5Library)4 ArrayList (java.util.ArrayList)2 List (java.util.List)2 OptionalInt (java.util.OptionalInt)2 IntStream (java.util.stream.IntStream)2 FileUtils (org.apache.commons.io.FileUtils)2 IntRange (org.apache.commons.lang.math.IntRange)2 DefaultRealMatrixChangingVisitor (org.apache.commons.math3.linear.DefaultRealMatrixChangingVisitor)2 DefaultRealMatrixPreservingVisitor (org.apache.commons.math3.linear.DefaultRealMatrixPreservingVisitor)2 CommandLineProgramTest (org.broadinstitute.hellbender.CommandLineProgramTest)2 StandardArgumentDefinitions (org.broadinstitute.hellbender.cmdline.StandardArgumentDefinitions)2 UserException (org.broadinstitute.hellbender.exceptions.UserException)2 PoNTestUtils (org.broadinstitute.hellbender.tools.pon.PoNTestUtils)2 PCATangentNormalizationResult (org.broadinstitute.hellbender.tools.pon.coverage.pca.PCATangentNormalizationResult)2