Identify and estimate a transfer function model for the period up until December 2014, using July 2008 as the date of intervention.
rfss <- window(RFSS, end=c(2008,6))
lrfss <- log(rfss)
dlrfss1 <-diff(lrfss)
dlrfss12 <- diff(lrfss, lag=12)
d2lrfss1_12 <- diff(diff(lrfss), lag=12)
par(mfrow=c(2,3))
plot(RFSS, main=expression(RFSS))
plot(lrfss, main=expression(log(rfss)))
plot.new()
plot(dlrfss1, main=expression(paste(Delta, "log(rfss)")))
plot(dlrfss12, main=expression(paste(Delta[12], "log(rfss)")))
plot(d2lrfss1_12, main=expression(paste(Delta, Delta[12], "log(rfss)")))
maxlag <- 20
par(mfrow=c(2,4), max=c(3,3,4,2))
## Warning in par(mfrow = c(2, 4), max = c(3, 3, 4, 2)): "max" is not a
## graphical parameter
Acf(lrfss, type='correlation', lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("ACF for lag(rfss)")))
Acf(dlrfss1, type='correlation', lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("ACF for",Delta," lag(rfss)")))
Acf(dlrfss12, type='correlation',lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("ACF for",Delta[12]," lag(rfss)")))
Acf(d2lrfss1_12, type='correlation',lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("ACF for",Delta, Delta[12],"lag(rfss)")))
Acf(lrfss, type='partial',lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("PACF for lag(rfss)")))
Acf(dlrfss1, type='partial',lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("PACF for",Delta," lag(rfss)")))
Acf(dlrfss12, type='partial',lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("PACF for",Delta[12]," lag(rfss)")))
Acf(d2lrfss1_12, type='partial',lag=maxlag, na.action=na.omit, ylab="", main=expression(paste("PACF for",Delta,Delta[12]," lag(rfss)")))
m1 <- Arima(dlrfss1, order = c(1,1,0), seasonal = list(order = c(1,0,0), period = 12))
m1
## Series: dlrfss1
## ARIMA(1,1,0)(1,0,0)[12]
##
## Coefficients:
## ar1 sar1
## -0.6250 0.9666
## s.e. 0.0557 0.0113
##
## sigma^2 estimated as 0.001257: log likelihood=360.86
## AIC=-715.73 AICc=-715.6 BIC=-705.89
tsdiag(m1, gof.lag=36)
m2 <- Arima(dlrfss1, order = c(1,1,0), seasonal = list(order = c(2,0,0), period = 12))
m2
## Series: dlrfss1
## ARIMA(1,1,0)(2,0,0)[12]
##
## Coefficients:
## ar1 sar1 sar2
## -0.6253 0.8858 0.0841
## s.e. 0.0557 0.0748 0.0767
##
## sigma^2 estimated as 0.001253: log likelihood=361.46
## AIC=-714.92 AICc=-714.71 BIC=-701.81
tsdiag(m2, gof.lag=36)
m3 <- Arima(dlrfss1, order = c(2,0,0), seasonal = list(order = c(1,0,0), period = 12))
m3
## Series: dlrfss1
## ARIMA(2,0,0)(1,0,0)[12] with non-zero mean
##
## Coefficients:
## ar1 ar2 sar1 mean
## -0.6777 -0.4288 0.9706 0.0028
## s.e. 0.0645 0.0642 0.0102 0.0107
##
## sigma^2 estimated as 0.0004412: log likelihood=466.07
## AIC=-922.14 AICc=-921.82 BIC=-905.72
tsdiag(m3, gof.lag=36)
m4 <- Arima(dlrfss1, order = c(2,0,0), seasonal = list(order = c(1,0,0), period = 12))
m4
## Series: dlrfss1
## ARIMA(2,0,0)(1,0,0)[12] with non-zero mean
##
## Coefficients:
## ar1 ar2 sar1 mean
## -0.6777 -0.4288 0.9706 0.0028
## s.e. 0.0645 0.0642 0.0102 0.0107
##
## sigma^2 estimated as 0.0004412: log likelihood=466.07
## AIC=-922.14 AICc=-921.82 BIC=-905.72
tsdiag(m4, gof.lag=36)
m5 <- Arima(dlrfss1, order = c(4,1,0), seasonal = list(order = c(1,0,0), period = 12))
m5
## Series: dlrfss1
## ARIMA(4,1,0)(1,0,0)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 sar1
## -1.3001 -1.2152 -0.6653 -0.2957 0.9717
## s.e. 0.0683 0.1050 0.1044 0.0685 0.0101
##
## sigma^2 estimated as 0.0006354: log likelihood=427.42
## AIC=-842.84 AICc=-842.4 BIC=-823.18
tsdiag(m5, gof.lag=36)
m6 <- Arima(dlrfss1, order = c(2,0,1), seasonal = list(order = c(2,0,0), period = 12))
m6
## Series: dlrfss1
## ARIMA(2,0,1)(2,0,0)[12] with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 sar1 sar2 mean
## -0.3899 -0.2998 -0.4412 0.7088 0.2667 0.0033
## s.e. 0.2078 0.1430 0.2218 0.0810 0.0840 0.0071
##
## sigma^2 estimated as 0.0004048: log likelihood=475.32
## AIC=-936.64 AICc=-936.05 BIC=-913.66
tsdiag(m6, gof.lag=36)
m7 <- Arima(dlrfss1, order = c(3,0,0), seasonal = list(order = c(1,0,0), period = 12))
m7
## Series: dlrfss1
## ARIMA(3,0,0)(1,0,0)[12] with non-zero mean
##
## Coefficients:
## ar1 ar2 ar3 sar1 mean
## -0.7098 -0.4814 -0.0768 0.9692 0.0029
## s.e. 0.0710 0.0809 0.0721 0.0107 0.0096
##
## sigma^2 estimated as 0.0004422: log likelihood=466.63
## AIC=-921.27 AICc=-920.83 BIC=-901.57
tsdiag(m7, gof.lag=36)
m8 <- Arima(dlrfss1, order = c(9,1,0), seasonal = list(order = c(1,0,0), period = 12))
m8
## Series: dlrfss1
## ARIMA(9,1,0)(1,0,0)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7
## -1.6051 -1.9386 -1.9608 -1.9519 -1.6934 -1.3168 -0.9821
## s.e. 0.0706 0.1271 0.1766 0.2062 0.2175 0.2059 0.1773
## ar8 ar9 sar1
## -0.6367 -0.1426 0.9721
## s.e. 0.1289 0.0711 0.0101
##
## sigma^2 estimated as 0.0004297: log likelihood=466.64
## AIC=-911.28 AICc=-909.84 BIC=-875.22
tsdiag(m8, gof.lag=36)
m9 <- Arima(dlrfss1, order = c(9,1,0), seasonal = list(order = c(0,1,0), period = 12))
m9
## Series: dlrfss1
## ARIMA(9,1,0)(0,1,0)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7
## -1.6045 -1.9247 -1.9245 -1.9111 -1.6666 -1.3007 -0.9877
## s.e. 0.0728 0.1308 0.1808 0.2108 0.2220 0.2099 0.1810
## ar8 ar9
## -0.6656 -0.1638
## s.e. 0.1327 0.0745
##
## sigma^2 estimated as 0.0004338: log likelihood=453.37
## AIC=-886.73 AICc=-885.46 BIC=-854.58
tsdiag(m9, gof.lag=36)
m10<- Arima(dlrfss1, order = c(9,1,0), seasonal = list(order = c(2,0,0), period = 12))
m10
## Series: dlrfss1
## ARIMA(9,1,0)(2,0,0)[12]
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7
## -1.6370 -2.0075 -1.9865 -1.9971 -1.7100 -1.3334 -0.9818
## s.e. 0.0707 0.1297 0.1834 0.2138 0.2267 0.2132 0.1843
## ar8 ar9 sar1 sar2
## -0.6449 -0.1156 0.5703 0.4140
## s.e. 0.1330 0.0724 0.0714 0.0719
##
## sigma^2 estimated as 0.0003634: log likelihood=480.98
## AIC=-937.96 AICc=-936.25 BIC=-898.62
tsdiag(m10, gof.lag=36)
rfss <- RFSS
m <- m10
m.f <- forecast(m, h=48)
m.f.err <- lrfss - m.f$mean
## Warning in .cbind.ts(list(e1, e2), c(deparse(substitute(e1))[1L],
## deparse(substitute(e2))[1L]), : non-intersecting series
