I'm trying to set up an AR(2) model in the dlm context. I've generated a
time series utilizing the code:
am = 800; #sample size
des = 200; #initial values to be discarded
V = 0.5
v = rnorm((am+des+1),0,sqrt(V))
W = 0.9
w = rnorm((am+des+1),0,sqrt(W))
U = 0.9
u = rnorm((am+des+1),0,sqrt(U))
phi1 = 0.6;
phi2 = -0.4;
mu=matrix(0,nrow=(am+des+1))
mu[2,1] = 10;
x=matrix(0,nrow=(am+des+1))
x[1,1] = 0;
x[2,1] = 0;
#----------------------------------------------------------
yg=NULL
for (i in 3:(am+des+1)) {
mu[i] = mu[i-1] + w[i]
x[i] = phi1*x[i-1] + phi2*x[i-2] + u[i]
yg[i] = mu[i] + x[i] + v[i]
}
y=NULL
for (i in 1:(am + 1)) {
y[i] = yg[(i+des)]
}
And obtained the estimates through:
buildfun = function(theta) {
dlm(FF=t(c(1,1,0)),GG=matrix(c(1,0,0,0,theta[1],theta[2],0,1,0),nrow=3,byrow=T),
V=exp(theta[3]),W=diag(c(exp(theta[4]),exp(theta[5]),0)),
m0=c(0,0,0),C0=diag(c(0.3,0.8,0.7)))
}
estMLE = dlmMLE (y, parm = c(0,0,0,0,0), build=buildfun)
phis=estMLE$par[1:2]
variances=exp(estMLE$par[3:5])
c(phis,variances)
The estimates are not very close to the values defined in the series generated.
I think there is some sort of problem in the model specified in
"buildfun", however, I cannot identify it. Any help is appreciated.
[[alternative HTML version deleted]]