Search in sources :

Example 31 with Array2DRowRealMatrix

use of org.apache.commons.math3.linear.Array2DRowRealMatrix in project gatk-protected by broadinstitute.

the class SNPSegmenter method writeSegmentFile.

/**
     * Write segment file based on maximum-likelihood estimates of the minor allele fraction at SNP sites,
     * assuming the specified allelic bias.  These estimates are converted to target coverages,
     * which are written to a temporary file and then passed to {@link RCBSSegmenter}.
     * @param snps                  TargetCollection of allelic counts at SNP sites
     * @param sampleName            sample name
     * @param outputFile            segment file to write to and return
     * @param allelicBias           allelic bias to use in estimate of minor allele fraction
     */
public static void writeSegmentFile(final TargetCollection<AllelicCount> snps, final String sampleName, final File outputFile, final double allelicBias) {
    Utils.validateArg(snps.totalSize() > 0, "Must have a positive number of SNPs to perform SNP segmentation.");
    try {
        final File targetsFromSNPCountsFile = File.createTempFile("targets-from-snps", ".tsv");
        final List<Target> targets = snps.targets().stream().map(ac -> new Target(name(ac), ac.getInterval())).collect(Collectors.toList());
        final RealMatrix minorAlleleFractions = new Array2DRowRealMatrix(snps.targetCount(), 1);
        minorAlleleFractions.setColumn(0, snps.targets().stream().mapToDouble(ac -> ac.estimateMinorAlleleFraction(allelicBias)).toArray());
        ReadCountCollectionUtils.write(targetsFromSNPCountsFile, new ReadCountCollection(targets, Collections.singletonList(sampleName), minorAlleleFractions));
        //segment SNPs based on observed log_2 minor allele fraction (log_2 is applied in CBS.R)
        RCBSSegmenter.writeSegmentFile(sampleName, targetsFromSNPCountsFile.getAbsolutePath(), outputFile.getAbsolutePath(), false);
    } catch (final IOException e) {
        throw new UserException.CouldNotCreateOutputFile("Could not create temporary output file during " + "SNP segmentation.", e);
    }
}
Also used : Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) List(java.util.List) UserException(org.broadinstitute.hellbender.exceptions.UserException) RCBSSegmenter(org.broadinstitute.hellbender.utils.segmenter.RCBSSegmenter) AllelicCount(org.broadinstitute.hellbender.tools.exome.alleliccount.AllelicCount) Utils(org.broadinstitute.hellbender.utils.Utils) RealMatrix(org.apache.commons.math3.linear.RealMatrix) IOException(java.io.IOException) Collections(java.util.Collections) Collectors(java.util.stream.Collectors) File(java.io.File) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) IOException(java.io.IOException) UserException(org.broadinstitute.hellbender.exceptions.UserException) File(java.io.File)

Example 32 with Array2DRowRealMatrix

use of org.apache.commons.math3.linear.Array2DRowRealMatrix in project gatk-protected by broadinstitute.

the class IntegerCopyNumberTransitionProbabilityCacheUnitTest method testBasicSoundness.

@Test
public void testBasicSoundness() {
    for (final RealMatrix transitionMatrix : TRANSITION_MATRICES) {
        final IntegerCopyNumberTransitionProbabilityCache cache = new IntegerCopyNumberTransitionProbabilityCache(new IntegerCopyNumberTransitionMatrix(transitionMatrix, 0));
        for (final int dist : DISTANCES) {
            final RealMatrix transitionMatrixExponentiated = cache.getTransitionProbabilityMatrix(dist);
            /* assert positivity */
            Assert.assertTrue(Arrays.stream(transitionMatrixExponentiated.getData()).flatMapToDouble(Arrays::stream).allMatch(d -> d >= 0));
            /* assert conservation of probability */
            for (int c = 0; c < transitionMatrix.getColumnDimension(); c++) {
                Assert.assertEquals(Arrays.stream(transitionMatrixExponentiated.getColumn(c)).sum(), 1.0, EPSILON);
            }
            /* assert correctness, T(2*d) = T(d)*T(d) */
            assertEqualMatrices(cache.getTransitionProbabilityMatrix(2 * dist), transitionMatrixExponentiated.multiply(transitionMatrixExponentiated));
        }
        /* assert loss of initial state over long distances, i.e. all columns must be equal */
        final RealMatrix longRangeTransitionMatrix = cache.getTransitionProbabilityMatrix(Integer.MAX_VALUE);
        final double[] firstColumn = longRangeTransitionMatrix.getColumn(0);
        final RealMatrix syntheticLongRangeTransitionMatrix = new Array2DRowRealMatrix(firstColumn.length, firstColumn.length);
        for (int i = 0; i < firstColumn.length; i++) {
            syntheticLongRangeTransitionMatrix.setColumn(i, firstColumn);
        }
        assertEqualMatrices(longRangeTransitionMatrix, syntheticLongRangeTransitionMatrix);
        final double[] stationary = cache.getStationaryProbabilityVector().toArray();
        ArrayAsserts.assertArrayEquals(stationary, firstColumn, EPSILON);
    }
}
Also used : Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) Arrays(java.util.Arrays) Assert(org.testng.Assert) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Test(org.testng.annotations.Test) ArrayAsserts(org.testng.internal.junit.ArrayAsserts) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) Arrays(java.util.Arrays) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 33 with Array2DRowRealMatrix

use of org.apache.commons.math3.linear.Array2DRowRealMatrix in project gatk-protected by broadinstitute.

the class IntegerCopyNumberTransitionMatrixUnitTest method testPadding.

@Test
public void testPadding() {
    final IntegerCopyNumberTransitionMatrix data = new IntegerCopyNumberTransitionMatrix(new Array2DRowRealMatrix(new double[][] { { 1, 2, 3 }, { 4, 5, 6 }, { 7, 8, 9 } }), 2);
    final RealMatrix expected = new Array2DRowRealMatrix(new double[][] { { 1.0 / 12, 2.0 / 15, 3.0 / 18, 0, 0 }, { 4.0 / 12, 5.0 / 15, 6.0 / 18, 0, 0 }, { 7.0 / 12, 8.0 / 15, 9.0 / 18, 0, 0 }, { 0, 0, 0, 1, 0 }, { 0, 0, 0, 0, 1 } });
    Assert.assertEquals(data.getTransitionMatrix().subtract(expected).getNorm(), 0, 1e-12);
}
Also used : Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 34 with Array2DRowRealMatrix

use of org.apache.commons.math3.linear.Array2DRowRealMatrix in project gatk-protected by broadinstitute.

the class IntegerCopyNumberTransitionMatrixUnitTest method testUnnormalizedProbability.

@Test
public void testUnnormalizedProbability() {
    /* it should normalize unnormalized transition matrices and give a warning */
    final IntegerCopyNumberTransitionMatrix transitionMatrix = new IntegerCopyNumberTransitionMatrix(new Array2DRowRealMatrix(new double[][] { { 1, 2, 3 }, { 4, 5, 6 }, { 7, 8, 9 } }), 0);
    for (int i = 0; i < 3; i++) {
        final double[] col = transitionMatrix.getTransitionMatrix().getColumn(i);
        Assert.assertEquals(Arrays.stream(col).sum(), 1.0, 1e-12);
    }
}
Also used : Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test)

Example 35 with Array2DRowRealMatrix

use of org.apache.commons.math3.linear.Array2DRowRealMatrix in project gatk by broadinstitute.

the class SomaticGenotypingEngine method getAsRealMatrix.

//convert a likelihood matrix of alleles x reads into a RealMatrix
public static RealMatrix getAsRealMatrix(final LikelihoodMatrix<Allele> matrix) {
    final RealMatrix result = new Array2DRowRealMatrix(matrix.numberOfAlleles(), matrix.numberOfReads());
    result.walkInOptimizedOrder(new DefaultRealMatrixChangingVisitor() {

        @Override
        public double visit(int row, int column, double value) {
            return matrix.get(row, column);
        }
    });
    return result;
}
Also used : Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Array2DRowRealMatrix(org.apache.commons.math3.linear.Array2DRowRealMatrix) DefaultRealMatrixChangingVisitor(org.apache.commons.math3.linear.DefaultRealMatrixChangingVisitor)

Aggregations

Array2DRowRealMatrix (org.apache.commons.math3.linear.Array2DRowRealMatrix)141 RealMatrix (org.apache.commons.math3.linear.RealMatrix)101 Test (org.testng.annotations.Test)60 IntStream (java.util.stream.IntStream)31 BaseTest (org.broadinstitute.hellbender.utils.test.BaseTest)28 File (java.io.File)27 Collectors (java.util.stream.Collectors)25 ArrayList (java.util.ArrayList)24 Assert (org.testng.Assert)24 List (java.util.List)22 SimpleInterval (org.broadinstitute.hellbender.utils.SimpleInterval)22 Target (org.broadinstitute.hellbender.tools.exome.Target)18 java.util (java.util)15 Random (java.util.Random)14 ReadCountCollection (org.broadinstitute.hellbender.tools.exome.ReadCountCollection)14 ParamUtils (org.broadinstitute.hellbender.utils.param.ParamUtils)14 DataProvider (org.testng.annotations.DataProvider)14 Stream (java.util.stream.Stream)13 Arrays (java.util.Arrays)12 DoubleStream (java.util.stream.DoubleStream)12