Examples with SingularMatrixException org.apache.commons.math3.linear.SingularMatrixException used on opensource projects

Example 6 with SingularMatrixException

use of org.apache.commons.math3.linear.SingularMatrixException in project knime-core by knime.

the class IrlsLearner method learn.

/**
 * {@inheritDoc}
 */
@Override
public LogRegLearnerResult learn(final TrainingData<ClassificationTrainingRow> trainingData, final ExecutionMonitor exec) throws CanceledExecutionException, InvalidSettingsException {
    exec.checkCanceled();
    int iter = 0;
    boolean converged = false;
    final int tcC = trainingData.getTargetDimension() + 1;
    final int rC = trainingData.getFeatureCount() - 1;
    final RealMatrix beta = MatrixUtils.createRealMatrix(1, (tcC - 1) * (rC + 1));
    Double loglike = 0.0;
    Double loglikeOld = 0.0;
    exec.setMessage("Iterative optimization. Processing iteration 1.");
    // main loop
    while (iter < m_maxIter && !converged) {
        RealMatrix betaOld = beta.copy();
        loglikeOld = loglike;
        // Do heavy work in a separate thread which allows to interrupt it
        // note the queue may block if no more threads are available (e.g. thread count = 1)
        // as soon as we stall in 'get' this thread reduces the number of running thread
        Future<Double> future = ThreadPool.currentPool().enqueue(new Callable<Double>() {

            @Override
            public Double call() throws Exception {
                final ExecutionMonitor progMon = exec.createSubProgress(1.0 / m_maxIter);
                irlsRls(trainingData, beta, rC, tcC, progMon);
                progMon.setProgress(1.0);
                return likelihood(trainingData.iterator(), beta, rC, tcC, exec);
            }
        });
        try {
            loglike = future.get();
        } catch (InterruptedException e) {
            future.cancel(true);
            exec.checkCanceled();
            throw new RuntimeException(e);
        } catch (ExecutionException e) {
            if (e.getCause() instanceof RuntimeException) {
                throw (RuntimeException) e.getCause();
            } else {
                throw new RuntimeException(e.getCause());
            }
        }
        if (Double.isInfinite(loglike) || Double.isNaN(loglike)) {
            throw new RuntimeException(FAILING_MSG);
        }
        exec.checkCanceled();
        // test for decreasing likelihood
        while ((Double.isInfinite(loglike) || Double.isNaN(loglike) || loglike < loglikeOld) && iter > 0) {
            converged = true;
            for (int k = 0; k < beta.getColumnDimension(); k++) {
                if (abs(beta.getEntry(0, k) - betaOld.getEntry(0, k)) > m_eps * abs(betaOld.getEntry(0, k))) {
                    converged = false;
                    break;
                }
            }
            if (converged) {
                break;
            }
            // half the step size of beta
            beta.setSubMatrix((beta.add(betaOld)).scalarMultiply(0.5).getData(), 0, 0);
            exec.checkCanceled();
            loglike = likelihood(trainingData.iterator(), beta, rC, tcC, exec);
            exec.checkCanceled();
        }
        // test for convergence
        converged = true;
        for (int k = 0; k < beta.getColumnDimension(); k++) {
            if (abs(beta.getEntry(0, k) - betaOld.getEntry(0, k)) > m_eps * abs(betaOld.getEntry(0, k))) {
                converged = false;
                break;
            }
        }
        iter++;
        LOGGER.debug("#Iterations: " + iter);
        LOGGER.debug("Log Likelihood: " + loglike);
        StringBuilder betaBuilder = new StringBuilder();
        for (int i = 0; i < beta.getColumnDimension() - 1; i++) {
            betaBuilder.append(Double.toString(beta.getEntry(0, i)));
            betaBuilder.append(", ");
        }
        if (beta.getColumnDimension() > 0) {
            betaBuilder.append(Double.toString(beta.getEntry(0, beta.getColumnDimension() - 1)));
        }
        LOGGER.debug("beta: " + betaBuilder.toString());
        exec.checkCanceled();
        exec.setMessage("Iterative optimization. #Iterations: " + iter + " | Log-likelihood: " + DoubleFormat.formatDouble(loglike) + ". Processing iteration " + (iter + 1) + ".");
    }
    StringBuilder warnBuilder = new StringBuilder();
    if (iter >= m_maxIter) {
        warnBuilder.append("The algorithm did not reach convergence after the specified number of epochs. " + "Setting the epoch limit higher might result in a better model.");
    }
    // The covariance matrix
    RealMatrix covMat = null;
    if (m_calcCovMatrix) {
        try {
            covMat = new QRDecomposition(A).getSolver().getInverse().scalarMultiply(-1);
        } catch (SingularMatrixException sme) {
            if (warnBuilder.length() > 0) {
                warnBuilder.append("\n");
            }
            warnBuilder.append("The covariance matrix could not be calculated because the" + " observed fisher information matrix was singular.");
        }
    }
    RealMatrix betaMat = MatrixUtils.createRealMatrix(tcC - 1, rC + 1);
    for (int i = 0; i < beta.getColumnDimension(); i++) {
        int r = i / (rC + 1);
        int c = i % (rC + 1);
        betaMat.setEntry(r, c, beta.getEntry(0, i));
    }
    m_warning = warnBuilder.length() > 0 ? warnBuilder.toString() : null;
    return new LogRegLearnerResult(betaMat, covMat, iter, loglike);
}

Example 7 with SingularMatrixException

use of org.apache.commons.math3.linear.SingularMatrixException in project imagingbook-common by imagingbook.

the class Matrix method inverse.

// Matrix inversion ---------------------------------------
/**
 * @param A a square matrix.
 * @return the inverse of A or null if A is non-square or singular.
 */
public static double[][] inverse(final double[][] A) {
    RealMatrix M = MatrixUtils.createRealMatrix(A);
    if (!M.isSquare())
        return null;
    else {
        double[][] Ai = null;
        try {
            // new LUDecomposition(M).getSolver().getInverse();
            RealMatrix Mi = MatrixUtils.inverse(M);
            Ai = Mi.getData();
        } catch (SingularMatrixException e) {
        }
        return Ai;
    }
}

Example 8 with SingularMatrixException

use of org.apache.commons.math3.linear.SingularMatrixException in project imagingbook-common by imagingbook.

the class Matrix method solve.

// use general method, i.e. double[][] inverse(double[][] A)
// @Deprecated
// public static double[][] inverse3x3(final double[][] A) {
// double[][] B = duplicate(A);
// final double det = determinant3x3(B);
// if (Math.abs(det) < Arithmetic.EPSILON_DOUBLE)
// return null;
// else {
// final double a00 = B[0][0];
// final double a01 = B[0][1];
// final double a02 = B[0][2];
// final double a10 = B[1][0];
// final double a11 = B[1][1];
// final double a12 = B[1][2];
// final double a20 = B[2][0];
// final double a21 = B[2][1];
// final double a22 = B[2][2];
// B[0][0] = (a11 * a22 - a12 * a21) / det;
// B[0][1] = (a02 * a21 - a01 * a22) / det;
// B[0][2] = (a01 * a12 - a02 * a11) / det;
// 
// B[1][0] = (a12 * a20 - a10 * a22) / det;
// B[1][1] = (a00 * a22 - a02 * a20) / det;
// B[1][2] = (a02 * a10 - a00 * a12) / det;
// 
// B[2][0] = (a10 * a21 - a11 * a20) / det;
// B[2][1] = (a01 * a20 - a00 * a21) / det;
// B[2][2] = (a00 * a11 - a01 * a10) / det;
// return B;
// }
// }
// numerically stable? should be replaced by standard inversion
// @Deprecated
// public static float[][] inverse3x3(final float[][] A) {
// final float[][] B = duplicate(A);
// final double det = determinant3x3(B);
// // IJ.log("   determinant = " + det);
// if (Math.abs(det) < Arithmetic.EPSILON_DOUBLE)
// return null;
// else {
// final double a00 = B[0][0];
// final double a01 = B[0][1];
// final double a02 = B[0][2];
// final double a10 = B[1][0];
// final double a11 = B[1][1];
// final double a12 = B[1][2];
// final double a20 = B[2][0];
// final double a21 = B[2][1];
// final double a22 = B[2][2];
// B[0][0] = (float) ((a11 * a22 - a12 * a21) / det);
// B[0][1] = (float) ((a02 * a21 - a01 * a22) / det);
// B[0][2] = (float) ((a01 * a12 - a02 * a11) / det);
// 
// B[1][0] = (float) ((a12 * a20 - a10 * a22) / det);
// B[1][1] = (float) ((a00 * a22 - a02 * a20) / det);
// B[1][2] = (float) ((a02 * a10 - a00 * a12) / det);
// 
// B[2][0] = (float) ((a10 * a21 - a11 * a20) / det);
// B[2][1] = (float) ((a01 * a20 - a00 * a21) / det);
// B[2][2] = (float) ((a00 * a11 - a01 * a10) / det);
// return B;
// }
// }
// ------------------------------------------------------------------------
// Finds the EXACT solution x for A.x = b
public static double[] solve(final double[][] A, double[] b) {
    RealMatrix AA = MatrixUtils.createRealMatrix(A);
    RealVector bb = MatrixUtils.createRealVector(b);
    DecompositionSolver solver = new LUDecomposition(AA).getSolver();
    double[] x = null;
    try {
        x = solver.solve(bb).toArray();
    } catch (SingularMatrixException e) {
    }
    return x;
}

Example 9 with SingularMatrixException

use of org.apache.commons.math3.linear.SingularMatrixException in project FSensor by KalebKE.

the class RotationKalmanFilter method correct.

/**
 * Correct the current state estimate with an actual measurement.
 *
 * @param z
 *            the measurement vector
 * @throws NullArgumentException
 *             if the measurement vector is {@code null}
 * @throws DimensionMismatchException
 *             if the dimension of the measurement vector does not fit
 * @throws SingularMatrixException
 *             if the covariance matrix could not be inverted
 */
public void correct(final RealVector z) throws NullArgumentException, DimensionMismatchException, SingularMatrixException {
    // sanity checks
    MathUtils.checkNotNull(z);
    if (z.getDimension() != measurementMatrix.getRowDimension()) {
        throw new DimensionMismatchException(z.getDimension(), measurementMatrix.getRowDimension());
    }
    // S = H * P(k) * H' + R
    RealMatrix s = measurementMatrix.multiply(errorCovariance).multiply(measurementMatrixT).add(measurementModel.getMeasurementNoise());
    // Inn = z(k) - H * xHat(k)-
    RealVector innovation = z.subtract(measurementMatrix.operate(stateEstimation));
    // calculate gain matrix
    // K(k) = P(k)- * H' * (H * P(k)- * H' + R)^-1
    // K(k) = P(k)- * H' * S^-1
    // instead of calculating the inverse of S we can rearrange the formula,
    // and then solve the linear equation A x X = B with A = S', X = K' and
    // B = (H * P)'
    // K(k) * S = P(k)- * H'
    // S' * K(k)' = H * P(k)-'
    RealMatrix kalmanGain = new CholeskyDecomposition(s).getSolver().solve(measurementMatrix.multiply(errorCovariance.transpose())).transpose();
    // update estimate with measurement z(k)
    // xHat(k) = xHat(k)- + K * Inn
    stateEstimation = stateEstimation.add(kalmanGain.operate(innovation));
    // update covariance of prediction error
    // P(k) = (I - K * H) * P(k)-
    RealMatrix identity = MatrixUtils.createRealIdentityMatrix(kalmanGain.getRowDimension());
    errorCovariance = identity.subtract(kalmanGain.multiply(measurementMatrix)).multiply(errorCovariance);
}

Example 10 with SingularMatrixException

use of org.apache.commons.math3.linear.SingularMatrixException in project tetrad by cmu-phil.

the class SemBicScore method localScoreDiff.

@Override
public double localScoreDiff(int x, int y, int[] z) {
    Node _x = variables.get(x);
    Node _y = variables.get(y);
    List<Node> _z = getVariableList(z);
    double r;
    try {
        r = partialCorrelation(_x, _y, _z);
    } catch (SingularMatrixException e) {
        // System.out.println(SearchLogUtils.determinismDetected(_z, _x));
        return Double.NaN;
    }
    int p = 2 + z.length;
    int N = covariances.getSampleSize();
    return -N * Math.log(1.0 - r * r) - p * getPenaltyDiscount() * Math.log(N);
// return localScore(y, append(z, x)) - localScore(y, z);
}