These functions run the iterative equations of the Kalman smoother
for a state space model upon the output from the Kalman filter.
Usage
KS(y, ss, kf)
KS.deriv(y, ss, kf)
Arguments
y
a numeric time series or vector.
ss
a list containing the matrices of the state space model.
kf
a list containing the output returned by the Kalman filter
KF.
Details
See the details section and the section ‘state space representation’
in KF.
Missing observations are allowed.
The input kf passed to KS.deriv must contain the
derivative terms related to the filter that are returned by
KF.deriv or KF.deriv.C.
When the Kalman filter was found to convergence at some iteration,
i.e., kf$convit is not null, these functions use steady state values for
N and varahat in the intermediate iterations of the smoother.
For example, if the filter converged at iteration 15 in a series of length n,
the equations of the smoother are run for the first iterations from observation n
to n-15; then the steady state values are used until there are 15 iterations
remaining. In the last iterations, from observation 15 to 1 the equations of
the smoother are evaluated again.
In practice, if the disturbance smoother is to be run as well, using the functions
described in KFKSDS will be slightly more efficient.
Value
A list containing the following elements:
ahat
smoothed state vector.
varhat
covariance matrix of ahat.
r
weighted sum of innovations used to obtain ahat.
N
intermediate matrix used to obtain varahat
The function KS.deriv returns also the derivatives referred to each
of the elements defined above, named respectively dahat, dvarahat,
dr and dN.
References
Durbin, J. and Koopman, S. J. (2001).
Time Series Analysis by State Space Methods.
Oxford University Press.
Harvey, A. C. (1989).
Forecasting, Structural Time Series Models and the Kalman Filter.
Cambridge University Press.
See Also
KF, KSDS;
char2numeric in package stsm.
Examples
# local level plus seasonal model with arbitrary parameter values
# for the 'JohnsonJohnson' time series
m <- stsm::stsm.model(model = "llm+seas", y = JohnsonJohnson,
pars = c("var1" = 2, "var2" = 15, "var3" = 30))
ss <- stsm::char2numeric(m)
kf <- KF(m@y, ss)
ks <- KS(m@y, ss, kf)
plot(ks$ahat[,1:2], main = "smoothed state vector")
kfd <- KF.deriv(m@y, ss)
ksd <- KS.deriv(m@y, ss, kfd)
all.equal(ks$ahat, ksd$ahat)
# extended output is required if 'KF.deriv.C' is used to obtain
# the necessary elements from the filter, set return.all = TRUE
kfdc <- KF.deriv.C(m@y, ss, return.all = TRUE)
ksd <- KS.deriv(m@y, ss, kfdc)
all.equal(ks$ahat, ksd$ahat)
# compare analytical and numerical derivatives
# yield same results up to a tolerance error
fcn <- function(x, model, type, i)
{
m <- stsm::set.pars(model, x)
ss <- stsm::char2numeric(m)
kf <- KF(m@y, ss)
ks <- KS(m@y, ss, kf)
switch(type, "ahat" = sum(ks$ahat[,i]), "r" = sum(ks$r[,i]))
}
dahat <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "ahat", i = 1)
all.equal(dahat, colSums(ksd$dahat[,1,]))
dahat <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "ahat", i = 2)
all.equal(dahat, colSums(ksd$dahat[,2,]))
dahat <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "ahat", i = 3)
all.equal(dahat, colSums(ksd$dahat[,3,]))
dr <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "r", i = 1)
all.equal(dr, colSums(ksd$dr[,1,]), check.attributes = FALSE)
dr <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "r", i = 2)
all.equal(dr, colSums(ksd$dr[,2,]), check.attributes = FALSE)
dr <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "r", i = 3)
all.equal(dr, colSums(ksd$dr[,3,]), check.attributes = FALSE)
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.
R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.
Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.
> library(KFKSDS)
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/KFKSDS/KS.Rd_%03d_medium.png", width=480, height=480)
> ### Name: KS
> ### Title: Kalman Smoother for State Space Models
> ### Aliases: KS KS.deriv
> ### Keywords: ts, model
>
> ### ** Examples
>
> # local level plus seasonal model with arbitrary parameter values
> # for the 'JohnsonJohnson' time series
> m <- stsm::stsm.model(model = "llm+seas", y = JohnsonJohnson,
+ pars = c("var1" = 2, "var2" = 15, "var3" = 30))
> ss <- stsm::char2numeric(m)
>
> kf <- KF(m@y, ss)
> ks <- KS(m@y, ss, kf)
>
> plot(ks$ahat[,1:2], main = "smoothed state vector")
>
> kfd <- KF.deriv(m@y, ss)
> ksd <- KS.deriv(m@y, ss, kfd)
> all.equal(ks$ahat, ksd$ahat)
[1] TRUE
>
> # extended output is required if 'KF.deriv.C' is used to obtain
> # the necessary elements from the filter, set return.all = TRUE
> kfdc <- KF.deriv.C(m@y, ss, return.all = TRUE)
> ksd <- KS.deriv(m@y, ss, kfdc)
> all.equal(ks$ahat, ksd$ahat)
[1] TRUE
>
> # compare analytical and numerical derivatives
> # yield same results up to a tolerance error
> fcn <- function(x, model, type, i)
+ {
+ m <- stsm::set.pars(model, x)
+ ss <- stsm::char2numeric(m)
+ kf <- KF(m@y, ss)
+ ks <- KS(m@y, ss, kf)
+ switch(type, "ahat" = sum(ks$ahat[,i]), "r" = sum(ks$r[,i]))
+ }
>
> dahat <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "ahat", i = 1)
> all.equal(dahat, colSums(ksd$dahat[,1,]))
[1] "Mean relative difference: 3.093305e-07"
> dahat <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "ahat", i = 2)
> all.equal(dahat, colSums(ksd$dahat[,2,]))
[1] "Mean relative difference: 8.430606e-08"
> dahat <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "ahat", i = 3)
> all.equal(dahat, colSums(ksd$dahat[,3,]))
[1] "Mean relative difference: 1.467548e-06"
> dr <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "r", i = 1)
> all.equal(dr, colSums(ksd$dr[,1,]), check.attributes = FALSE)
[1] TRUE
> dr <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "r", i = 2)
> all.equal(dr, colSums(ksd$dr[,2,]), check.attributes = FALSE)
[1] TRUE
> dr <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "r", i = 3)
> all.equal(dr, colSums(ksd$dr[,3,]), check.attributes = FALSE)
[1] TRUE
>
>
>
>
>
> dev.off()
null device
1
>