Модель функционального смешанного эффекта с использованием KFAS из-за низкой точности фиксированных смешанных эффектов

Question

Я попытался подобрать функциональную модель смешанных эффектов (с функциональными фиксированными эффектами и функциональными случайными эффектами), используя представление модели пространства состояний, которое было установлено фильтром Калмана. Я думаю, что функция KFS в r позаботится об одномерном подходе фильтрации (предназначенном для многомерной модели пространства состояний). Результаты показывают, что модель не очень подходит. На рисунке прикрепленной ссылки красные линии на графиках представляют собой смешанные эффекты, а черные линии - истинные смешанные эффекты. Следующий код в основном реализовал методологию в Guo (2002). Может кто-нибудь сказать мне, что я указал неправильно для модели? Большое спасибо!!!

p=2 #two functional fixed effects
q=1 #one functional random intercept
m=30 #number of time points 
n=20 #n subjects
ymm = t(ym) #ymm is an m-by-n matrix specified for state space model in KFS
fj.list = sapply(1:m, function(i){
 if(i == 1) dt = t[i] 
 else dt = t[i] - t[i-1]
 fjj = getfj(Xj = X, Zj = z, p=2, q=1, n = n) 
 return(fjj)
 }
)
hj.list = sapply(1:m, function(i){
  if(i == 1) dt = t[i] 
  else dt = t[i] - t[i-1]
  hjj = gethj(Xj = ymm[m,], Zj = z, p=2, q=1, n = n, dt = dt) 
  return(hjj)
}
)

wblock = matrix(c(dt^3/3, dt^2/2,dt^2/2, dt),ncol = 2)

w = array(0, c(2*(p+q*n), 2*(p+q*n),m))

w[1:(2*(p+q*n)),1:(2*(p+q*n)),] =as.matrix(bdiag(rep(list(wblock), p+q*n)))

ww = 2*p+2*q*n

Ht = array(0, c(2*(p+q*n), 2*(p+q*n),m))

Ht[,,] = hj.list[[1]]


fj.list2 = sapply(1:m, function(i){
  if(i == 1) dt = t[i] 
  else dt = t[i] - t[i-1]
  fjj = getfj(Xj = X, Zj = z, p=2, q=1, n = n) 
  return(fjj)
}
)
oi = seq(7, 43, by=2)
fj.list2[[1]][,oi] = fj.list2[[1]][,5]

Zt =  array(0, c(n, 2*(p+q*n),m))

Zt[,,] = fj.list[[1]]
Rj = array(0, c(2*(p+q*n), 2*(p+q*n), m))
Rj[,,] = diag(2*(p+q*n))

hvar = array(0, c(n,n,1))
hvar[,,] = diag(c(1,rep(0,(n-1))))

model<-SSModel(ymm~-1+ SSMcustom(Z=Zt,
                                 T=Ht,
                                 Q = w,
                                 R = Rj),H = diag(n))



updatefn<-function(pars,model){
  model["H"] <- diag(pars[1],n)

  for(i in 1:p){ 
    model$Q[(2*(i-1)+1):(2*i), (2*(i-1)+1):(2*i),]<-model$Q[(2*(i-1)+1):(2*i), (2*(i-1)+1):(2*i), ]*pars[(i+1)]
  }
  model$Q[(2*p+1):(2*p+2*q*n), (2*p+1):(2*p+2*q*n),]<-model$Q[(2*p+1):(2*p+2*q*n), (2*p+1):(2*p+2*q*n),]*pars[4]

  model
}
fit<-fitSSM(model,inits=c(1,1,1,500),updatefn=updatefn,method="L-BFGS-B",lower = rep(10e-6,4))
out1 = KFS(fit$model,nsim = 100)

Оценки против истинной картины смешанных эффектов