Class Web Page



1. Forecasting with ARMA

This is using Model 1 (direct fit with ARMA) from the last lecture.

  library(forecast)
  source('https://nmimoto.github.io/R/TS-00.txt')

  D = read.csv("https://nmimoto.github.io/datasets/lake.csv")
  Lake  = ts(D[,1], start=1875, freq=1)
  length(Lake)
## [1] 97
  plot(Lake, type="o")

  Fit1 = auto.arima(Lake, d=0, stepwise=FALSE, approximation=FALSE)  
  Fit1
## Series: Lake 
## ARIMA(1,0,1) with non-zero mean 
## 
## Coefficients:
##          ar1     ma1    mean
##       0.7665  0.3393  9.1290
## s.e.  0.0773  0.1123  0.3861
## 
## sigma^2 = 0.4784:  log likelihood = -101.09
## AIC=210.18   AICc=210.62   BIC=220.48
  Randomness.tests(Fit1$residuals)

##   B-L test H0: the series is uncorrelated
##   M-L test H0: the square of the series is uncorrelated
##   J-B test H0: the series came from Normal distribution
##   SD         : Standard Deviation of the series
##       BL15  BL20  BL25  ML15  ML20   JB    SD
## [1,] 0.963 0.952 0.934 0.567 0.641 0.89 0.684


10-step Prediction

  plot(forecast(Fit1, 10))
  abline(h=mean(Lake), col="red")

  forecast(Fit1, 10)
##      Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 1972       9.752287 8.865860 10.63871 8.396614 11.10796
## 1973       9.606775 8.285205 10.92834 7.585609 11.62794
## 1974       9.495234 7.975010 11.01546 7.170252 11.82022
## 1975       9.409735 7.784069 11.03540 6.923494 11.89598
## 1976       9.344197 7.659652 11.02874 6.767908 11.92049
## 1977       9.293960 7.575760 11.01216 6.666200 11.92172
## 1978       9.255451 7.517780 10.99312 6.597913 11.91299
## 1979       9.225933 7.476923 10.97494 6.551053 11.90081
## 1980       9.203307 7.447668 10.95895 6.518289 11.88832
## 1981       9.185963 7.426440 10.94548 6.495006 11.87692


Rolling 1-step Prediciton

  source('https://nmimoto.github.io/R/TS-00.txt')

  Pred = Rolling1step.forecast(Lake,
                               window.size=60,
                               Arima.order=c(1,0,1),
                               include.mean=TRUE)
## 
##   Total length 97 , window size 60 .
##   Last 37 obs retrospectively forecasted with Rolling 1-step
##       prediction using same order and fized window size.
## 
##   Average Prediction Error:  0.0058
##   root Mean Squared Error of Prediction:   0.7303

  # Pred



2. Foreasting with ARMA with Linear Trend

This is using Model 2 (linear trend + ARMA).

  Fit2 = Arima(Lake, order=c(1,0,1), xreg=time(Lake) )
  Fit2
## Series: Lake 
## Regression with ARIMA(1,0,1) errors 
## 
## Coefficients:
##          ar1     ma1  intercept     xreg
##       0.6682  0.3817    55.5443  -0.0242
## s.e.  0.0936  0.1136    18.0324   0.0094
## 
## sigma^2 = 0.4612:  log likelihood = -98.66
## AIC=207.33   AICc=207.98   BIC=220.2


10-step Prediction

  Reg2=lm(Lake~time(Lake))
  h=10
  plot(forecast(Fit2, h, xreg=last(time(Lake))+(1:h)/frequency(Lake)))
  abline(Reg2, col="red")

  forecast(Fit2, h, xreg=last(time(Lake))+(1:h)/frequency(Lake))
##      Point Forecast    Lo 80     Hi 80    Lo 95    Hi 95
## 1972       9.369630 8.499292 10.239968 8.038563 10.70070
## 1973       8.857415 7.595493 10.119337 6.927472 10.78736
## 1974       8.507123 7.105238  9.909007 6.363126 10.65112
## 1975       8.265032 6.804977  9.725086 6.032071 10.49799
## 1976       8.095244 6.609951  9.580538 5.823684 10.36680
## 1977       7.973771 6.477345  9.470197 5.685185 10.26236
## 1978       7.884584 6.383214  9.385954 5.588436 10.18073
## 1979       7.816971 6.313398  9.320543 5.517454 10.11649
## 1980       7.763773 6.259218  9.268329 5.462755 10.06479
## 1981       7.720210 6.215216  9.225203 5.418520 10.02190


Rolling 1-step Prediciton

  source('https://nmimoto.github.io/R/TS-00.txt')

  Pred = Rolling1step.forecast(Lake,
                               window.size=60,
                               Arima.order=c(1,0,1),
                               include.mean=TRUE,     # need this for intercept
                               xreg = TRUE)
## 
##   Total length 97 , window size 60 .
##   Last 37 obs retrospectively forecasted with Rolling 1-step
##       prediction using same order and fized window size.
## 
##   Average Prediction Error:  0.2225
##   root Mean Squared Error of Prediction:   0.7817

  # Pred



3. Forecasting with ARIMA

  Fit3 = Arima(Lake, order=c(1,1,2), include.drift=TRUE)
  Fit3
## Series: Lake 
## ARIMA(1,1,2) with drift 
## 
## Coefficients:
##          ar1      ma1      ma2    drift
##       0.6911  -0.6271  -0.3729  -0.0241
## s.e.  0.0956   0.1251   0.1157   0.0100
## 
## sigma^2 = 0.4663:  log likelihood = -99
## AIC=208.01   AICc=208.67   BIC=220.83


10-step Prediction

  plot(forecast(Fit3, 10))

  forecast(Fit3, 10)
##      Point Forecast    Lo 80     Hi 80    Lo 95    Hi 95
## 1972       9.406472 8.527087 10.285857 8.061569 10.75137
## 1973       8.919379 7.629143 10.209615 6.946133 10.89262
## 1974       8.575314 7.124264 10.026365 6.356124 10.79450
## 1975       8.330089 6.804788  9.855390 5.997342 10.66284
## 1976       8.153168 6.591336  9.714999 5.764552 10.54178
## 1977       8.023447 6.442886  9.604009 5.606187 10.44071
## 1978       7.926345 6.335830  9.516860 5.493862 10.35883
## 1979       7.851783 6.255788  9.447778 5.410920 10.29265
## 1980       7.792799 6.193677  9.391920 5.347154 10.23844
## 1981       7.744578 6.143610  9.345547 5.296109 10.19305


Rolling 1-step Prediction

  source('https://nmimoto.github.io/R/TS-00.txt')

  Pred = Rolling1step.forecast(Lake,
                               window.size=60,
                               Arima.order=c(1,1,2),
                               include.drift=TRUE)
## 
##   Total length 97 , window size 60 .
##   Last 37 obs retrospectively forecasted with Rolling 1-step
##       prediction using same order and fized window size.
## 
##   Average Prediction Error:  0.2124
##   root Mean Squared Error of Prediction:   0.7781

  # Pred



4. Forecasting with ARIMA

  Fit4 = Arima(Lake, order=c(0,1,0))
  Fit4
## Series: Lake 
## ARIMA(0,1,0) 
## 
## sigma^2 = 0.5383:  log likelihood = -106.49
## AIC=214.98   AICc=215.02   BIC=217.54


10-step Prediction

  plot(forecast(Fit4, 10))

  forecast(Fit4, 10)
##      Point Forecast    Lo 80    Hi 80    Lo 95    Hi 95
## 1972           9.96 9.019758 10.90024 8.522024 11.39798
## 1973           9.96 8.630297 11.28970 7.926395 11.99361
## 1974           9.96 8.331453 11.58855 7.469353 12.45065
## 1975           9.96 8.079516 11.84048 7.084048 12.83595
## 1976           9.96 7.857555 12.06244 6.744588 13.17541
## 1977           9.96 7.656887 12.26311 6.437693 13.48231
## 1978           9.96 7.472354 12.44765 6.155473 13.76453
## 1979           9.96 7.300594 12.61941 5.892790 14.02721
## 1980           9.96 7.139274 12.78073 5.646072 14.27393
## 1981           9.96 6.986694 12.93331 5.412721 14.50728


Rolling 1-step Prediction

  source('https://nmimoto.github.io/R/TS-00.txt')

  Pred = Rolling1step.forecast(Lake,
                               window.size=60,
                               Arima.order=c(0,1,0),
                               include.drift=FALSE)
## 
##   Total length 97 , window size 60 .
##   Last 37 obs retrospectively forecasted with Rolling 1-step
##       prediction using same order and fized window size.
## 
##   Average Prediction Error:  0.0843
##   root Mean Squared Error of Prediction:   0.7908

  # Pred



5. Comparing MSE

  # From Model1, 2, 3, estimate for $\sigma^2$
  Sigma_sq = c(0.4784, 0.4612, 0.4663, 0.5383)

  # Estimates of sigma
  sqrt(Sigma_sq)
## [1] 0.6916647 0.6791171 0.6828616 0.7336893
  # 0.6916  0.6791  0.6828  0.7336

  # Rolling 1-step prediction rMSE for the three models were
  # 0.7303, 0.7817, 0.7781, 0.7908