MultivariateNormalDistribution throws "matrix is singular"

Question

I want to generate a random vector with N(0, C) distribution, i.e. normal distribution with 0 mean and given covariance matrix C.

I'm using MultivariateNormalDistribution from Apache Commons:

double[] means = new double[2];
double[][] C = {
{ 5.1455479254351755, -2.0050191427617987 },
{ -2.0050191427617987, 0.7812776833676598 } };
new MultivariateNormalDistribution(new JDKRandomGenerator(), means, C);

And getting an exception matrix is singular:

Exception in thread "main" org.apache.commons.math3.linear.SingularMatrixException: matrix is singular
    at org.apache.commons.math3.linear.EigenDecomposition$Solver.getInverse(EigenDecomposition.java:533)
    at org.apache.commons.math3.distribution.MultivariateNormalDistribution.<init>(MultivariateNormalDistribution.java:125)
    at javabbob.Experiment.main(Experiment.java:52)

I've read here that this means, that the matrix is not invertible. Ok.

But, all I want is a random vector with N(0, C) distribution. I can use whatever method.

In Multivariate normal distribution Wikipedia article it's written:

The covariance matrix is allowed to be singular (in which case the corresponding distribution has no density). This case arises frequently in statistics (...)

How can I generate such a random vector in Java?

I've also tried CholeskyDecomposition with the same C array:

RealMatrix covMatrix = new Array2DRowRealMatrix(C);
CholeskyDecomposition choleskyDecomposition = new CholeskyDecomposition(covMatrix);

And it also doesn't work, throwing NonPositiveDefiniteMatrixException:

Exception in thread "main" org.apache.commons.math3.linear.NonPositiveDefiniteMatrixException: 0 is smaller than, or equal to, the minimum (0): not positive definite matrix: value 0 at index 1
    at org.apache.commons.math3.linear.CholeskyDecomposition.<init>(CholeskyDecomposition.java:142)
    at org.apache.commons.math3.linear.CholeskyDecomposition.<init>(CholeskyDecomposition.java:85)
    at javabbob.Experiment.main(Experiment.java:59)

For people with similar problem:

Not Positive Definite Matrices - Causes and Cures gave me some insight.

pjs · Answer 1 · 2013-10-28T14:35:42.873

The online matrix calculator at bluebit.gr gives the lower triangular Cholesky decomposition of your C matrix as:

L = { { 2.268380  0.000000},
      {-0.883899  0.000000} };

This means you should be able to generate a single standard normal, Z, and transform to obtain your two correlated normals as X₁ = 2.268380*Z and X₂ = -0.883899*Z. Formally it would be X₂ = -0.883899*Z₁ + 0.000000*Z₂, but since the lower right entry in the decomposition is zero you don't need a second independent Z. Add suitable means if they should be non-zero.

Addendum

Sorry, I thought you were looking for the solution to this specific problem. In general, decompose the variance/covariance matrix C using Cholesky decomposition to get the lower triangular "root" L. If M is a vector of means and Z is a vector of independent normals, then X = M + LZ will be a vector of correlated normals with the desired mean and variance/covariance structure. I'm not current in Java, so you'll have to find or implement a suitable matrix library. A quick web search showed both Apache Commons and NIST have Cholesky decomposition implementations in Java.

I think this should work, but Apache library throws `NonPositiveDefiniteMatrixException` and doesn't give L, please check my updated question. — Adam Stelmaszczyk, Nov 01 '13 at 09:48
I tried NIST and it's working! From the docs: *If the matrix is not symmetric or positive definite, the constructor returns a partial decomposition* :) Thank you! — Adam Stelmaszczyk, Nov 01 '13 at 10:20

MultivariateNormalDistribution throws "matrix is singular"

1 Answers1