kalmanFilter {pomp}R Documentation

Kalman filter

Description

The basic Kalman filter for multivariate, linear, Gaussian processes.

Usage

kalmanFilter(object, X0 = rinit(object), A, Q, C, R, tol = 1e-06)

Arguments

object

a pomp object containing data;

X0

length-m vector containing initial state. This is assumed known without uncertainty.

A

m×mm\times m latent state-process transition matrix. E[Xt+1Xt]=A.XtE[X_{t+1} | X_t] = A.X_t.

Q

m×mm\times m latent state-process covariance matrix. Var[Xt+1Xt]=QVar[X_{t+1} | X_t] = Q

C

n×mn\times m link matrix. E[YtXt]=C.XtE[Y_t | X_t] = C.X_t.

R

n×nn\times n observation process covariance matrix. Var[YtXt]=RVar[Y_t | X_t] = R

tol

numeric; the tolerance to be used in computing matrix pseudoinverses via singular-value decomposition. Singular values smaller than tol are set to zero.

Details

If the latent state is XX, the observed variable is YY, XtRmX_t \in R^m, YtRnY_t \in R^n, and

XtMVN(AXt1,Q)X_t \sim \mathrm{MVN}(A X_{t-1}, Q)

YtMVN(CXt,R)Y_t \sim \mathrm{MVN}(C X_t, R)

where MVN(M,V)\mathrm{MVN}(M,V) denotes the multivariate normal distribution with mean MM and variance VV. Then the Kalman filter computes the exact likelihood of YY given AA, CC, QQ, and RR.

Value

A named list containing the following elements:

object

the ‘pomp’ object

A, Q, C, R

as in the call

filter.mean

E[Xty1,,yt]E[X_t|y^*_1,\dots,y^*_t]

pred.mean

E[Xty1,,yt1]E[X_t|y^*_1,\dots,y^*_{t-1}]

forecast

E[Yty1,,yt1]E[Y_t|y^*_1,\dots,y^*_{t-1}]

cond.logLik

f(yty1,,yt1)f(y^*_t|y^*_1,\dots,y^*_{t-1})

logLik

f(y1,,yT)f(y^*_1,\dots,y^*_T)

See Also

enkf, eakf

Examples

 # takes too long for R CMD check

  if (require(dplyr)) {

    gompertz() -> po

    po |>
      as.data.frame() |>
      mutate(
        logY=log(Y)
      ) |>
      select(time,logY) |>
      pomp(times="time",t0=0) |>
      kalmanFilter(
        X0=c(logX=0),
        A=matrix(exp(-0.1),1,1),
        Q=matrix(0.01,1,1),
        C=matrix(1,1,1),
        R=matrix(0.01,1,1)
      ) -> kf

    po |>
      pfilter(Np=1000) -> pf

    kf$logLik
    logLik(pf) + sum(log(obs(pf)))
    
  }


[Package pomp version 5.11.0.0 Index]