KFAS (version 1.5.1)

logLik.SSModel: Log-likelihood of the State Space Model.

Description

Function logLik.SSmodel computes the log-likelihood value of a state space model.

Usage

# S3 method for SSModel
logLik(
  object,
  marginal = FALSE,
  nsim = 0,
  antithetics = TRUE,
  theta,
  check.model = TRUE,
  transform = c("ldl", "augment"),
  maxiter = 50,
  seed,
  convtol = 1e-08,
  expected = FALSE,
  H_tol = 1e+15,
  transform_tol,
  ...
)

Value

Log-likelihood of the model.

Arguments

object

State space model of class SSModel.

marginal

Logical. Compute marginal instead of diffuse likelihood (see details). Default is FALSE.

nsim

Number of independent samples used in estimating the log-likelihood of the non-Gaussian state space model. Default is 0, which gives good starting value for optimization. Only used for non-Gaussian model.

antithetics

Logical. If TRUE, two antithetic variables are used in simulations, one for location and another for scale. Default is TRUE. Only used for non-Gaussian model.

theta

Initial values for conditional mode theta. Only used for non-Gaussian model.

check.model

Logical. If TRUE, function is.SSModel is called before computing the likelihood. Default is TRUE.

transform

How to transform the model in case of non-diagonal covariance matrix \(H\). Defaults to "ldl". See function transformSSM for details.

maxiter

The maximum number of iterations used in linearisation. Default is 50. Only used for non-Gaussian model.

seed

The value is used as a seed via set.seed function. Only used for non-Gaussian model.

convtol

Tolerance parameter for convergence checks for Gaussian approximation. Only used for non-Gaussian model.

expected

Logical value defining the approximation of H_t in case of Gamma and negative binomial distribution. Default is FALSE which matches the algorithm of Durbin & Koopman (1997), whereas TRUE uses the expected value of observations in the equations, leading to results which match with glm (where applicable). The latter case was the default behaviour of KFAS before version 1.3.8. Essentially this is the difference between observed and expected information. Only used for non-Gaussian model.

H_tol

Tolerance parameter for check max(H) > tol_H, which suggests that the approximation converged to degenerate case with near zero signal-to-noise ratio. Default is very generous 1e15. Only used for non-Gaussian model.

transform_tol

Tolerance parameter for LDL decomposition in case of a non-diagonal H and transform = "ldl". See transformSSM and ldl for details.

...

Ignored.

Details

As KFAS is based on diffuse initialization, the log-likelihood is also diffuse, which coincides with restricted likelihood (REML) in an appropriate (mixed) models. However, in typical REML estimation constant term \(log|X'X|\) is omitted from the log-likelihood formula. Similar term is also missing in diffuse log-likelihood formulations of state space models, but unlike in simpler linear models this term is not necessarily constant. Therefore omitting this term can lead to suboptimal results in model estimation if there is unknown parameters in diffuse parts of Zt or Tt (Francke et al. 2011). Therefore so called marginal log-likelihood (diffuse likelihood + extra term) is recommended. See also Gurka (2006) for model comparison in mixed model settings using REML with and without the additional (constant) term. The marginal log-likelihood can be computed by setting marginal = TRUE.

Note that for non-Gaussian models with importance sampling derivative-free optimization methods such as Nelder-Mead might be more reliable than methods which use finite difference approximations. This is due to noise caused by the relative stopping criterion used for finding approximating Gaussian model.

References

Francke, M. K., Koopman, S. J. and De Vos, A. F. (2010), Likelihood functions for state space models with diffuse initial conditions. Journal of Time Series Analysis, 31: 407--414.

Gurka, M. J (2006), Selecting the Best Linear Mixed Model Under REML. The American Statistician, Vol. 60.

Casals, J., Sotoca, S., Jerez, M. (2014), Minimally conditioned likelihood for a nonstationary state space model. Mathematics and Computers in Simulation, Vol. 100.

Examples

# Example of estimating AR model with covariates, and how to deal with possible
# non-stationarity in optimization.

set.seed(1)
x <- rnorm(100)
y <- 2 * x + arima.sim(n = 100, model = list(ar = c(0.5, -0.3)))

model<- SSModel(y ~ SSMarima(ar = c(0.5, -0.3),  Q = 1) + x, H = 0)

obj <- function(pars, model, estimate = TRUE) {
  #guard against stationarity
  armamod <- try(SSMarima(ar = artransform(pars[1:2]), Q = exp(pars[3])), silent = TRUE)
  if(class(armamod) == "try-error") {
    return(Inf)
  } else {
    # use advanced subsetting method for SSModels, see ?`[.SSModel`
    model["T", states = "arima"] <- armamod$T
    model["Q", eta = "arima"]  <- armamod$Q
    model["P1", states = "arima"]  <- armamod$P1
    if(estimate) {
      -logLik(model)
    } else {
      model
    }
  }
}
fit <- optim(p = c(0.5,-0.5,1), fn = obj, model = model, method ="BFGS")

model <- obj(fit$par, model, FALSE)
model$T
model$Q
coef(KFS(model), last = TRUE)