MARSS() outputs marssMLE objects. fitted(MLEobj), where MLEobj is the output from a MARSS() call, will return the modeled value of \(\mathbf{y}_t\) or \(\mathbf{x}_t\). For \(\mathbf{y}_t\), this is \(\mathbf{Z}_t\tilde{\mathbf{x}}_t+\mathbf{a}_t\). For \(\mathbf{x}_t\), this is \(\mathbf{B}_t\tilde{\mathbf{x}}_{t-1}+\mathbf{u}_t\). If you want the estimate of \(\mathbf{x}_t\), then use tidy.marssMLE.
# S3 method for marssMLE
fitted(object, ...,
type = c("ytT", "xtT", "ytt", "ytt1", "xtt1"),
interval = c("none", "confidence", "prediction"),
conf.level = 0.95,
output = c("tibble", "matrix"))A marssMLE object.
Other arguments. Used for backwards compatibility with old arguments.
Fitted values for the observations (y) or the states (x). If 'tT', then the estimate at time 't' is conditioned on all the data. If 'tt', then the estimate is conditioned on data up to time t. If 'tt1', estimate is conditioned on data up to time t-1. This is also known as one-step-ahead estimate or for y, the innovations.
If interval="confidence", then the standard error and confidence interval of the fitted value is returned.
Level for the intervals if interval != "none".
tibble or list of matrices
If interval="none" (the default), a T column matrix of the fitted values with one row for each observation (or state) time series is returned.
If interval = "confidence", the following are returned in a list:
Fitted values of observations (y) or states (x). See details.
Standard errors of fitted values
Lower confidence level at alpha = 1-conf.level. The interval is approximated using qnorm(alpha/2)*.se.fit + .fitted.
Upper confidence level. The interval is approximated using qnorm(1-alpha/2)*.se.fit + .fitted.
The confidence interval is for the fitted value, i.e. \mathbf{Z}_t\tilde{\mathbf{x}}_t+\mathbf{a}_tZ(t) tildex(t)+a(t) or \mathbf{B}_t\tilde{\mathbf{x}}_{t-1}+\mathbf{u}_tB(t) tildex(t-1)+u(t).
If interval = "prediction", the following are returned in a list:
Fitted values of observations (y) or states (x). See details.
Standard deviation of new \(\mathbf{x}_t\) or \(\mathbf{y}_t\) iid values.
Lower range at alpha = 1-conf.level. The interval is approximated using qnorm(alpha/2)*.sd + .fitted.
Upper range at level. The interval is approximated using qnorm(1-alpha/2)*.sd + .fitted.
The prediction interval is for new \mathbf{x}_tx(t) or \mathbf{y}_ty(t), not estimated \mathbf{x}_tx(t) nor \mathbf{y}_ty(t) used in the model. Do not plot the observed data nor the states estimates on these intervals. For that you need the residuals intervals provided by augment.marssMLE.
observation fitted values
The model predicted (fitted) \(\hat{\mathbf{y}}_t\) is \(\mathbf{Z}_t\tilde{\mathbf{x}}_t+\mathbf{a}_t\), where the model is written in marss form. See MARSS.marss for a discussion of the conversion of MARSS models with covariates (\(\mathbf{c}\) and \(\mathbf{d}\)) into marss form which is how models are written in the internal MARSS algorithms).
\(\tilde{\mathbf{x}}_t\) is the expected value of the states at time \(t\). If type="ytT", \(\tilde{\mathbf{x}}_t\) is the expected value conditioned on all the data, i.e. xtT returned by MARSSkf(). If type="ytt1", then expected value uses only the data up to time \(t-1\), i.e. xtt1 returned by MARSSkf(). These are commonly known as the one step ahead predictions for a state-space model. If type="ytt", then the expected value uses the data up to time \(t\), i.e. xtt returned by MARSSkf().
If interval="confidence", the se and interval is for the fitted \(\mathbf{y}\). The standard error of the fitted values is \(\mathbf{Z}_t \tilde{\mathbf{V}}_t \mathbf{Z}_t^\top\). If interval="prediction", the standard deviation of new iid \(\mathbf{y}\) datasets is returned. The standard deviation of new \(\mathbf{y}\) is \(\mathbf{Z}_t \tilde{\mathbf{V}}_t \mathbf{Z}_t^\top\). \(\tilde{\mathbf{V}}_t\) is either conditioned on 1:T, 1:t, or 1:t-1 depending on type. Do not plot the data used in the model on these intervals. For that you want the conditional model residuals (and se's). Use augment.marssMLE.
state fitted values
The model predicted \(\mathbf{x}(t)\) given \(\mathbf{x}_{t-1}\) is \(\mathbf{B}(t)\tilde{\mathbf{x}}(t-1)+\mathbf{u}(t)\), where the model is written in "marss" form (MARSS.marss). This type of state fitted value is used for process outlier detection and shock detection. See residuals.marssMLE and read the references cited.
If you want estimates of the states, rather than the model predictions based on \(\mathbf{x}_{t-1}\) then you'll want either the states estimate conditioned on all the data (or conditioned on the data up to time \(t-1\) or up to time \(t\)). These are returned by MARSSkf() in xtT, xtt1 and xtt respectively. Which one you want depends on your objective and application. You can also use the tidy.marssMLE() function to return a data.frame (tibble) with the estimated states with standard errors and intervals.
\(\tilde{\mathbf{x}}_{t-1}\) used in the prediction is the expected value of the states at time \(t-1\). If type="xtT", this is the expected value at time \(t-1\) conditioned on all the data, i.e. xtT[,t-1] returned by MARSSkf(). If type="xtt1", it is the expected value conditioned on the data up to time \(t-1\), i.e. xtt[,t-1] returned by MARSSkf(). The fitted state values conditioned on data up to \(t\) is not provided. This would require the expected value of states at time \(t\) conditioned on data up to time \(t+1\), and this is not output by the Kalman filter. Only conditioning on data up to \(t-1\) and T are output.
The intervals returned by fitted.marssMLE for the fitted states, are not typically what one uses or needs--however might be useful for simulation work. If you are doing outlier detection or shock detection, you need the intervals on the smoothed state residuals which are returned by residuals.marssMLE or augment.marssMLE. If you want intervals on the states estimates, use tidy.marssMLE.
If interval="confidence", the standard error of the fitted values (meaning the 'expected value of \(\mathbf{X}_t\)' given \(\mathbf{X}_{t-1}\)) is returned. The standard error of the fitted value is \(\mathbf{B}_t \tilde{\mathbf{V}}_{t-1} \mathbf{B}_t^\top\). If interval="prediction", the standard deviation of \(\mathbf{X}_t\) given \(\mathbf{X}_{t-1}\) is output. The latter is \(\mathbf{B}_t \tilde{\mathbf{V}}_{t-1} \mathbf{B}_t^\top + \mathbf{Q}\) (notice it includes Q). \(\tilde{\mathbf{V}}_{t-1}\) is either conditioned on 1:T or 1:t-1 depending on type. Do not plot (or compare) the estimate of \(\mathbf{x}_t\) to the intervals for the fitted values. You need the conditional states residuals intervals in that case; use augment.marssMLE.
# NOT RUN {
dat <- t(harborSeal)
dat <- dat[c(2, 11, 12), ]
fit <- MARSS(dat, model = list(Z = factor(c("WA", "OR", "OR"))))
fitted(fit)
# Example of fitting a stochastic level model to the Nile River flow data
# red line is smooothed level estimate
# grey line is one-step-ahead prediction using xtt1
nile <- as.vector(datasets::Nile)
mod.list <- list(
Z = matrix(1), A = matrix(0), R = matrix("r"),
B = matrix(1), U = matrix(0), Q = matrix("q"),
x0 = matrix("pi")
)
fit <- MARSS(nile, model = mod.list, silent = TRUE)
plot(nile, type = "p", pch = 16, col = "blue")
lines(fitted(fit, type="ytT")[1, ], col = "red", lwd = 2)
lines(fitted(fit, type="ytt1")[1, ], col = "grey", lwd = 2)
# }
Run the code above in your browser using DataLab