STLF function in the FORECAST package

Question

I am trying to forecast a yearly time series on a weekly bases (52 weeks a year and I have 164 weeks data). As the frequency is larger than 24, R advices me to use "stlf" rather than "ets" to avoid seasonality being ignored. The "stlf" function works perfectly well and I got the following:

> WR.ets<-stlf(WeeklyReferral,method="ets")
> summary(WR.ets)

Forecast method: STL +  ETS(A,A,N)

Model Information:
ETS(A,A,N) 

Call:
 ets(y = x.sa, model = etsmodel) 

  Smoothing parameters:
    alpha = 0.0262 
    beta  = 1e-04 

  Initial states:
    l = 93.1548 
    b = 0.1159 

  sigma:  12.6201

     AIC     AICc      BIC 
1675.954 1676.205 1688.353 

Error measures:
                     ME     RMSE      MAE       MPE     MAPE      MASE
Training set -0.1869514 12.62011 9.790321 -2.589141 11.12905 0.5990874

Forecasts:
         Point Forecast     Lo 80     Hi 80     Lo 95     Hi 95
2013.423       95.39869  79.22537 111.57201  70.66373 120.13364
2013.442       95.03434  78.85538 111.21330  70.29075 119.77793
...............................................................

The point forecast gives the mean of the predicted value. However, what I want is the actual forecast value rather than the mean. Thus I am trying to understand how it works and break down the steps. I use "stl" decomposition firstly on the time series

temp<-stl(WeeklyReferral,s.window="periodic", robust=TRUE)
> temp
 Call:
 stl(x = WeeklyReferral, s.window = "periodic", robust = TRUE)

Components
Time Series:
Start = c(2010, 15) 
End = c(2013, 22) 
Frequency = 52 
            seasonal     trend   remainder
2010.269   7.1597729  82.33453  -0.4943046
2010.288  -1.4283001  82.69446   5.7338358
..........................................
2013.404   8.0046803 117.74388  -0.7485615

Then I use "trend+remainder" as the new time series to forecast for 3 months (12 periods). I use the last state vector obtained by "stlf" function as the initial state vector in my following formulas. And add the seasonal values at the same week last year back to the forecasted values as the "stlf" function shows the model is ETS(A,A,N).

y<-c(rep(NA,13))
l<-c(rep(NA,13))
b<-c(rep(NA,13))
e<-c(rep(NA,12))
alpha<-0.0262
beta<-0.0001

y[1]<-117.74388-0.7485615
l[1]<-109.66913
b[1]<-0.11284923

for (j in 1:1000){
  for(i in 2:13){
e[i-1]=rnorm(sd=12.6201,n=1)
b[i]<-b[i-1]+beta*e[i-1]
l[i]<-l[i-1]+b[i-1]+alpha*e[i-1]
y[i]<-l[i-1]+b[i-1]+e[i-1]+temp$time.series[i+164-52,1]
}}

Am I right?

I tried to use "ets" function on the new decomposed time series and it gave different parameters (alpha, beta, l,b, sigma) and it didn't give any forecasted values.

Any opinions are appreciated.

Solution

As far as I can tell from the comments above, you actually want to simulate future sample paths from the model rather than obtain point forecasts or interval forecasts. The following code will do it.

# STL decomposition
temp <- stl(WeeklyReferral, s.window="periodic", robust=TRUE)
# Seasonally adjusted data
sa <- seasadj(temp)
seascomp <- tail(temp$time.series,52)[,1]
# ETS model
fit <- ets(sa, "ZZN")
# Simulations from ETS model with re-seasonalization
sim <- matrix(0, nrow=52, ncol=1000)
for(i in 1:1000)
  sim[,i] <- simulate(fit, nsim=52) + seascomp

The matrix sim contains 1000 future sample paths each of length 52.

Comments

[MichelleX]

Thank you so much Professor!

[MichelleX]

One more question, Professor. When I apply "stlf" function, it shows that "Forecast method: STL + ETS(A,N,N)". However, when I apply STL decomposition manually and apply ETS, the fitted model I got from ETS is "ETS(A,A,N)". Any possible reason for this? And also how do I know whether the seasonal components need to be added or multiplied to the forecast value. In other words, when I apply the STL decomposition, how to test if the model is additive or multiplicative?