@Override public SlopeCoefficients estimateCoefficients(final DerivationEquation eq) throws EstimationException { final double[][] sourceTriangleMatrix = eq.getCovarianceLowerTriangularMatrix(); // Copy matrix and enhance it to a full matrix as expected by CholeskyDecomposition // FIXME: Avoid copy job to speed-up the solving process e.g. by extending the CholeskyDecomposition constructor final int length = sourceTriangleMatrix.length; final double[][] matrix = new double[length][]; for (int i = 0; i < length; i++) { matrix[i] = new double[length]; final double[] s = sourceTriangleMatrix[i]; final double[] t = matrix[i]; for (int j = 0; j <= i; j++) { t[j] = s[j]; } for (int j = i + 1; j < length; j++) { t[j] = sourceTriangleMatrix[j][i]; } } final RealMatrix coefficients = new Array2DRowRealMatrix(matrix, false); try { final DecompositionSolver solver = new CholeskyDecomposition(coefficients).getSolver(); final RealVector constants = new ArrayRealVector(eq.getConstraints(), true); final RealVector solution = solver.solve(constants); return new DefaultSlopeCoefficients(solution.toArray()); } catch (final NonPositiveDefiniteMatrixException e) { throw new EstimationException("Matrix inversion error due to data is linearly dependent", e); } }
@Override protected RealVector solve(final RealMatrix jacobian, final RealVector residuals) { try { final Pair<RealMatrix, RealVector> normalEquation = computeNormalMatrix(jacobian, residuals); final RealMatrix normal = normalEquation.getFirst(); final RealVector jTr = normalEquation.getSecond(); return new CholeskyDecomposition( normal, SINGULARITY_THRESHOLD, SINGULARITY_THRESHOLD) .getSolver() .solve(jTr); } catch (NonPositiveDefiniteMatrixException e) { throw new ConvergenceException(LocalizedFormats.UNABLE_TO_SOLVE_SINGULAR_PROBLEM, e); } }
/** * Creates a multivariate normal distribution with the given mean vector and * covariance matrix. * <br/> * The number of dimensions is equal to the length of the mean vector * and to the number of rows and columns of the covariance matrix. * It is frequently written as "p" in formulae. * * @param rng Random Number Generator. * @param means Vector of means. * @param covariances Covariance matrix. * @throws DimensionMismatchException if the arrays length are * inconsistent. * @throws SingularMatrixException if the eigenvalue decomposition cannot * be performed on the provided covariance matrix. * @throws NonPositiveDefiniteMatrixException if any of the eigenvalues is * negative. */ public MultivariateNormalDistribution(RandomGenerator rng, final double[] means, final double[][] covariances) throws SingularMatrixException, DimensionMismatchException, NonPositiveDefiniteMatrixException { super(rng, means.length); final int dim = means.length; if (covariances.length != dim) { throw new DimensionMismatchException(covariances.length, dim); } for (int i = 0; i < dim; i++) { if (dim != covariances[i].length) { throw new DimensionMismatchException(covariances[i].length, dim); } } this.means = MathArrays.copyOf(means); covarianceMatrix = new Array2DRowRealMatrix(covariances); // Covariance matrix eigen decomposition. final EigenDecomposition covMatDec = new EigenDecomposition(covarianceMatrix); // Compute and store the inverse. covarianceMatrixInverse = covMatDec.getSolver().getInverse(); // Compute and store the determinant. covarianceMatrixDeterminant = covMatDec.getDeterminant(); // Eigenvalues of the covariance matrix. final double[] covMatEigenvalues = covMatDec.getRealEigenvalues(); for (int i = 0; i < covMatEigenvalues.length; i++) { if (covMatEigenvalues[i] < 0) { throw new NonPositiveDefiniteMatrixException(covMatEigenvalues[i], i, 0); } } // Matrix where each column is an eigenvector of the covariance matrix. final Array2DRowRealMatrix covMatEigenvectors = new Array2DRowRealMatrix(dim, dim); for (int v = 0; v < dim; v++) { final double[] evec = covMatDec.getEigenvector(v).toArray(); covMatEigenvectors.setColumn(v, evec); } final RealMatrix tmpMatrix = covMatEigenvectors.transpose(); // Scale each eigenvector by the square root of its eigenvalue. for (int row = 0; row < dim; row++) { final double factor = FastMath.sqrt(covMatEigenvalues[row]); for (int col = 0; col < dim; col++) { tmpMatrix.multiplyEntry(row, col, factor); } } samplingMatrix = covMatEigenvectors.multiply(tmpMatrix); }
/** * Creates a multivariate normal distribution with the given mean vector and * covariance matrix. * <br/> * The number of dimensions is equal to the length of the mean vector * and to the number of rows and columns of the covariance matrix. * It is frequently written as "p" in formulae. * <p> * <b>Note:</b> this constructor will implicitly create an instance of * {@link Well19937c} as random generator to be used for sampling only (see * {@link #sample()} and {@link #sample(int)}). In case no sampling is * needed for the created distribution, it is advised to pass {@code null} * as random generator via the appropriate constructors to avoid the * additional initialisation overhead. * * @param means Vector of means. * @param covariances Covariance matrix. * @throws DimensionMismatchException if the arrays length are * inconsistent. * @throws SingularMatrixException if the eigenvalue decomposition cannot * be performed on the provided covariance matrix. * @throws NonPositiveDefiniteMatrixException if any of the eigenvalues is * negative. */ public MultivariateNormalDistribution(final double[] means, final double[][] covariances) throws SingularMatrixException, DimensionMismatchException, NonPositiveDefiniteMatrixException { this(new Well19937c(), means, covariances); }
/** * Creates a multivariate normal distribution with the given mean vector and * covariance matrix. * <br/> * The number of dimensions is equal to the length of the mean vector * and to the number of rows and columns of the covariance matrix. * It is frequently written as "p" in formulae. * * @param means Vector of means. * @param covariances Covariance matrix. * @throws DimensionMismatchException if the arrays length are * inconsistent. * @throws SingularMatrixException if the eigenvalue decomposition cannot * be performed on the provided covariance matrix. * @throws NonPositiveDefiniteMatrixException if any of the eigenvalues is * negative. */ public MultivariateNormalDistribution(final double[] means, final double[][] covariances) throws SingularMatrixException, DimensionMismatchException, NonPositiveDefiniteMatrixException { this(new Well19937c(), means, covariances); }
