[R] Package dlm generates unstable results?

[Taste of R]

 
Hi, All,
 
This is the first time I seriously use this package. However, I am confused 
that the result is quite unstable. Maybe I wrote something wrong in the code? 
So could anybody give me some hint? Many thanks.
 
My test model is really simple.
Y_t = X_t * a_t + noise(V)，(no Intercept here)
a_t = a_{t-1} + noise(W)
 
I first run the following code: (I shall provide data at the end of the mail)
 
BuildMod <- function(x){
 return(dlm(
  m0  = x[1],
  C0  = x[2],
  FF  = 1,
  GG  = 1,
  V   = x[3],
  W   = x[4],
  JFF = 1,
  X   = X
  ))
}
ModFit  <- dlmMLE(Y,rep(1,4),BuildMod,debug=T)
dlmMod  <- BuildMod(ModFit$par)
V <- dlmMod$V
W <- dlmMod$W
ModFilt <- dlmFilter(Y,dlmMod)
v <- tail(dlmSvd2var(ModFilt$U.C,ModFilt$D.C),1)
m <- tail(ModFilt$m,1)
 
The results are:
V = 5.945003e-05
W =  0.0003086623

v =  9.850526e-05 (the estimated variance for a_t  after we observe the last 
pair of observations)
m = -0.02965614 (the estimated mean for a_t after we observe the last pair of 
observations)

Since I observe that all the numbers are quite small, I re-estimate the whole 
thing using a smaller starting value, ie.,
 
ModFit <- dlmMLE(Y,rep(0.1,4),BuildMod,debug=T)
 
I was hoping to get a better estimate with a staring value closer to the final 
estimates and checkk the robustness of the result. However, I got the following 
message:
 
Error in dlm(m0 = x[1], C0 = x[2], FF = 1, GG = 1, V = x[3], W = x[4],  : 
  C0 is not a valid variance matrix
 
I do not know exactly how to debug this. It seems to me that during MLE, 
somehow the C0 value was turned into a negative number. This is quite scary, So 
would anybody give me some hints? Thanks a lot. Here is the data for Y,
 
0.0125678739370109
-0.00241285475528163
0.00386919876129071
-0.00352839097011217
0.00285344714211614
0.00374266510625097
0.00797807743013259
-0.00543459628953192
-0.0138447399853609
-0.0102614592879934
-0.0225111772602310
-0.0127304918143123
-0.00730849659351113
-0.0206703167742092
0.0228898867615212
-0.000489089315662759
0.00760340725359960
 
And here is the data for X
0.481797086735748
-0.336701996049702
-0.403677908907445
-0.660006432637389
-0.598885119922226
-0.386026586737966
-0.861884498592061
-0.575467614543903
-1.09697129563504
0.156970856187597
0.816563280464663
0.223472270913202
-0.901251445288375
0.402507298840520
0.655744153537491
1.97576880567968
1.05416962424537
0.128426774782304
0.0930179549205073
-0.120974956861827
0.446329808606182
1.15541338624045
0.658459632511477



        [[alternative HTML version deleted]]


______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.