Search in sources :

Example 36 with Covariance

use of org.apache.commons.math3.stat.correlation.Covariance in project narchy by automenta.

the class MyCMAESOptimizer method updateBD.

/**
 * Update B and D from C.
 *
 * @param negccov Negative covariance factor.
 */
private void updateBD(double negccov) {
    if (ccov1 + ccovmu + negccov > 0 && (iterations % 1.0 / (ccov1 + ccovmu + negccov) / dimension / dimensionDivisorWTF) < 1) {
        // to achieve O(N^2)
        C = triu(C, 0).add(triu(C, 1).transpose());
        // enforce symmetry to prevent complex numbers
        final EigenDecomposition eig = new EigenDecomposition(C);
        // eigen decomposition, B==normalized eigenvectors
        B = eig.getV();
        D = eig.getD();
        diagD = diag(D);
        if (min(diagD) <= 0) {
            for (int i = 0; i < dimension; i++) {
                if (diagD.getEntry(i, 0) < 0) {
                    diagD.setEntry(i, 0, 0);
                }
            }
            final double tfac = max(diagD) / big_magic_number_WTF;
            C = C.add(eye(dimension, dimension).scalarMultiply(tfac));
            diagD = diagD.add(ones(dimension, 1).scalarMultiply(tfac));
        }
        if (max(diagD) > big_magic_number_WTF * min(diagD)) {
            final double tfac = max(diagD) / big_magic_number_WTF - min(diagD);
            C = C.add(eye(dimension, dimension).scalarMultiply(tfac));
            diagD = diagD.add(ones(dimension, 1).scalarMultiply(tfac));
        }
        diagC = diag(C);
        // D contains standard deviations now
        diagD = sqrt(diagD);
        // O(n^2)
        BD = times(B, repmat(diagD.transpose(), dimension, 1));
    }
}
Also used : EigenDecomposition(org.apache.commons.math3.linear.EigenDecomposition)

Example 37 with Covariance

use of org.apache.commons.math3.stat.correlation.Covariance in project MindsEye by SimiaCryptus.

the class PCAUtil method pcaFeatures.

/**
 * Pca features inv tensor [ ].
 *
 * @param covariance        the covariance
 * @param components        the components
 * @param featureDimensions the feature dimensions
 * @param power             the power
 * @return the tensor [ ]
 */
public static Tensor[] pcaFeatures(final RealMatrix covariance, final int components, final int[] featureDimensions, final double power) {
    @Nonnull final EigenDecomposition decomposition = new EigenDecomposition(covariance);
    final int[] orderedVectors = IntStream.range(0, components).mapToObj(x -> x).sorted(Comparator.comparing(x -> -decomposition.getRealEigenvalue(x))).mapToInt(x -> x).toArray();
    return IntStream.range(0, orderedVectors.length).mapToObj(i -> {
        @Nonnull final Tensor src = new Tensor(decomposition.getEigenvector(orderedVectors[i]).toArray(), featureDimensions).copy();
        return src.scale(1.0 / src.rms()).scale((Math.pow(decomposition.getRealEigenvalue(orderedVectors[i]) / decomposition.getRealEigenvalue(orderedVectors[0]), power)));
    }).toArray(i -> new Tensor[i]);
}
Also used : IntStream(java.util.stream.IntStream) BlockRealMatrix(org.apache.commons.math3.linear.BlockRealMatrix) DoubleStatistics(com.simiacryptus.util.data.DoubleStatistics) Tensor(com.simiacryptus.mindseye.lang.Tensor) Supplier(java.util.function.Supplier) Collectors(java.util.stream.Collectors) RecycleBin(com.simiacryptus.mindseye.lang.RecycleBin) List(java.util.List) Stream(java.util.stream.Stream) EigenDecomposition(org.apache.commons.math3.linear.EigenDecomposition) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Comparator(java.util.Comparator) Nonnull(javax.annotation.Nonnull) EigenDecomposition(org.apache.commons.math3.linear.EigenDecomposition) Tensor(com.simiacryptus.mindseye.lang.Tensor) Nonnull(javax.annotation.Nonnull)

Example 38 with Covariance

use of org.apache.commons.math3.stat.correlation.Covariance in project MindsEye by SimiaCryptus.

the class PCAUtil method getCovariance.

/**
 * Forked from Apache Commons Math
 *
 * @param stream the stream
 * @return covariance covariance
 */
@Nonnull
public static RealMatrix getCovariance(@Nonnull final Supplier<Stream<double[]>> stream) {
    final int dimension = stream.get().findAny().get().length;
    final List<DoubleStatistics> statList = IntStream.range(0, dimension * dimension).mapToObj(i -> new DoubleStatistics()).collect(Collectors.toList());
    stream.get().forEach(array -> {
        for (int i = 0; i < dimension; i++) {
            for (int j = 0; j <= i; j++) {
                statList.get(i * dimension + j).accept(array[i] * array[j]);
            }
        }
        RecycleBin.DOUBLES.recycle(array, array.length);
    });
    @Nonnull final RealMatrix covariance = new BlockRealMatrix(dimension, dimension);
    for (int i = 0; i < dimension; i++) {
        for (int j = 0; j <= i; j++) {
            final double v = statList.get(i + dimension * j).getAverage();
            covariance.setEntry(i, j, v);
            covariance.setEntry(j, i, v);
        }
    }
    return covariance;
}
Also used : IntStream(java.util.stream.IntStream) BlockRealMatrix(org.apache.commons.math3.linear.BlockRealMatrix) DoubleStatistics(com.simiacryptus.util.data.DoubleStatistics) Tensor(com.simiacryptus.mindseye.lang.Tensor) Supplier(java.util.function.Supplier) Collectors(java.util.stream.Collectors) RecycleBin(com.simiacryptus.mindseye.lang.RecycleBin) List(java.util.List) Stream(java.util.stream.Stream) EigenDecomposition(org.apache.commons.math3.linear.EigenDecomposition) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Comparator(java.util.Comparator) Nonnull(javax.annotation.Nonnull) BlockRealMatrix(org.apache.commons.math3.linear.BlockRealMatrix) RealMatrix(org.apache.commons.math3.linear.RealMatrix) Nonnull(javax.annotation.Nonnull) DoubleStatistics(com.simiacryptus.util.data.DoubleStatistics) BlockRealMatrix(org.apache.commons.math3.linear.BlockRealMatrix) Nonnull(javax.annotation.Nonnull)

Example 39 with Covariance

use of org.apache.commons.math3.stat.correlation.Covariance in project GDSC-SMLM by aherbert.

the class MultivariateGaussianMixtureExpectationMaximizationTest method canComputeCovariance.

@SeededTest
void canComputeCovariance(RandomSeed seed) {
    final UniformRandomProvider rng = RngUtils.create(seed.getSeed());
    final int rows = 20;
    final int cols = 3;
    final double[][] data = new double[rows][cols];
    for (int i = 0; i < rows; i++) {
        for (int j = 0; j < cols; j++) {
            data[i][j] = rng.nextDouble();
        }
    }
    // Compute the same mean as the Covariance class
    final double[] means = getColumnMeans(data);
    final double[][] obs = MultivariateGaussianMixtureExpectationMaximization.covariance(means, data);
    final double[][] exp = new Covariance(data).getCovarianceMatrix().getData();
    Assertions.assertArrayEquals(exp, obs);
}
Also used : Covariance(org.apache.commons.math3.stat.correlation.Covariance) UniformRandomProvider(org.apache.commons.rng.UniformRandomProvider) SeededTest(uk.ac.sussex.gdsc.test.junit5.SeededTest)

Aggregations

RealMatrix (org.apache.commons.math3.linear.RealMatrix)27 Array2DRowRealMatrix (org.apache.commons.math3.linear.Array2DRowRealMatrix)10 IntStream (java.util.stream.IntStream)9 ArrayList (java.util.ArrayList)8 List (java.util.List)7 Nonnull (javax.annotation.Nonnull)7 Collectors (java.util.stream.Collectors)6 Stream (java.util.stream.Stream)5 Covariance (org.apache.commons.math3.stat.correlation.Covariance)5 java.util (java.util)4 Supplier (java.util.function.Supplier)4 Nullable (javax.annotation.Nullable)4 ImmutablePair (org.apache.commons.lang3.tuple.ImmutablePair)4 TooManyEvaluationsException (org.apache.commons.math3.exception.TooManyEvaluationsException)4 EigenDecomposition (org.apache.commons.math3.linear.EigenDecomposition)4 Logger (org.apache.logging.log4j.Logger)4 UserException (org.broadinstitute.hellbender.exceptions.UserException)4 Nd4jIOUtils (org.broadinstitute.hellbender.tools.coveragemodel.nd4jutils.Nd4jIOUtils)4 Utils (org.broadinstitute.hellbender.utils.Utils)4 ParamUtils (org.broadinstitute.hellbender.utils.param.ParamUtils)4