approxSSM: Linear Gaussian Approximation for Exponential Family State Space Model

Description

Function approxSMM performs a linear Gaussian approximation of an exponential family state space model.

Usage

approxSSM(model, theta, maxiter = 50, tol = 1e-08)

Arguments

model

A non-Gaussian state space model object of class SSModel.

theta

Initial values for conditional mode theta.

maxiter

The maximum number of iterations used in approximation. Default is 50.

tol

Tolerance parameter for convergence checks.

Value

An object of class SSModel which contains the approximating Gaussian state space model with following additional components:

thetahat

Mode of $p(\theta|y)$.

iterations

Number of iterations used.

difference

Relative difference in the last step of approximation algorithm.

Details

This function is rarely needed itself, it is mainly available for illustrative and debugging purposes. The underlying Fortran code is used by other functions of KFAS for non-Gaussian state space modelling.

The linear Gaussian approximating model is defined by $$\tilde y_t = Z_t \alpha_t + \epsilon_t, \quad \epsilon_t \sim N(0,\tilde H_t),$$ $$\alpha_{t+1} = T_t \alpha_t + R_t \eta_t, \quad \eta_t \sim N(0,Q_t),$$

and $\alpha_1 \sim N(a_1,P_1)$,

where $\tilde y$ and $\tilde H$ are chosen in a way that the linear Gaussian approximating model has the same conditional mode of $\theta=Z\alpha$ given the observations $y$ as the original non-Gaussian model. Models also have a same curvature at the mode.

The approximation of the exponential family state space model is based on iterative weighted least squares method, see McCullagh and Nelder (1983) p.31 and Durbin Koopman (2012) p. 243.

References

McCullagh, P. and Nelder, J. A. (1983). Generalized linear models. Chapman and Hall.
Koopman, S.J. and Durbin, J. (2012). Time Series Analysis by State Space Methods. Second edition. Oxford University Press.

Examples

Run this code

# NOT RUN {
# A Gamma example modified from ?glm (with log-link)
clotting <- data.frame(
  u = c(5,10,15,20,30,40,60,80,100),
  lot1 = c(118,58,42,35,27,25,21,19,18),
  lot2 = c(69,35,26,21,18,16,13,12,12))

glmfit1 <- glm(lot1 ~ log(u), data = clotting, family = Gamma(link = "log"))
glmfit2 <- glm(lot2 ~ log(u), data = clotting, family = Gamma(link = "log"))

# Model lot1 and lot2 together (they are still assumed independent)
# Note that Gamma distribution is parameterized by 1/dispersion i.e. shape parameter
model <- SSModel(cbind(lot1, lot2) ~ log(u),
                u = 1/c(summary(glmfit1)$dispersion, summary(glmfit2)$dispersion),
                data = clotting, distribution = "gamma")
approxmodel <- approxSSM(model)

# Conditional modes of linear predictor:
approxmodel$thetahat
cbind(glmfit1$linear.predictor, glmfit2$linear.predictor)

KFS(approxmodel)
summary(glmfit1)
summary(glmfit2)

# approxSSM uses modified step-halving for more robust convergence than glm:
y <- rep (0:1, c(15, 10))
suppressWarnings(glm(formula = y ~ 1, family = binomial(link = "logit"), start = 2))
model <- SSModel(y~1, dist = "binomial")
KFS(model, theta = 2)
KFS(model, theta = 7)
# }

Run the code above in your browser using DataCamp Workspace