# Covariance matrix error because of float precision

Hi all.

I am trying to implement the BEKK model(https://arxiv.org/pdf/math/0702815.pdf) using Stan and I face the following problem. When adding terms like A * A^T it returns some exceptions (probably because of float precision). I would appreciate any help to avoid this problem. The following is that part of my code and the exception message.

Best,
Edgar

``````parameters {
matrix[N,N] A;
matrix[N,N] B;
cov_matrix[N] S0;
}
transformed parameters {
cov_matrix[N] S[T];t
for(t in 1:T)
{
if (t == 1){
S[t] = S0;
}
else{
S[t] = S0 + A * y[t-1] * y[t-1]' * A' + B * S[t-1] * B';  // the line that
// throws the
// exception
}
}
}
``````

You would to use functions like `crossprod`, `tcrossprod`, and `quad_form` to enforce symmetry. Something like

``````S[t] = S0 + tcrossprod(A * y[t-1]) + quad_form(S[t - 1], B');
``````
1 Like

I’ve generally had to resort to passing my matrices through something like:

``````  matrix makesym(matrix mat, int verbose, int pd){
matrix[rows(mat),cols(mat)] out;
for(coli in 1:cols(mat)){
if(pd ==1 && mat[coli,coli] < 1e-5){
if(verbose > 0) print("diagonal too low (",mat[coli,coli],") during makesym row ", coli, " col ", coli);
out[coli,coli] = 1e-5;
} else out[coli,coli] = mat[coli,coli];
for(rowi in coli:rows(mat)){
if(rowi > coli) {
out[rowi,coli] = mat[rowi,coli];
out[coli,rowi] = mat[rowi,coli];
}
}
}
return out;
}``````
1 Like