run KFAS

testrun_KFAS.R

library(KFAS)
set.seed(123456)
alpha <- beta <- numeric(100)
alpha[1:2] <- rnorm(2)
for(i in 3:100){
  alpha[i] <- rnorm(1,2*alpha[i-1]-alpha[i-2],0.2)
}
beta <- rnorm(100,alpha,0.2)
y1<-rnbinom(100,size=10,mu=exp(alpha))
y2<-rnbinom(100,size=10,mu=exp(beta))
plot(y2,type="l")
lines(y1,type="l",col="blue")

Zt=matrix(c(1,0,0,
            0,0,1),byrow=TRUE,2,3)
Ht=matrix(0,3,3)
Tt=matrix(c(2,-1,0,
            1,0,0,
            1,0,0),3,3,byrow = TRUE)
Qt=matrix(c(NA,0,0,
            0,0,0,
            0,0,NA),byrow = TRUE,3,3)
a1 = c(0,0,0)
P1 <- matrix(0, 3, 3)
P1inf <- diag(3)
mod1 <-SSModel(cbind(c(y1[1:90],rep(NA,10)),c(y2[1:90],rep(NA,10)))~-1+SSMcustom(Z=Zt,T=Tt,Q=Qt,a1=a1),distribution = rep("negative binomial",2))

diag(mod1$P1inf)
diag(mod1$P1) <- 100^2

updatefn <- function(pars,model){
  model$Q[1,1,] <- exp(pars[1])
  model$Q[3,3,] <- exp(pars[2])
  model$u[,] <- exp(pars[3:4])
  return(model)
}

fit1 <-fitSSM(mod1,inits = c(1,1,1,1),updatefn = updatefn,nsim=1000)
pre1 <-predict(fit1$model,nsim = 2000)
plot(y1)
lines(pre1$y1,col="blue",lwd=2)
plot(y2)
lines(pre1$y2,col="blue",lwd=2)