PREDICT 413 - Homework 2

The homework assignment here works through the following problems in the Hynman text:

• 7.8 Exercises 1 through 4
• 8.11 Exercises 5 through 8

Chapter 7 - Exponential Smoothing

7.8.1

Consider the pigs series - the number of pigs slaughtered in Victoria each month. * Use the ses function in R to find the optimal values of \(\alpha\) and \(\ell_0\) and generate forecasts for the next four months. * Compute a 95% prediction interval for the first forecast using \(\hat{y} = 1.96s\) where \(s\) is the standard deviation of the residuals. Compare your interval with the interval produced by R.

The intervals produced in R expand as they get farther from the observed data whereas the manually produced intervals are static.

ses.fit.sum <- summary(ses.fit)

Forecast method: Simple exponential smoothing

Model Information:
Simple exponential smoothing 

Call:
 ses(y = pigs) 

  Smoothing parameters:
    alpha = 0.2971 

  Initial states:
    l = 77258.8888 

  sigma:  10253.6

     AIC     AICc      BIC 
4462.955 4463.086 4472.665 

Error measures:
                   ME    RMSE      MAE        MPE     MAPE      MASE       ACF1
Training set 385.9193 10253.6 7961.387 -0.9226055 9.274012 0.7966253 0.01284174

Forecasts:
         Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
Sep 1995       98816.43 85675.91 111956.9 78719.75 118913.1
Oct 1995       98816.43 85108.12 112524.7 77851.38 119781.5
Nov 1995       98816.43 84562.93 113069.9 77017.58 120615.3
Dec 1995       98816.43 84037.84 113595.0 76214.52 121418.3
Jan 1996       98816.43 83530.77 114102.1 75439.03 122193.8
Feb 1996       98816.43 83039.99 114592.9 74688.46 122944.4
Mar 1996       98816.43 82564.03 115068.8 73960.53 123672.3
Apr 1996       98816.43 82101.62 115531.2 73253.33 124379.5
May 1996       98816.43 81651.65 115981.2 72565.17 125067.7
Jun 1996       98816.43 81213.19 116419.7 71894.60 125738.3

7.8.2

Write your own function to implement simple exponential smoothing. The function should take arguments y (the time series), \(\alpha\) (the smoothing parameter) and \(\ell_0\) (the initial level). It should return the forecast of the next observation in the series. Does it give the same forecast as ses?

I didn’t get the same output as the book here. My model was significantly more conservative than the SES model output in the R package.

ts(
  bind_cols(data.frame(pigs), 
            data.frame(ses = yaz.ses(data.frame(pigs), 0.2971, 77258.8888)[[1]],
                       r_model = ses.fit$fitted))
)%>%
  autoplot()+
  theme_yaz()+
  labs(y = 'Number of Pigs Slaughtered',
       title = 'Predicted vs. Actual Pig Slaughtering Rates')+
  scale_color_manual(name = 'Series', values = yaz_cols[3:5])

[[1]]
          Jan      Feb      Mar      Apr      May      Jun      Jul      Aug      Sep      Oct      Nov      Dec
1980 76997.18 75680.73 64368.94 82954.03 85525.73 82749.92 87178.50 84113.61 82272.19 84922.91 81221.05 84448.15
1981 77148.99 78456.83 81532.41 82894.61 84828.14 84093.71 85052.45 83131.99 82601.08 81543.40 84389.32 86856.74
1982 77149.89 79788.43 82592.16 82164.64 83481.09 83393.15 84106.19 82259.11 84811.50 77455.01 82102.25 89084.39
1983 78437.22 80556.14 85849.87 81615.30 85237.54 86974.09 83385.42 85615.45 82969.18 83296.29 86834.75 87043.02
1984 81333.65 83712.53 86150.83 82285.26 88500.59 88084.35 88332.73 89989.95 82115.91 87261.09 84999.27 90011.94
1985 84927.07 85010.86 87381.71 85845.71 87458.66 83850.98 84605.61 79657.71 80144.06 84181.95 80908.50 80825.90
1986 78881.09 77924.13 76542.02 76236.90 77190.59 77611.58 80146.44 76848.63 74974.82 76802.28 73373.75 77288.04
1987 73109.33 71947.07 77577.71 75809.67 70934.26 75028.00 79699.00 75141.79 77813.31 79797.94 78601.52 80109.00
1988 74825.67 77941.36 80501.77 74121.25 75773.12 76125.78 76923.79 79880.24 76719.39 76281.46 77520.37 77671.00
1989 73993.79 76224.12 78083.37 75308.46 78614.89 76777.92 76748.21 78735.51 76690.27 77192.07 77654.06 78004.64
1990 76879.53 77501.06 77356.67 79291.38 83410.97 80762.33 82085.01 79483.61 76448.43 81608.76 78559.03 80983.96
1991 78376.91 77557.80 77262.49 79571.55 82653.96 77944.93 84990.65 82758.54 81429.61 84384.87 88330.65 84353.67
1992 82192.86 82580.58 84069.64 85181.69 84800.51 86507.05 83146.85 81311.96 81144.10 80667.55 79190.07 83797.50
1993 76080.33 77759.24 82351.22 81916.86 81082.90 81357.71 85816.29 85072.95 85225.96 84545.00 82170.58 84702.47
1994 78790.18 79123.22 86692.44 82975.13 84737.23 84907.17 81583.51 83738.08 87006.18 84490.34 83315.90 85479.98
1995 80718.95 81025.26 86012.68 79352.88 88440.87 86020.40 80417.99 84165.61                                    

[[2]]
              Jan          Feb          Mar          Apr          May          Jun          Jul          Aug
1980  -183.957409 -1109.290184 -9060.344142  4003.115971  5810.762214  3859.648668  6972.492349  4818.185667
1981   -77.244466   842.032193  3003.856813  3961.349653  5320.425641  4804.193950  5478.093491  4128.206093
1982   -76.617971  1778.015379  3748.759094  3448.250436  4373.583211  4311.769061  4812.964877  3514.658882
1983   828.249308  2317.636208  6038.597479  3062.120826  5608.195572  6828.816215  4306.339439  5873.829354
1984  2864.148479  4536.263020  6250.143879  3533.036062  7901.792924  7609.219867  7783.803076  8948.665685
1985  5389.966560  5448.857068  7115.333157  6035.673836  7169.420538  4633.578541  5164.010780  1686.129480
1986  1140.243703   467.597152  -503.887405  -718.357448   -48.008044   247.906319  2029.657445  -288.373204
1987 -2916.727595 -3733.676776   224.099518 -1018.657274 -4445.583666 -1568.093187  1715.157071 -1488.110688
1988 -1710.307500   479.709384  2279.420027 -2205.447200 -1044.343560  -796.460462  -235.538811  1842.544341
1989 -2295.035952  -727.337206   579.530884 -1370.956166   953.130599  -338.075122  -358.958281  1037.916224
1990  -266.654718   170.220968    68.728815  1428.640129  4324.298956  2462.565331  3392.283570  1563.754168
1991   785.856495   210.107802     2.529201  1625.568319  3792.196065   482.215363  5434.656520  3865.704785
1992  3468.089437  3740.614662  4787.278591  5568.935233  5301.004303  6500.532956  4138.647673  2848.903773
1993  -828.411695   351.695620  3579.396675  3274.084890  2687.894617  2881.063838  6014.999509  5492.502871
1994  1076.341237  1310.441449  6630.843868  4017.943014  5256.523174  5375.974843  3039.775846  4554.222537
1995  2432.075919  2647.381289  6153.037190  1471.868268  7859.817775  6158.466811  2220.529519  4854.731195
              Sep          Oct          Nov          Dec
1980  3523.847472  5387.042918  2785.001306  5053.330037
1981  3755.024042  3011.583581  5011.981382  6746.327737
1982  5308.731071   137.852071  3404.395802  8312.146999
1983  4013.766382  4243.689962  6730.874199  6877.265144
1984  3414.002055  7030.547531  5440.712636  8964.119223
1985  2027.986793  4866.216932  2565.310474  2507.255292
1986 -1605.474042  -320.950932 -2730.867480    20.488718
1987   389.702969  1784.697990   943.733177  2003.344665
1988  -379.214945  -687.032709   183.795021   289.672637
1989  -399.680441   -46.963886   277.769237   524.190513
1990  -569.669355  3057.526531   913.870260  2618.353697
1991  2931.601082  5008.848908  7782.341255  4986.921591
1992  2730.913925  2395.948054  1357.428557  4595.988855
1993  5600.051140  5121.409135  3452.427068  5232.089878
1994  6851.370027  5082.984123  4257.472847  5778.602149
1995                                                    

7.8.3

I didn’t get the optim() function to work, but here is what I’ve attempted…

yaz.ses(pigs, 0.2971, 77258.8888)[[1]]

Error in formals(sys.function(sys.parent())) : node stack overflow
Error during wrapup: node stack overflow

7.8.4

Cannot do this without the last step working….

Chapter 8: ARIMA Models

8.11,5

The KPSS unit roor test is used for optimizing the amount of differencing used. The appropriate amount of differencing is 6.73.

tseries::kpss.test(ts(retaildata%>%select(A3349414R)))

p-value smaller than printed p-value

    KPSS Test for Level Stationarity

data:  ts(retaildata %>% select(A3349414R))
KPSS Level = 6.7273, Truncation lag parameter = 4, p-value = 0.01

8.11.6

Using the reduced value for \(\theta\) completely smoothed the curve to a flat line.

8.11.7

Fitting and plotting an appropriate ARIMA model for wmurders using the auto.arima() function yields an ARIMA (1,2,1) model.

auto.arima(fpp::wmurders)

Series: fpp::wmurders 
ARIMA(1,2,1)                    

Coefficients:
          ar1      ma1
      -0.2434  -0.8261
s.e.   0.1553   0.1143

sigma^2 estimated as 0.04632:  log likelihood=6.44
AIC=-6.88   AICc=-6.39   BIC=-0.97

The residuals are mostly normally distributed and non-autocorrelated.

checkresiduals(arima.fit)

    Ljung-Box test

data:  Residuals from ARIMA(1,2,1)
Q* = 12.419, df = 8, p-value = 0.1335

Model df: 2.   Total lags used: 10

Plot the forecast and the next 3 values of the time series.

8.11.8

Forecasting USGDP using auto.arima vs. the ETS model yields slightly different outputs. The ETS model continues upward in an almost linear manner while the ARIMA model is slightly more conservative.

autoplot(usgdp)+
  autolayer(forecast(arima.fit, h = 3), series = 'ARIMA')+
  autolayer(forecast(ets(usgdp), h = 3), series = 'ETS')+
  xlim(2004,2007)+
  ylim(7500,12500)

Scale for 'x' is already present. Adding another scale for 'x', which will replace the existing scale.