PROC IML RANDNORMAL Error

[sasuser54321]

When running the code:

proc iml;
mean=%str({0.193,0.156,1.271});
cov=%str({0.003 0.001 0.002, 0.001 0.0003 0.001, 0.002 0.001 0.008});
print cov;
print mean;
x=randnormal(1, mean, cov);
print x;
quit iml;

My log shows:

968 proc iml;
NOTE: IML Ready
969 mean=%str({0.193,0.156,1.271});
970 cov=%str({0.003 0.001 0.002, 0.001 0.0003 0.001, 0.002 0.001 0.008});
971 x=randnormal(1, mean, cov);
NOTE: Module RANDNORMAL loaded from the storage SASHELP.IMLMLIB.
NOTE: Module ROWVEC loaded from the storage SASHELP.IMLMLIB.
ERROR: The covariance matrix is not symmetric positive definite
972 print x;
ERROR: Matrix x has not been set to a value.

statement : PRINT at line 972 column 5
973 quit iml;
NOTE: Exiting IML.
NOTE: The SAS System stopped processing this step because of errors.
NOTE: PROCEDURE IML used (Total process time):
real time 0.01 seconds
cpu time 0.00 seconds

However, when running the same code with a minor modification, the code runs with success:

proc iml;
mean=%str({0.193,0.156,1.271});
cov=%str({0.003 0.001 0.002, 0.001 0.0004 0.001, 0.002 0.001 0.008});
print cov;
print mean;
x=randnormal(1, mean, cov);
print x;
quit iml;

It appears to be due to the fact that the 0.0003 in the first code is too small for randnormal because 0.0004 allows the code to run successfully.

Any thoughts on this?

[Rick_SAS · Accepted solution]

The reason that the RANDNORMAL function is complaining is that it is trying to simulate data from a multivariate normal distribution whose (population) covariance is equal to the matrix that you supply. But you are not supplying a valid covariance matrix. Every covariance matrix is positive definite, but your matrix is not.

A real symmetric matrix is "positive definite" when all of its eigenvalues are positive. If you compute the eigenvalues of this matrix when the diagonal element is 0.0003, you will see that the smallest eigenvalue is negative:

proc iml;
cov={0.003 0.001  0.002, 
     0.001 0.0003 0.001, 
     0.002 0.001  0.008};
eval = eigval(cov);
print eval;

In contrast, when you increase the size of the diagonal element to 0.0004, the eigenvalues are all positive:

cov={0.003 0.001  0.002, 
     0.001 0.0004 0.001, 
     0.002 0.001  0.008};
eval = eigval(cov);
print eval;

[sasuser54321]

Your right. I was able to resolve the issue by using the actual values in the matrix as opposed to rounded values. When using the actual values as opposed to rounded values, the covariance matrix is now valid and becomes positive definite with all eigen values of the matrix now positive.

Actual:

proc iml;
cov={.002537487 .000632960 .002179278,
.000632960 .000264712 .000680551,
.002179278 .000680551 .007732319};
eval = eigval(cov);
print eval;

Rounded:

proc iml;
cov={0.003 0.001 0.002,
0.001 0.0003 0.001,
0.002 0.001 0.008};
eval = eigval(cov);
print eval;

Thank you!