-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
157 lines (120 loc) · 3.58 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
---
output: github_document
---
<!-- README.md is generated from README.Rmd. Please edit that file -->
```{r, include = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "man/figures/README-",
out.width = "100%"
)
```
# wex
<!-- badges: start -->
<!-- badges: end -->
`wex` is an R package designed to compute the exact observation weights for the Kalman filter and smoother using the method described in Koopman and Harvey (2003). Built on top of the FKF package, `wex` enhances the functionality of the existing packages and allows to get further insights from the state space models, as is illustrated in a number of motivating examples below.
## Installation
You can install the development version of wex from [GitHub](https://github.com/) with:
``` r
# install.packages("devtools")
devtools::install_github("timginker/wex")
```
## Example 1: Local level model
In this example, we fit the local level model to the Nile dataset and compute the associated smoothed and filtered values.
The resulting estimates are presented in the plot below:
```{r echo=F, eval=T}
library(FKF)
library(wex)
data(Nile)
y <- Nile
y[c(3, 10)] <- NA # NA values can be handled
# Computing smoothed and filtered values of the local level
fit_kf <- fkf(a0 = y[1],
P0 = matrix(100),
dt = matrix(0),
ct = matrix(0),
Tt = matrix(1),
Zt = matrix(1),
HHt = matrix(1385.066),
GGt = matrix(15124.13),
yt = t(y))
fit_ks <- fks(fit_kf)
# Filtered values of the local level
mu_t <- fit_kf$att[1,]
mu_t <- ts(mu_t)
tsp(mu_t) <- tsp(y)
# Smoothed values of the local level
mu_T <- fit_ks$ahatt[1,]
mu_T <- ts(mu_T)
tsp(mu_T) <- tsp(y)
```
```{r eval=T, echo=F}
plot.ts(
y,
col = "darkgrey",
xlab = "",
ylab = "",
lwd = 1.5,
main = "Flow of the River Nile and its Estimated Smoothed and Filtered Levels"
)
lines(mu_t, col = "black")
lines(mu_T, col = "blue")
legend(
"bottomleft",
legend = c("Flow", "Filtered Level", "Smoothed Level"),
lwd = c(2, 1),
col = c("darkgrey", "black", "blue"),
bty = "n"
)
```
Now, w.l.o.g., let's consider the 50th value of the estimated local level
```{r}
cat("smoothed level[50] = ",mu_T[50])
```
It is computed as (insert formula here).
Now, We can compute the weight of each observation using the `wex` function, and compare the local level estimates obtained from the weighted average of the observed data with the associated estimates obtained from the Kalman filter and smoother.
```{r}
wts=wex(Tt=matrix(1),
Zt=matrix(1),
HHt = matrix(1385.066),
GGt = matrix(15124.13),
yt = t(y),
t=50)
```
We can also visualize the weights assigned to each observation:
```{r}
par(mfrow = c(2, 1),
mar = c(2.2, 2.2, 1, 1),
cex = 0.8)
plot(
wts$Wt,
col = "darkgrey",
xlab = "",
ylab = "",
lwd = 1.5,
type="l",
main="Filtering weights"
)
plot(
wts$WtT,
col = "blue",
xlab = "",
ylab = " ",
lwd = 1.5,
type="l",
main="Smoothing weights"
)
```
It is also easy to verify the identity between the smoothed and filtered levels obtained from the Kalman filter and the corresponding estimates computed using the weights.
```{r echo=F}
cat("\n Smoothed level computed using the weights = ",
sum(y*as.numeric(wts$WtT),na.rm = T),
" \n Smoothed level from the Kalman Filter = ",fit_ks$ahatt[50])
```
```{r echo=F}
cat("\n Filtered level computed using the weights = ",
sum(y*as.numeric(wts$Wt),na.rm = T),
" \n Filtered level from the Kalman Filter = ",fit_kf$att[50])
```
## Example 2: