Skip to content

Commit 4ba8d8a

Browse files
javlacallegaborcsardi
authored andcommitted
version 1.4
1 parent cbac082 commit 4ba8d8a

File tree

7 files changed

+210
-13
lines changed

7 files changed

+210
-13
lines changed

DESCRIPTION

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
Package: KFKSDS
2-
Version: 1.3
3-
Date: 2014-01-26
2+
Version: 1.4
3+
Date: 2014-06-15
44
Title: Kalman Filter, Smoother and Disturbance Smoother
55
Description: Naive implementation of the Kalman filter, smoother and disturbance
66
smoother for state space models. An interface to other implementations of the
@@ -12,6 +12,6 @@ Suggests: dlm, dse, FKF, numDeriv, stsm.class
1212
NeedsCompilation: yes
1313
Encoding: UTF-8
1414
License: GPL-2
15-
Packaged: 2014-01-26 19:55:17 UTC; javlacalle
1615
Repository: CRAN
17-
Date/Publication: 2014-01-27 10:58:23
16+
Packaged: 2014-06-15 12:18:11 UTC; javlacalle
17+
Date/Publication: 2014-06-15 16:50:44

MD5

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,20 +1,23 @@
1-
03d1cb1d69a86bd4c516af0aa5ed9504 *DESCRIPTION
2-
d063863027eed469ca5cd54e8a391817 *NAMESPACE
1+
8a71ef5f35a27c5e7891de83b40fc3f7 *DESCRIPTION
2+
617483d584aecafedd9068bb8fb39f9d *NAMESPACE
33
9dba3cc633ea2a111e55844c4f981c11 *R/DS.R
44
37c2ebbf9911d53c2ebdb898fec54186 *R/KF-interfaces.R
55
c89965ee7739de1739216746c2c53e06 *R/KF.R
66
bae5a673da6acd61df2cfaca27648f07 *R/KFKSDS.R
77
c0544b4399abf02856d7c40fed3b9edf *R/KS.R
88
34d3874118c8807f699f17f4c3b052a3 *R/KSDS.R
9+
076178d1487e60e2fdef263c9cd3ffed *R/predict-stsmSS.R
10+
d55097663c287c802838746a525ce556 *inst/CHANGES
911
24d96fbf6f632e9eb90872985feffc92 *inst/tests/DS-derivatives.R
1012
69813e915d4209c76f8772b0f9f6d616 *inst/tests/KF-derivatives.R
1113
9d20b6238ebf1805c862640bf60c72a9 *inst/tests/KS-derivatives.R
1214
6f947dafc1ce4c3cd2a5edc1b1ffc5c0 *man/DS.Rd
1315
78f0abcd966f7cd50d3eeaa0af774808 *man/KF-interfaces.Rd
14-
3833d7097a0ef69456619b66004a01fc *man/KF.Rd
16+
ee5e5ee1725a7411ae806e35817ef1af *man/KF.Rd
1517
6a3cc26861a11821bf759e49eacd3061 *man/KFKSDS-package.Rd
1618
d918bfdcb9941a84c0cfbc7d0337ee99 *man/KFKSDS.Rd
1719
8b904ee1fa46d75ff96944c89137f208 *man/KS.Rd
20+
b0613b67bc9f2eef053573055fe2f173 *man/predict-stsmSS.Rd
1821
3c40f594bc53834ac537ffe254436c51 *src/KF-deriv.cpp
1922
882373f8b41d314460c268382df0a415 *src/KF.c
2023
87ac4254dacb55c8f1600b23e2d669f2 *src/KF.h

NAMESPACE

Lines changed: 8 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -2,17 +2,19 @@
22
useDynLib(KFKSDS)
33

44
export(
5+
DS,
6+
DS.deriv,
7+
KalmanFilter,
58
KF,
69
KF.C,
710
KF.deriv,
811
KF.deriv.C,
9-
KS,
10-
KS.deriv,
11-
DS,
12-
DS.deriv,
1312
KSDS,
1413
KSDS.deriv,
1514
KFKSDS.deriv.C,
1615
KFKSDS.deriv.steady.C,
17-
make.KF.args,
18-
KalmanFilter)
16+
KS,
17+
KS.deriv,
18+
make.KF.args)
19+
20+
S3method(predict, stsmSS)

R/predict-stsmSS.R

Lines changed: 50 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,50 @@
1+
2+
predict.stsmSS <- function(object, y, n.ahead = 12L, ...)
3+
{
4+
kf <- KF(y, object)
5+
6+
n <- length(y)
7+
Z <- object$Z
8+
mT <- object$T
9+
H <- object$H
10+
Q <- object$Q
11+
tmT <- t(mT)
12+
13+
m <- ncol(Z)
14+
a <- matrix(nrow = n.ahead + 1, ncol = m)
15+
P <- array(NA, dim = c(m, m, n.ahead + 1))
16+
pred <- var <- rep(0, n.ahead)
17+
18+
a[1,] <- if (m == 1) kf$a.upd[n] else kf$a.upd[n,]
19+
P[,,1] <- if (m == 1) kf$P.upd[n] else kf$P.upd[,,n]
20+
21+
for (i in seq_len(n.ahead))
22+
{
23+
ip1 <- i + 1
24+
a[ip1,] <- mT %*% a[i,]
25+
P[,,ip1] <- mT %*% P[,,i] %*% tmT + Q
26+
27+
pred[i] <- Z %*% a[ip1,]
28+
var[i] <- H + sum(diag(crossprod(Z) %*% P[,,ip1]))
29+
}
30+
31+
s <- frequency(y)
32+
ytsp <- tsp(y)
33+
t0 <- ytsp[2L] + 1 / s
34+
35+
a <- ts(a[-1,], start = t0, frequency = s)
36+
Pres <- ts(matrix(nrow = n.ahead, ncol = m), start = t0, frequency = s)
37+
if (m == 1)
38+
{
39+
Pres <- ts(P[-1], start = t0, frequency = s)
40+
} else {
41+
P <- P[,,-1]
42+
for (i in seq_len(n.ahead))
43+
Pres[i,] <- diag(P[,,i])
44+
}
45+
46+
pred <- ts(pred, start = t0, frequency = s)
47+
var <- ts(var, start = t0, frequency = s)
48+
49+
list(pred = pred, se = sqrt(var), a = a, P = Pres)
50+
}

inst/CHANGES

Lines changed: 27 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,27 @@
1+
2+
Version: 1.4 2014-06-15
3+
------------------------
4+
5+
o Added the method "predict" for objects of class "stsmSS"
6+
Those objects are returned by "stsm.class::char2numeric".
7+
The method is equivalent to "predict.StructTS" but the forecasts
8+
of the components are also returned. In addition, it is straightforward
9+
to use for a model fitted by any of the maximum likelihood procedures
10+
available in package "stsm".
11+
12+
Version: 1.3 2014-01-26
13+
------------------------
14+
15+
o Fixed the following warning returned by the compiler at some points
16+
in the code:
17+
"warning: ISO C++ forbids variable length array ‘df’ [-Wvla]"
18+
19+
o Removed suggested package "sspir" since it is no longer available-
20+
21+
o Removed suggested package "KFAS" since the examples using this
22+
package relied on a previous version of "KFAS".
23+
24+
Version: 1.1 2014-01-25
25+
------------------------
26+
27+
o First version submitted to CRAN.

man/KF.Rd

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -194,11 +194,13 @@ all.equal(kf$mll, kfdc$mll)
194194
# as expected the versions that use compiled C code
195195
# are faster that the versions written fully in R
196196
# e.g. not including derivatives
197+
\dontrun{
197198
system.time(for(i in seq_len(10)) kf <- KF(m@y, ss))
198199
system.time(for(i in seq_len(10)) kfc <- KF.C(m@y, ss))
199200
# e.g. including derivatives
200201
system.time(for(i in seq_len(10)) kfd <- KF.deriv(m@y, ss))
201202
system.time(for(i in seq_len(10)) kfdc <- KF.deriv.C(m@y, ss, return.all = TRUE))
203+
}
202204

203205
# compare analytical and numerical derivatives
204206
# they give same results up to a tolerance error
@@ -221,18 +223,22 @@ all.equal(df, colSums(kfdc$df), check.attributes = FALSE)
221223
# no calls to compiled C code in either case
222224
# looking at these timings, using analytical derivatives is
223225
# expected to be useful in optimization algorithms
226+
\dontrun{
224227
system.time(for (i in seq_len(10))
225228
numdv <- numDeriv::grad(func = fcn, x = m@pars, model = m, type = "v"))
226229
system.time(for(i in seq_len(10)) kfdv <- colSums(KF.deriv(m@y, ss)$dv))
230+
}
227231

228232
# compare timings when convergence is not checked with the case
229233
# when steady state values are used after convergence is observed
230234
# computation time is reduced substantially
235+
\dontrun{
231236
n <- length(m@y)
232237
system.time(for(i in seq_len(20)) a <- KF.deriv(m@y, ss, convergence = c(0.001, n)))
233238
system.time(for(i in seq_len(20)) b <- KF.deriv(m@y, ss, convergence = c(0.001, 10)))
234239
# the results are the same up to a tolerance error
235240
all.equal(colSums(a$dv), colSums(b$dv))
236241
}
242+
}
237243

238244
\keyword{ts, model}

man/predict-stsmSS.Rd

Lines changed: 109 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,109 @@
1+
\name{predict.stsmSS}
2+
\alias{predict.stsmSS}
3+
4+
\title{Kalman Filter for State Space Models}
5+
6+
\description{
7+
These functions run the iterative equations of the Kalman filter
8+
for a state space model.
9+
}
10+
11+
\usage{
12+
\S3method{predict}{stsmSS}(object, y, n.ahead = 12L, ...)
13+
}
14+
15+
\arguments{
16+
\item{object}{a list containing the matrices of the state space model.}
17+
\item{y}{a numeric time series.}
18+
\item{n.ahead}{a numeric. The number of steps ahead to predict.}
19+
\item{...}{further arguments. Currently omitted.}
20+
}
21+
22+
\value{
23+
A list containing the following elements:
24+
item{pred}{a time series containing \code{n.ahead} predictions.}
25+
item{se}{a time series containing the standard errors of \code{pred}.}
26+
item{a}{a univariate or multivariate time series object containing \code{n.ahead}
27+
predictions for the state vector.}
28+
item{P}{a univariate or multivariate time series object containing the square of
29+
the standard errors of \code{a}.}
30+
}
31+
32+
\details{
33+
This function computes the same values as the
34+
function \link[stats]{predict.StructTS} from the \pkg{stats} package but
35+
the predictions of the components are also returned.
36+
}
37+
38+
\references{
39+
Harvey, A. C. (1989).
40+
\emph{Forecasting, Structural Time Series Models and the Kalman Filter}.
41+
Cambridge University Press.
42+
}
43+
44+
\examples{
45+
## local level model
46+
## Nile time series
47+
y <- Nile
48+
m <- stsm.class::stsm.model(model = "local-level", y = y, transPars = "StructTS")
49+
fit <- StructTS(y, "level")
50+
m <- stsm.class::set.pars(m, as.vector(fit$coef[c(2,1)]) * 100 / var(y))
51+
ss <- stsm.class::char2numeric(m, P0cov = TRUE)
52+
res <- predict(ss, y, 5)
53+
54+
# display forecasts and confidence intervals
55+
plot(cbind(y, res$pred), type = "n", plot.type = "single")
56+
lines(y)
57+
lines(res$pred, col = "blue")
58+
lines(res$pred + 2 * res$se, col = "red", lty = 2)
59+
lines(res$pred - 2 * res$se, col = "red", lty = 2)
60+
61+
# for the whole series, the above is the same as "predict.StructTS"
62+
all.equal(res$pred, predict(fit, 5)$pred)
63+
all.equal(res$se, predict(fit, 5)$se)
64+
65+
## basic Structural model
66+
## AirPassengers time series (in logarithms)
67+
y <- log(AirPassengers)
68+
m <- stsm.class::stsm.model(model = "BSM", y = y, transPars = "StructTS")
69+
fit <- StructTS(y, "BSM")
70+
m <- stsm.class::set.pars(m, as.vector(fit$coef[c(4,1:3)]) * 100 / var(y))
71+
ss <- stsm.class::char2numeric(m, P0cov = TRUE)
72+
res <- predict(ss, y, 12)
73+
74+
all.equal(res$pred, predict(fit, 12)$pred)
75+
all.equal(res$se, predict(fit, 12)$se)
76+
77+
# forecasts and confidence intervals for the series
78+
# scaled back to original scale
79+
expy <- exp(y)
80+
plot(cbind(expy, exp(res$pred + 2 * res$se)), type = "n", plot.type = "single")
81+
lines(expy)
82+
lines(exp(res$pred), col = "blue")
83+
lines(exp(res$pred + 2 * res$se), col = "red", lty = 2)
84+
lines(exp(res$pred - 2 * res$se), col = "red", lty = 2)
85+
86+
# forecasts for the trend component
87+
# the aproach in StructTS() seems to seasonal fluctuations in the trend
88+
# see the "stsm" package for a more flexible interface for maximum likelihood
89+
# procedures to fit a structural time series model
90+
trend <- exp(fitted(fit)[,1])
91+
plot(cbind(trend, AirPassengers),
92+
type = "n", plot.type = "single")
93+
lines(AirPassengers, col = "gray")
94+
lines(trend)
95+
lines(exp(res$a[,1]), col = "blue")
96+
lines(exp(res$a[,1] + 2 * sqrt(res$P[,1])), col = "red", lty = 2)
97+
lines(exp(res$a[,1] - 2 * sqrt(res$P[,1])), col = "red", lty = 2)
98+
99+
# forecasts for the seasonal component
100+
seas <- exp(fitted(fit)[,3])
101+
plot(cbind(seas, exp(res$a[,3]) + 2 * sqrt(res$P[,3])),
102+
type = "n", plot.type = "single")
103+
lines(seas)
104+
lines(exp(res$a[,3]), col = "blue")
105+
lines(exp(res$a[,3] + 2 * sqrt(res$P[,3])), col = "red", lty = 2)
106+
lines(exp(res$a[,3] - 2 * sqrt(res$P[,3])), col = "red", lty = 2)
107+
}
108+
109+
\keyword{ts, model}

0 commit comments

Comments
 (0)