In your Stan program, you have two transformed parameters that require symmetry:
the only transformed parameter you have that requires symmetry is this one:
transformed parameters {
cov_matrix[nt] H[T];
...
cov_matrix[nt] Qr[T];
So that means that one of the entries for H[t] has not been defined (nan is the default initialization precisely to catch these kinds of errors).
The definitions are:
Qr[1,] = Qr1_init;
H[1] = Qr[1,];
Even better would be to not have three copies of Qr1_init in the program.
This should be fine, because Qr1_init is defined to be of type cov_matrix.
So the problem is probably this:
for (t in 2:T) {
...
Qr[t,] = (1 - a_q - b_q) * S + a_q * (u[t-1,] * u[t-1,]') + b_q * Qr[t-1,];
...
H[t,] = quad_form_diag(R[t, ], D[t, ]);
I would recommend writing a little test for symmetry here and reporting an error using reject("t = ", t) if it’s not symmetric to let you know which iteration failed. I’d do it with a function (we’ve been meaning to add these to the language for years, but haven’t yet):
functions {
void test_symmetry(matrix x) {
if (rows(x) != cols(x)) reject("columns different sizes");
for (m in 1:rows(x)) {
for (n in m + 1:cols(x)) {
if (x[m, n] != x[n, m]) reject("x", {m, n}, " != x", {n, m});
}
and then insert
print('testing symmetry for H", t)`
test_symmetry(H[t]);
after each H[t] is defined and the same for the other one.
Style comment: Everywhere (left and right sides of equations) you have A[t, ] , it’s cleaner and probably more efficient to just use A[t]. If you wanted to use BUGS style or R style, that’d have to be A[t, , ] in all the cases you used it. Again, that has the same meaning as just A[t] if A has at least three indexes.