atsa-es / MARSS

Multivariate Autoregressive State-Space Modeling with R
https://atsa-es.github.io/MARSS/
GNU General Public License v3.0
51 stars 28 forks source link

MARSSkfas when V0 not 0 #96

Open eeholmes opened 4 years ago

eeholmes commented 4 years ago

Version 3.11.5 Info 9/27/2021

This is completed but before closing completely do the following.

-Bug. Only made it happen for Z as column vector with V0 non-zero. Guessing it is in MARSSkfss() and seems to be in the smoother part.

dat <- t(harborSealWA)
dat <- dat[2:4, ] # remove the year row
# fit a model with 1 hidden state and 3 observation time series
kemfit <- MARSS(dat, model = list(
  Z=matrix(1,3,1),
  x0=matrix(0),
  V0=matrix(10000)
))
MARSSkfas(kemfit)$VtT[,,1:2]
MARSSkfss(kemfit)$VtT[,,1:2]

Gives

 MARSSkfas(kemfit)$VtT[,,1:2]
[1] 0.006583545 0.012637968
 MARSSkfss(kemfit)$VtT[,,1:2]
[1] 0.00658371 0.01263808

With Z identity, ok

kemfit <- MARSS(dat, model = list(
  x0=matrix(0,3,1),
  V0=diag(10000,3)
))
 diag(MARSSkfas(kemfit)$VtT[,,2])
[1] 0.04039229 0.01485954 0.01070933
 diag(MARSSkfss(kemfit)$VtT[,,2])
[1] 0.04039229 0.01485954 0.01070933
eeholmes commented 4 years ago

MARSSkfas(kemfit, return.lag.one=FALSE)$VtT[,,1:2] is the same so not in the stack part of code.

Happens with tinitx=1 or tinitx=0

eeholmes commented 4 years ago

Not a bug per se but numerical inaccuracy in the classic Kalman filter. When Z is a column vector, an inversion of Z %*% Vtt1 %*% t(Z) + R is done. If Vtt1 is much larger than R then that is ill-conditioned. With V0 = 10000, Vtt1 at t=1 is much larger than R and this is ill-conditioned.

MARSSkfss is not used in the fitting functions. Addressed by changing the default kf function in the residuals to use MARSSkfas.

eeholmes commented 4 years ago

Good enough. Part of the numerical problems with the classic filter/smoother. Avoid making the V0 so bad by not trying to make it diffuse.

eeholmes commented 4 years ago

Try this for Kt[,,1] which worked for the V0T problem