org.apache.commons.math3.linear.RealMatrix Java Examples

/**
 * Update of the covariance matrix C for diagonalOnly > 0
 *
 * @param hsig Flag indicating a small correction.
 * @param bestArz Fitness-sorted matrix of the gaussian random values of the
 * current offspring.
 * @param xold xmean matrix of the previous generation.
 */
private void updateCovarianceDiagonalOnly(boolean hsig,
                                          final RealMatrix bestArz,
                                          final RealMatrix xold) {
    // minor correction if hsig==false
    double oldFac = hsig ? 0 : ccov1Sep * cc * (2. - cc);
    oldFac += 1. - ccov1Sep - ccovmuSep;
    diagC = diagC.scalarMultiply(oldFac) // regard old matrix
            // plus rank one update
            .add(square(pc).scalarMultiply(ccov1Sep))
            // plus rank mu update
            .add((times(diagC, square(bestArz).multiply(weights)))
                    .scalarMultiply(ccovmuSep));
    diagD = sqrt(diagC); // replaces eig(C)
    if (diagonalOnly > 1 && iterations > diagonalOnly) {
        // full covariance matrix from now on
        diagonalOnly = 0;
        B = eye(dimension, dimension);
        BD = diag(diagD);
        C = diag(diagC);
    }
}
 

/**
 * Test R Swiss fertility dataset against R.
 */
@Test
public void testSwissFertility() throws Exception {
     RealMatrix matrix = createRealMatrix(swissData, 47, 5);
     PearsonsCorrelation corrInstance = new PearsonsCorrelation(matrix);
     RealMatrix correlationMatrix = corrInstance.getCorrelationMatrix();
     double[] rData = new double[] {
           1.0000000000000000, 0.3530791836199747, -0.6458827064572875, -0.6637888570350691,  0.4636847006517939,
             0.3530791836199747, 1.0000000000000000,-0.6865422086171366, -0.6395225189483201, 0.4010950530487398,
            -0.6458827064572875, -0.6865422086171366, 1.0000000000000000, 0.6984152962884830, -0.5727418060641666,
            -0.6637888570350691, -0.6395225189483201, 0.6984152962884830, 1.0000000000000000, -0.1538589170909148,
             0.4636847006517939, 0.4010950530487398, -0.5727418060641666, -0.1538589170909148, 1.0000000000000000
     };
     TestUtils.assertEquals("correlation matrix", createRealMatrix(rData, 5, 5), correlationMatrix, 10E-15);

     double[] rPvalues = new double[] {
             0.01491720061472623,
             9.45043734069043e-07, 9.95151527133974e-08,
             3.658616965962355e-07, 1.304590105694471e-06, 4.811397236181847e-08,
             0.001028523190118147, 0.005204433539191644, 2.588307925380906e-05, 0.301807756132683
     };
     RealMatrix rPMatrix = createLowerTriangularRealMatrix(rPvalues, 5);
     fillUpper(rPMatrix, 0d);
     TestUtils.assertEquals("correlation p values", rPMatrix, corrInstance.getCorrelationPValues(), 10E-15);
}
 

public MultivariateTDistribution(RealVector mean, RealMatrix covarianceMatrix, double degreesOfFreedom) {
    this.mean = mean;
    if (mean.getDimension() > 1) {
        this.precisionMatrix = MatrixUtils.blockInverse(covarianceMatrix, (-1 + covarianceMatrix.getColumnDimension()) / 2);
    } else {
        this.precisionMatrix = MatrixUtils.createRealIdentityMatrix(1).scalarMultiply(1. / covarianceMatrix.getEntry(0, 0));
    }
    this.dof = degreesOfFreedom;

    this.D = mean.getDimension();

    double determinant = new LUDecomposition(covarianceMatrix).getDeterminant();

    this.multiplier = Math.exp(Gamma.logGamma(0.5 * (D + dof)) - Gamma.logGamma(0.5 * dof)) /
            Math.pow(Math.PI * dof, 0.5 * D) /
            Math.pow(determinant, 0.5);
}
 

@Test
public void testBasicStats() {

    VectorialCovariance stat = new VectorialCovariance(points[0].length, true);
    for (int i = 0; i < points.length; ++i) {
        stat.increment(points[i]);
    }

    Assert.assertEquals(points.length, stat.getN());

    RealMatrix c = stat.getResult();
    double[][] refC    = new double[][] {
            { 8.0470, -1.9195, -3.4445},
            {-1.9195,  1.0470,  3.2795},
            {-3.4445,  3.2795, 12.2070}
    };

    for (int i = 0; i < c.getRowDimension(); ++i) {
        for (int j = 0; j <= i; ++j) {
            Assert.assertEquals(refC[i][j], c.getEntry(i, j), 1.0e-12);
        }
    }

}
 

/**
 * Builds a correlated random vector generator from its mean
 * vector and covariance matrix.
 *
 * @param mean Expected mean values for all components.
 * @param covariance Covariance matrix.
 * @param small Diagonal elements threshold under which  column are
 * considered to be dependent on previous ones and are discarded
 * @param generator underlying generator for uncorrelated normalized
 * components.
 * @throws org.apache.commons.math3.linear.NonPositiveDefiniteMatrixException
 * if the covariance matrix is not strictly positive definite.
 * @throws DimensionMismatchException if the mean and covariance
 * arrays dimensions do not match.
 */
public CorrelatedRandomVectorGenerator(double[] mean,
                                       RealMatrix covariance, double small,
                                       NormalizedRandomGenerator generator) {
    int order = covariance.getRowDimension();
    if (mean.length != order) {
        throw new DimensionMismatchException(mean.length, order);
    }
    this.mean = mean.clone();

    final RectangularCholeskyDecomposition decomposition =
        new RectangularCholeskyDecomposition(covariance, small);
    root = decomposition.getRootMatrix();

    this.generator = generator;
    normalized = new double[decomposition.getRank()];

}
 

/**
 * Update of the evolution paths ps and pc.
 *
 * @param zmean Weighted row matrix of the gaussian random numbers generating
 * the current offspring.
 * @param xold xmean matrix of the previous generation.
 * @return hsig flag indicating a small correction.
 */
private boolean updateEvolutionPaths(RealMatrix zmean, RealMatrix xold) {
    ps = ps.scalarMultiply(1. - cs).add(
            B.multiply(zmean).scalarMultiply(
                    Math.sqrt(cs * (2. - cs) * mueff)));
    normps = ps.getFrobeniusNorm();
    boolean hsig = normps /
        Math.sqrt(1. - Math.pow(1. - cs, 2. * iterations)) /
            chiN < 1.4 + 2. / (dimension + 1.);
    pc = pc.scalarMultiply(1. - cc);
    if (hsig) {
        pc = pc.add(xmean.subtract(xold).scalarMultiply(
                Math.sqrt(cc * (2. - cc) * mueff) / sigma));
    }
    return hsig;
}
 

@Test
public void testTransposeTimesSelf() {
  FastByIDMap<float[]> M = new FastByIDMap<float[]>();
  M.put(1L, new float[] {4.0f, -1.0f, -5.0f});
  M.put(2L, new float[] {2.0f, 0.0f, 3.0f});
  RealMatrix MTM = MatrixUtils.transposeTimesSelf(M);
  assertArrayEquals(new double[]{20.0, -4.0, -14.0}, MTM.getRow(0));
  assertArrayEquals(new double[]{-4.0, 1.0, 5.0}, MTM.getRow(1));
  assertArrayEquals(new double[]{-14.0, 5.0, 34.0}, MTM.getRow(2));
}
 

/**
 * Constructor
 * @param matrix    the matrix to decompose
 */
private XLUD(RealMatrix matrix) {
    this.lud = new LUDecomposition(matrix);
    this.l = LazyValue.of(() -> toDataFrame(lud.getL()));
    this.u = LazyValue.of(() -> toDataFrame(lud.getU()));
    this.p = LazyValue.of(() -> toDataFrame(lud.getP()));
}
 

/**
 * @param m Input matrix.
 * @return Row matrix representing the sums of the rows.
 */
private static RealMatrix sumRows(final RealMatrix m) {
    double[][] d = new double[1][m.getColumnDimension()];
    for (int c = 0; c < m.getColumnDimension(); c++) {
        double sum = 0;
        for (int r = 0; r < m.getRowDimension(); r++) {
            sum += m.getEntry(r, c);
        }
        d[0][c] = sum;
    }
    return new Array2DRowRealMatrix(d, false);
}
 

protected void fillUpper(RealMatrix matrix, double diagonalValue) {
    int dimension = matrix.getColumnDimension();
    for (int i = 0; i < dimension; i++) {
        matrix.setEntry(i, i, diagonalValue);
        for (int j = i+1; j < dimension; j++) {
            matrix.setEntry(i, j, matrix.getEntry(j, i));
        }
    }
}
 

private void checkAEqualUSVt(RealMatrix matrix) {
    BiDiagonalTransformer transformer = new BiDiagonalTransformer(matrix);
    RealMatrix u = transformer.getU();
    RealMatrix b = transformer.getB();
    RealMatrix v = transformer.getV();
    double norm = u.multiply(b).multiply(v.transpose()).subtract(matrix).getNorm();
    Assert.assertEquals(0, norm, 1.0e-14);
}
 

/**
 * @param start Start value.
 * @param end End value.
 * @param step Step size.
 * @return a sequence as column matrix.
 */
private static RealMatrix sequence(double start, double end, double step) {
    int size = (int) ((end - start) / step + 1);
    double[][] d = new double[size][1];
    double value = start;
    for (int r = 0; r < size; r++) {
        d[r][0] = value;
        value += step;
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * @param m Input matrix.
 * @return Row matrix representing the sums of the rows.
 */
private static RealMatrix sumRows(final RealMatrix m) {
    double[][] d = new double[1][m.getColumnDimension()];
    for (int c = 0; c < m.getColumnDimension(); c++) {
        double sum = 0;
        for (int r = 0; r < m.getRowDimension(); r++) {
            sum += m.getEntry(r, c);
        }
        d[0][c] = sum;
    }
    return new Array2DRowRealMatrix(d, false);
}
 

protected RealMatrix createRealMatrix(double[] data, int nRows, int nCols) {
    double[][] matrixData = new double[nRows][nCols];
    int ptr = 0;
    for (int i = 0; i < nRows; i++) {
        System.arraycopy(data, ptr, matrixData[i], 0, nCols);
        ptr += nCols;
    }
    return new Array2DRowRealMatrix(matrixData);
}
 

/**
 * @param size Number of rows.
 * @param popSize Population size.
 * @return a 2-dimensional matrix of Gaussian random numbers.
 */
private RealMatrix randn1(int size, int popSize) {
    double[][] d = new double[size][popSize];
    for (int r = 0; r < size; r++) {
        for (int c = 0; c < popSize; c++) {
            d[r][c] = random.nextGaussian();
        }
    }
    return new Array2DRowRealMatrix(d, false);
}
 

@Test
public void testLoad() {
  RealMatrix symmetric = randomSymmetricMatrix(500);
  Stopwatch stopwatch = new Stopwatch().start();
  int iterations = 100;
  for (int i = 0; i < iterations; i++) {
    MatrixUtils.getSolver(symmetric);
  }
  stopwatch.stop();
  long elapsedMS = stopwatch.elapsed(TimeUnit.MILLISECONDS);
  log.info("{}ms elapsed", elapsedMS);
  assertTrue(elapsedMS < 300 * iterations);
}
 

/**
 * Create a PearsonsCorrelation from a {@link Covariance}.  The correlation
 * matrix is computed by scaling the Covariance's covariance matrix.
 * The Covariance instance must have been created from a data matrix with
 * columns representing variable values.
 *
 * @param covariance Covariance instance
 */
public PearsonsCorrelation(Covariance covariance) {
    RealMatrix covarianceMatrix = covariance.getCovarianceMatrix();
    if (covarianceMatrix == null) {
        throw new NullArgumentException(LocalizedFormats.COVARIANCE_MATRIX);
    }
    nObs = covariance.getN();
    correlationMatrix = covarianceToCorrelation(covarianceMatrix);
}
 

/**
 * @param m Input matrix
 * @return Matrix representing the element-wise square (^2) of m.
 */
private static RealMatrix square(final RealMatrix m) {
    double[][] d = new double[m.getRowDimension()][m.getColumnDimension()];
    for (int r = 0; r < m.getRowDimension(); r++) {
        for (int c = 0; c < m.getColumnDimension(); c++) {
            double e = m.getEntry(r, c);
            d[r][c] = e * e;
        }
    }
    return new Array2DRowRealMatrix(d, false);
}
 

@Test(dataProvider = "tooManyZerosData")
public void testRemoveColumnsWithTooManyZeros(final ReadCountCollection readCount) {
    final RealMatrix counts = readCount.counts();
    final int[] numberOfZeros = IntStream.range(0, counts.getColumnDimension())
            .map(i -> (int) DoubleStream.of(counts.getColumn(i)).filter(d -> d == 0.0).count()).toArray();
    final int maximumNumberOfZeros = IntStream.of(numberOfZeros).max().getAsInt();
    for (int maxZeros = 0; maxZeros < maximumNumberOfZeros; maxZeros++) {
        final int maxZerosThres = maxZeros;
        final int expectedRemainingCount = (int) IntStream.of(numberOfZeros).filter(i -> i <= maxZerosThres).count();
        if (expectedRemainingCount == 0) {
            try {
                ReadCountCollectionUtils.removeColumnsWithTooManyZeros(readCount, maxZeros, false, NULL_LOGGER);
            } catch (final UserException.BadInput ex) {
                // expected.
                continue;
            }
            Assert.fail("expects an exception");
        }
        final ReadCountCollection rc = ReadCountCollectionUtils.removeColumnsWithTooManyZeros(readCount, maxZeros, false, NULL_LOGGER);
        Assert.assertEquals(rc.columnNames().size(), expectedRemainingCount);
        final int[] newIndices = new int[expectedRemainingCount];
        int nextIndex = 0;
        for (int i = 0; i < readCount.columnNames().size(); i++) {
            final String name = readCount.columnNames().get(i);
            final int newIndex = rc.columnNames().indexOf(name);
            if (numberOfZeros[i] <= maxZeros) {
                Assert.assertTrue(newIndex >= 0);
                newIndices[nextIndex++] = i;
            } else {
                Assert.assertEquals(newIndex, -1);
            }
        }
        Assert.assertEquals(nextIndex, expectedRemainingCount);
        for (int i = 1; i < newIndices.length; i++) {
            Assert.assertTrue(newIndices[i - 1] < newIndices[i]);
        }
    }
}
 

@Test
public void testSimplistic() {
    VectorialCovariance stat = new VectorialCovariance(2, true);
    stat.increment(new double[] {-1.0,  1.0});
    stat.increment(new double[] { 1.0, -1.0});
    RealMatrix c = stat.getResult();
    Assert.assertEquals( 2.0, c.getEntry(0, 0), 1.0e-12);
    Assert.assertEquals(-2.0, c.getEntry(1, 0), 1.0e-12);
    Assert.assertEquals( 2.0, c.getEntry(1, 1), 1.0e-12);
}
 

/**
 * @param m Input matrix.
 * @return Row matrix representing the sums of the rows.
 */
private static RealMatrix sumRows(final RealMatrix m) {
    double[][] d = new double[1][m.getColumnDimension()];
    for (int c = 0; c < m.getColumnDimension(); c++) {
        double sum = 0;
        for (int r = 0; r < m.getRowDimension(); r++) {
            sum += m.getEntry(r, c);
        }
        d[0][c] = sum;
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * @param m Input matrix.
 * @param k Diagonal position.
 * @return Upper triangular part of matrix.
 */
private static RealMatrix triu(final RealMatrix m, int k) {
    double[][] d = new double[m.getRowDimension()][m.getColumnDimension()];
    for (int r = 0; r < m.getRowDimension(); r++) {
        for (int c = 0; c < m.getColumnDimension(); c++) {
            d[r][c] = r <= c - k ? m.getEntry(r, c) : 0;
        }
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * @param n Number of rows.
 * @param m Number of columns.
 * @return n-by-m matrix filled with 1.
 */
private static RealMatrix ones(int n, int m) {
    double[][] d = new double[n][m];
    for (int r = 0; r < n; r++) {
        Arrays.fill(d[r], 1.0);
    }
    return new Array2DRowRealMatrix(d, false);
}
 

private void assertOutputContents(final File output, final List<Target> targets, final List<String> sampleNames, final double[][] counts) throws IOException {
    Assert.assertTrue(output.length() > 10);
    final ReadCountCollection readCounts = ReadCountCollectionUtils.parse(output);
    final List<String> sortedSampleNames = new ArrayList<>(sampleNames);
    Collections.sort(sortedSampleNames);
    Assert.assertEquals(readCounts.targets(), targets);
    Assert.assertEquals(readCounts.columnNames(), sortedSampleNames);
    final RealMatrix actualCounts = readCounts.counts();
    for (int i = 0; i < actualCounts.getColumnDimension(); i++) {
        final int actualColumn = readCounts.columnNames().indexOf(sampleNames.get(i));
        for (int j = 0; j < actualCounts.getRowDimension(); j++) {
            Assert.assertEquals(actualCounts.getEntry(j, actualColumn), counts[i][j], 0.00001);
        }
    }
}
 

/**
 * @param size Number of rows.
 * @param popSize Population size.
 * @return a 2-dimensional matrix of Gaussian random numbers.
 */
private RealMatrix randn1(int size, int popSize) {
    double[][] d = new double[size][popSize];
    for (int r = 0; r < size; r++) {
        for (int c = 0; c < popSize; c++) {
            d[r][c] = random.nextGaussian();
        }
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * @param mat Input matrix.
 * @param n Number of row replicates.
 * @param m Number of column replicates.
 * @return a matrix which replicates the input matrix in both directions.
 */
private static RealMatrix repmat(final RealMatrix mat, int n, int m) {
    int rd = mat.getRowDimension();
    int cd = mat.getColumnDimension();
    double[][] d = new double[n * rd][m * cd];
    for (int r = 0; r < n * rd; r++) {
        for (int c = 0; c < m * cd; c++) {
            d[r][c] = mat.getEntry(r % rd, c % cd);
        }
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * @param m Input matrix.
 * @param cols Columns to select.
 * @return Matrix representing the selected columns.
 */
private static RealMatrix selectColumns(final RealMatrix m, final int[] cols) {
    double[][] d = new double[m.getRowDimension()][cols.length];
    for (int r = 0; r < m.getRowDimension(); r++) {
        for (int c = 0; c < cols.length; c++) {
            d[r][c] = m.getEntry(r, cols[c]);
        }
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * Apply the pipeline to input data
 * @param data
 * @return
 */
protected final RealMatrix pipelineFitTransform(RealMatrix data) {
	RealMatrix operated = data;
	
	// Push through pipeline... fits the models in place
	for(PreProcessor pre: pipe)
		operated = pre.fit(operated).transform(operated);
	
	return operated;
}
 

/**
 * @param m Input matrix.
 * @return Row matrix representing the sums of the rows.
 */
private static RealMatrix sumRows(final RealMatrix m) {
    double[][] d = new double[1][m.getColumnDimension()];
    for (int c = 0; c < m.getColumnDimension(); c++) {
        double sum = 0;
        for (int r = 0; r < m.getRowDimension(); r++) {
            sum += m.getEntry(r, c);
        }
        d[0][c] = sum;
    }
    return new Array2DRowRealMatrix(d, false);
}
 

/**
 * @param m Input matrix 1.
 * @param n Input matrix 2.
 * @return the matrix where the elements of m and n are element-wise multiplied.
 */
private static RealMatrix times(final RealMatrix m, final RealMatrix n) {
    double[][] d = new double[m.getRowDimension()][m.getColumnDimension()];
    for (int r = 0; r < m.getRowDimension(); r++) {
        for (int c = 0; c < m.getColumnDimension(); c++) {
            d[r][c] = m.getEntry(r, c) * n.getEntry(r, c);
        }
    }
    return new Array2DRowRealMatrix(d, false);
}