Cholesky_decompose: m is not symmetric. m[1,2] = nan, but m[2,1] = nan

Modeling
,
1

I want to use state space smoother in my model, but I have got the error message. my model is

data {
  int<lower=1> T; //the number of sample
  int<lower=1> K; //number of K years
  
  matrix[11,T] y; //data2
 
  real h;
  vector[K] tau;  // represent K different year 
}


parameters {
  //parameters for mu
  
  vector<lower=0>[5] kappa;
  vector[4] kappaq;
  vector[5] mu;
  vector<lower=0>[5] sigma;
  vector<lower=0>[11] eps;
  vector[5] aa[T];
  real cr;
  real ca;
  real cg;
  real gr;
 
}
transformed parameters{
  
  //matrix[5,5] sigr=diag_matrix(sigma);
  matrix[5,5] phi;
  matrix[5,5] Q;
 
  vector[5] Mu;
  matrix[5,5] p0;
  vector[5] a1;
  matrix[5,5] p1;
  matrix[5,5] P[T];
  vector[5] a[T];
  matrix[5,5] sP[T];
  vector[5] sa[T];
 
  matrix[11,5] C1;
  matrix[5,4] C1a;
  vector[K] Ca;
  vector[K] Cg;
  vector[11] C0;
  matrix[11,T] my;
  vector[11] yita;
  matrix[11,11] f;
  matrix[11,11] R = diag_matrix(2*log(eps));
  for (i in 1:11){
    for(j in 1:5){
      C1[i,j] = 0;
    }
  }
   for(i in 1:5){
    for (j in 1:4){
      C1a[i,j] = (1-exp(-kappaq[j]*tau[i]))/(kappaq[j]*tau[i]);
    }
  }
  for (j in 1:K){
     Ca[j]=(cr+ca) * tau[j] - sigma[1]^2 *(tau[j]-C1[j,1])/(2*kappaq[1]^2)+sigma[1]^2 *C1[j,1]^2/(4*kappaq[1])-
          sigma[2]^2 *(tau[j]-C1[j,2])/(2*kappaq[2]^2)+sigma[2]^2 *C1[j,2]^2/(4*kappaq[2])-
          sigma[3]^2 *(tau[j]-C1[j,3])/(2*kappaq[3]^2)+sigma[3]^2 *C1[j,3]^2/(4*kappaq[3]);
  }
  for(j in 1:K){
     Cg[j]=(cr*gr+cg) * tau[j] - gr^2*sigma[1]^2 *(tau[j]-C1[j,1])/(2*kappaq[1]^2)+gr^2*sigma[1]^2 *C1[j,1]^2/(4*kappaq[1])-
           gr^2*sigma[2]^2 *(tau[j]-C1[j,2])/(2*kappaq[2]^2)+gr^2*sigma[2]^2 *C1[j,2]^2/(4*kappaq[2])-
           sigma[4]^2 *(tau[j]-C1[j,4])/(2*kappaq[4]^2)+sigma[4]^2 *C1[j,4]^2/(4*kappaq[4]);
  }
  C0[1]=cr;
  C0[2:6]=Ca;
  C0[7:11]=Cg;
  
  C1[1,1] = 1;
  C1[1,2] = 1;
  C1[1,5] = 1;
  for (i in 2:6){
    for(j in 1:3){
    C1[i,j] = C1a[i-1,j];
    }
  }
  for (i in 1:5){
    C1[i+6,1] = C1a[i,1]*gr;
    C1[i+6,2] = C1a[i,2]*gr;
    C1[i+6,4] = C1a[i,4];
  }
  
  
  phi = diag_matrix(exp(-(kappa)*h));
  Q = diag_matrix(exp(2*log(sigma)).*(1-exp(-2*kappa*h))./(2*kappa));

  p0 = diag_matrix(exp(2*log(sigma))./(2*kappa));
  Mu =  (diag_matrix(rep_vector(1,5))-phi) *mu;
  a1 = Mu + phi*mu;
  p1 = phi * p0 * phi';
  for (i in 1:T) {
    yita = y[,i] - C0 - C1 *  a1;
    f = C1 * p1 * C1' + R;
    a[i] =  a1 + p1 * C1' * inverse(f) * yita;
    P[i] = p1 - p1 * C1' * inverse(f) * C1 * p1;
    a1 = Mu + phi * a[i];
    p1 = phi * P[i] * phi' + Q;
    
  }
  sa[T] = a[T];
  sP[T] = P[T];
  for (j  in (T-1):1) {
    
    sa[j]  = a[j] + P[j] * phi' * inverse(phi * P[j] * phi' +Q) * (aa[j+1] - Mu - phi * a[j]);
    sP[j]  = P[j] - P[j] * phi' * inverse(phi * P[j] * phi'+ Q) * phi * P[j];
    
   
    
  }
 
 
  for (t in 1:T){
    my[,t] = C0 +C1*aa[t];

  }
}
  
model {
  kappa~normal(0,5);
  kappaq~normal(0,5);
  mu ~normal(0,5);
  sigma ~normal(0,1);

 for (t in 1:T){
    aa[t] ~ multi_normal_cholesky(sa[t],cholesky_decompose(0.5*(sP[t]+sP[t]')));
    y[,t] ~ multi_normal_cholesky(my[,t],diag_matrix(eps));
 }
 
}

And the wrong messages are:

Chain 1: Rejecting initial value:
Chain 1: Error evaluating the log probability at the initial value.
Chain 1: Exception: cholesky_decompose: m is not symmetric. m[1,2] = nan, but m[2,1] = nan (in ‘string’, line 126, column 4 to column 80)

I have no idea to solve this problem.

some suggestions?

2

Does cholesky decompose accept diagonal matrices?

3

sorry, I don’t understand your question. But actually I make cholesky decompose in diagonal matrix can really work in stan.If you want me to use cholesky decompose in

this line. I try it, but same problem exist.