inv_wishart resulting in asymmetric matrices

[ecod3r]

I am trying to use an inverse Wishart prior on a covariance matrix, however that is failing miserably no matter how easy the model I try. I have had problems with something as basic as

data {
    int<lower=0> K;
    matrix[K, K] W;
    vector[K] y;
    vector[K] mu;
}
parameters {
    matrix[K, K] Sigma;
}
model {
    Sigma ~ inv_wishart(2.0, W);
    y ~ multi_normal(mu, Sigma);
}

Where K = 2, W = [[1, 0], [0,1]], mu = [0, 0], y = [1, -0.5]. Note that I have fixed nu to 2.0. The errors I get are either look like the following:

Error in function stan::prob::inv_wishart_log(d): LDLT_Factor of random variableunderlying matrix is not positive definite. LDLT_Factor of random variablelast conditional variance is -1.5017314906815256.Rejecting proposed initial value with zero density.

and

Error in function stan::prob::multi_normal_log(N4stan5agrad3varE): Covariance matrixCovariance matrix is not symmetric. Covariance matrix[0,1] is -1.2426441399694403:0, but Covariance matrix[1,0] element is -0.488805:0Rejecting proposed initial value with zero density.

Eventually, Stan gives up printing

Initialization failed after 100 attempts. Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.

[randommm]

Hi, try instead:

parameters {
    cov_matrix[K] Sigma;
}

[ecod3r]

Thank you very much, that solved the problem!

It may save you (and some other clueless users such as myself) some time if it were stressed the need to use that type either in a note where the function is listed, or perhaps via a warning when the model is compiled.

[bob-carpenter]

On May 22, 2014, at 6:45 AM, ecoskian notifications@github.com wrote:

Thank you very much, that solved the problem!

It may save you (and some other clueless users such as myself) some time if it were stressed the need to use that type either in a note where the function is listed, or perhaps via a warning when the model is compiled.

That would be nice, but it's not so easy because there's no way for us to know at compile time if a user's constructed a proper symmetric, positive-definite matrix.

In general, we're trying to more heavily stress that in general the constraints on a parameter need to match its support. That's why if you have

Sigma ~ inv_wishart(...),

you need to ensure Sigma is symmetric, and positive definite, which can be achieved by declaring Sigma to be a cov_matrix parameter, or by constructing it as such. For instance, you can construct a 2 x 2 covariance matrix a out of a correlation and two scales:

data { cov_matrix[2] Omega; }

parameters { vector[2] sigma; real rho; }

model { matrix[2,2] a;

a[1,1] <- square(sigma1); a[2,2] <- square(sigma2); a[1,2] <- sigma1 * sigma2 * rho; a[2,1] <- a[1,2];

a ~ inv_wishart(Omega);

}

and then a is OK to be on the left-hand side of an inv_wishart sample.

Note that the scale and correlation parameters themselves need to be constrained here.

• Bob

— Reply to this email directly or view it on GitHub.