MARSSkem()
performs maximum-likelihood estimation, using an EM algorithm for constrained and unconstrained MARSS models. This is one of the base functions in the MARSS-package
.MARSSkem(MLEobj)
marssMLE
.marssMLE
object which was passed in, with additional components: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.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.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
:
[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object]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. RShowDoc("EMDerivation",package="MARSS")
to open a copy.MARSSmcinit
, MARSSkf
, marssMLE
, MARSSoptim
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