-
Notifications
You must be signed in to change notification settings - Fork 25
/
add.random.walk.holiday.Rd
166 lines (130 loc) · 4.85 KB
/
add.random.walk.holiday.Rd
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
154
155
156
157
158
159
160
161
162
163
164
165
166
\name{add.random.walk.holiday}
\alias{AddRandomWalkHoliday}
\alias{RandomWalkHolidayStateModel}
\Rdversion{1.1}
\title{
Random Walk Holiday State Model
}
\description{
Adds a random walk holiday state model to the state specification.
This model says
\deqn{%
y_t = \alpha_{d(t), t} + \epsilon_t
}{%
y[t] = alpha[d(t), t] + observation_error,
}
where there is one element in \eqn{\alpha_t}{alpha[, t]} for each day
in the holiday influence window. The transition equation is
\deqn{ %
\alpha_{d(t+1), t+1} = \alpha_{d(t+1), t} + \epsilon_{t+1}
}{ %
alpha[d(t+1), t+1] = alpha[d(t+1), t] + state_error
}
if t+1 occurs on day d(t+1) of the influence window, and
\deqn{
\alpha_{d(t+1), t+1} = \alpha_{d(t+1), t} %
}{%
alpha[d(t+1), t+1] = alpha[d(t+1), t]%
}
otherwise.
}
\usage{
AddRandomWalkHoliday(state.specification = NULL,
y,
holiday,
time0 = NULL,
sigma.prior = NULL,
initial.state.prior = NULL,
sdy = sd(as.numeric(y), na.rm = TRUE))
}
\arguments{
\item{state.specification}{A list of state components that you wish
augment. If omitted, an empty list will be assumed. }
\item{y}{The time series to be modeled, as a numeric vector
convertible to \code{\link[xts]{xts}}. This state model assumes \code{y}
contains daily data.}
\item{holiday}{An object of class \code{\link{Holiday}} describing the
influence window of the holiday being modeled.}
\item{time0}{An object convertible to \code{\link{Date}} containing
the date of the initial observation in the training data. If
omitted and \code{y} is a \code{\link[zoo]{zoo}} or
\code{\link[xts]{xts}} object, then \code{time0} will be obtained
from the index of \code{y[1]}.}
\item{sigma.prior}{An object created by \code{\link[Boom]{SdPrior}}
describing the prior distribution for the standard deviation of the
random walk increments.}
\item{initial.state.prior}{An object created using
\code{\link[Boom]{NormalPrior}}, describing the prior distribution
of the the initial state vector (at time 1).}
\item{sdy}{The standard deviation of the series to be modeled. This
will be ignored if \code{y} is provided, or if all the required
prior distributions are supplied directly. }
}
\value{
A list describing the specification of the random walk holiday state
model, formatted as expected by the underlying C++ code.
}
\references{
Harvey (1990), "Forecasting, structural time series, and the Kalman
filter", Cambridge University Press.
Durbin and Koopman (2001), "Time series analysis by state space
methods", Oxford University Press.
}
\author{
Steven L. Scott \email{[email protected]}
}
\seealso{
\code{\link{bsts}}.
\code{\link{RegressionHolidayStateModel}}
\code{\link{HierarchicalRegressionHolidayStateModel}}
}
\examples{
trend <- cumsum(rnorm(730, 0, .1))
dates <- seq.Date(from = as.Date("2014-01-01"), length = length(trend),
by = "day")
y <- zoo(trend + rnorm(length(trend), 0, .2), dates)
AddHolidayEffect <- function(y, dates, effect) {
## Adds a holiday effect to simulated data.
## Args:
## y: A zoo time series, with Dates for indices.
## dates: The dates of the holidays.
## effect: A vector of holiday effects of odd length. The central effect is
## the main holiday, with a symmetric influence window on either side.
## Returns:
## y, with the holiday effects added.
time <- dates - (length(effect) - 1) / 2
for (i in 1:length(effect)) {
y[time] <- y[time] + effect[i]
time <- time + 1
}
return(y)
}
## Define some holidays.
memorial.day <- NamedHoliday("MemorialDay")
memorial.day.effect <- c(.3, 3, .5)
memorial.day.dates <- as.Date(c("2014-05-26", "2015-05-25"))
y <- AddHolidayEffect(y, memorial.day.dates, memorial.day.effect)
presidents.day <- NamedHoliday("PresidentsDay")
presidents.day.effect <- c(.5, 2, .25)
presidents.day.dates <- as.Date(c("2014-02-17", "2015-02-16"))
y <- AddHolidayEffect(y, presidents.day.dates, presidents.day.effect)
labor.day <- NamedHoliday("LaborDay")
labor.day.effect <- c(1, 2, 1)
labor.day.dates <- as.Date(c("2014-09-01", "2015-09-07"))
y <- AddHolidayEffect(y, labor.day.dates, labor.day.effect)
## The holidays can be in any order.
holiday.list <- list(memorial.day, labor.day, presidents.day)
number.of.holidays <- length(holiday.list)
## In a real example you'd want more than 100 MCMC iterations.
niter <- 100
ss <- AddLocalLevel(list(), y)
ss <- AddRandomWalkHoliday(ss, y, memorial.day)
ss <- AddRandomWalkHoliday(ss, y, labor.day)
ss <- AddRandomWalkHoliday(ss, y, presidents.day)
model <- bsts(y, state.specification = ss, niter = niter, seed = 8675309)
## Plot model components.
plot(model, "comp")
## Plot the effect of the specific state component.
plot(ss[[2]], model)
}
\keyword{models}