MARSSkem() performs maximum-likelihood estimation, using an EM algorithm for constrained and unconstrained MARSS models. Users would not call this function directly normally. The function MARSS calls MARsskem. However users might want to used MARSSkem directly if they need to avoid some of the error-checking overhead associated with the MARSS function.
MARSSkem(MLEobj)An object of class marssMLE.
The marssMLE object which was passed in, with additional components:
String "kem".
Kalman filter output.
If MLEobj$control$trace = TRUE, a list with par = a record of each estimated parameter over all EM iterations and logLik = a record of the log likelikelihood at each iteration.
Number of iterations needed for convergence.
Did estimation converge successfully?
Converged in both the abstol test and the log-log plot test.
Some of the parameter estimates did not converge (based on the log-log plot test AND abstol tests) before MLEobj$control$maxit was reached. This is not an error per se.
No convergence diagnostics were computed because the MLE object had problems and was not fit. This isn't a convergence error just information.
No convergence diagnostics were computed because the MLE object was not fit. This isn't a convergence error just information.
Abstol convergence only. Some of the parameter estimates did not converge (based on the log-log plot test) before MLEobj$control$maxit was reached. However MLEobj$control$abstol was reached.
Log-log convergence only. Some of the parameter estimates did not converge (based on the abstol test) before MLEobj$control$maxit was reached. However the log-log convergence test was passed.
Abstol convergence only. Log-log convergence test was not computed because MLEobj$control$maxit was set to less than control$min.iter.conv.test.
Lack of convergence info. Parameter estimates did not converge based on the abstol test before MLEobj$control$maxit was reached. No log-log information since control$min.iter.conv.test is less than MLEobj$control$maxit so no log-log plot test could be done.
MLEobj$control$abstol was reached but the log-log plot test returned NAs. This is an odd error and you should set control$trace=TRUE and look at the outputted $iter.record to see what is wrong.
The EM algorithm was abandoned due to numerical errors. Usually this means one of the variances either went to zero or to all elements being equal. This is not an error per se. Most likely it means that your model is not very good for your data (too inflexible or too many parameters). Try setting control$trace=1 to view a detailed error report.
The algorithm was abandoned due to errors in the log-log convergence test. You should not get this error (it is included for debugging purposes to catch improper arguments passed into the log-log convergence test).
The algorithm was run for control$maxit iterations, control$abstol not reached, and the log-log convergence test returned errors. You should not get this error (it is included for debugging purposes to catch improper arguments passed into the log-log convergence test).
Other convergence errors. This is included for debugging purposes to catch misc. errors.
Log-likelihood.
State estimates from the Kalman filter.
Confidence intervals based on state standard errors, see caption of Fig 6.3 (p. 337) Shumway & Stoffer.
Any error messages.
To ensure that the global maximum-likelihood values are found, it is recommended that for your final fit initial parameter values be set using Monte Carlo initialization (MLEobj$control$MCInit = TRUE), particularly if the model is not a good fit to the data. This requires more compuation time, but reduces the chance of the algorithm terminating at a local maximum and not reaching the true MLEs. For many models and for draft analyses, this is unnecessary, but answers should be checked using an initial conditions search before reporting final values.
MARSSkem() calls a Kalman filter/smoother (MARSSkf) for hidden state estimation. The algorithm allows two options for the initial state conditions: fixed but unknown or a prior. In the first case, x0 (whether at t=0 or t=1) is treated as fixed but unknown (estimated); in this case, fixed$V0=0 and x0 is estimated. This is the default behavior. In the second case, the initial conditions are specified with a prior and V0!=0. In the later case, x0 or V0 may be estimated. MARSS will allow you to try to estimate both, but many researchers have noted that this is not robust so you should fix one or the other.
If you get errors, it generally means that the solution involves an ill-conditioned matrix. For example, your Q or R matrix is going to a value in which all elements have the same value, for example zero. If for example, you tried to fit a model with fixed and high R matrix and the variance in that R matrix was much higher than what is actually in the data, then you might drive Q to zero. Also if you try to fit a structurally inadequate model, then it is not unusual that Q will be driven to zero. For example, if you fit a model with 1 hidden state trajectory to data that clearly have 2 quite different hidden state trajectories, you might have this problem. Comparing the likelihood of this model to a model with more structural flexibility should reveal that the structually inflexible model is inadequate (much lower likelihood).
Convergence testing is done via a combination of two tests. The first test (abstol test) is the test that the change in the absolute value of the log-likelihood from one iteration to another is less than some tolerance value (abstol). The second test (log-log test) is that the slope of a plot of the log of the parameter value or log-likelihood versus the log of the iteration number is less than some tolerance. Both of these must be met to generate the Success! parameters converged output. If you want to circumvent one of these tests, then set the tolerance for the unwanted test to be high. That will guarantee that that test is met before the convergence test you want to use is met. The tolerance for the abstol test is set by control$abstol and the tolerance for the log-log test is set by control$conv.test.slope.tol. Anything over 1 is huge for both of these.
Objects of class marssMLE may be built from scatch but are easier to construct using MARSS with MARSS(..., fit=FALSE).
Options for MARSSkem() may be set using MLEobj$control. The commonly used elements of control are follows (see marssMLE:
minitMinimum number of EM iterations. You can use this to force the algorithm to do a certain number of iterations. This is helpful if your soln is not converging.
maxitMaximum number of EM iterations.
min.iter.conv.testThe minimum number of iterations before the log-log convergence test will be computed. If maxit is set less than this, then convergence will not be computed (and the algorithm will just run for maxit iterations).
kf.x0Whether to set the prior at t=0 ("x00") or at t=1 ("x10"). The default is "x00".
conv.test.deltaTThe number of iterations to use in the log-log convergence test. This defaults to 9.
abstolTolerance for log-likelihood change for the delta logLik convergence test. If log-likelihood changes less than this amount relative to the previous iteration, the EM algorithm exits. This is normally (default) set to NULL and the log-log convergence test is used instead.
allow.degenWhether to try setting Q or R elements to zero if they appear to be going to zero.
traceA positive integer. If not 0, a record will be created of each variable over all EM iterations and detailed warning messages (if appropriate) will be printed.
safeIf TRUE, MARSSkem will rerun MARSSkf after each individual parameter update rather than only after all parameters are updated. The latter is slower and unnecessary for many models, but in some cases, the safer and slower algorithm is needed because the ML parameter matrices have high condition numbers.
MCInitIf TRUE, Monte Carlo initialization will be performed by MARSSmcinit.
numInitsNumber of random initial value draws to be used with MARSSmcinit. Ignored if MCInit=FALSE.
numInitStepsMaximum number of EM iterations for each random initial value draw to be used with MARSSmcinit. Ignored if MCInit=FALSE.
boundsInitsLength 6 list. Each component is a length 2 vector of bounds on the uniform distributions (for B, U, and A) from which initial values will be drawn to be used with MARSSmcinit(). For R and Q, the 2 bounds specify the df and S(=diag(bound[2],m)) for a wishart distribution. Ignored if MCInit=FALSE. See Examples.
silentSuppresses printing of progress bars, error messages, warnings and convergence information.
R. H. Shumway and D. S. Stoffer (2006). Chapter 6 in Time Series Analysis and its Applications. Springer-Verlag, New York.
Ghahramani, Z. and Hinton, G. E. (1996) Parameter estimation for linear dynamical systems. Technical Report CRG-TR-96-2, University of Totronto, Dept. of Computer Science.
Harvey, A. C. (1989) Chapter 5 in Forecasting, structural time series models and the Kalman filter. Cambridge University Press, Cambridge, UK.
The user guide: Holmes, E. E., E. J. Ward, and M. D. Scheuerell (2012) Analysis of multivariate time-series using the MARSS package. NOAA Fisheries, Northwest Fisheries Science
Center, 2725 Montlake Blvd E., Seattle, WA 98112 Type RShowDoc("UserGuide",package="MARSS") to open a copy.
Holmes, E. E. (2012). Derivation of the EM algorithm for constrained and unconstrained multivariate autoregressive
state-space (MARSS) models. Technical Report. arXiv:1302.3919 [stat.ME]RShowDoc("EMDerivation",package="MARSS") to open a copy.
# NOT RUN {
dat = t(harborSeal)
dat = dat[2:4,]
#you can use MARSS to construct a proper marssMLE object.
MLEobj = MARSS(dat, model=list(Q="diagonal and equal", U="equal"), fit=FALSE)
#Pass this MLEobj to MARSSkem to do the fit.
kemfit = MARSSkem(MLEobj)
# }
Run the code above in your browser using DataLab