We are using a build in data set from package “forecast”. The original data is from Jan 1956 to Aug 1995.
library(forecast)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
gas
## Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 1956 1709 1646 1794 1878 2173 2321 2468 2416 2184 2121 1962 1825
## 1957 1751 1688 1920 1941 2311 2279 2638 2448 2279 2163 1941 1878
## 1958 1773 1688 1783 1984 2290 2511 2712 2522 2342 2195 1931 1910
## 1959 1730 1688 1899 1994 2342 2553 2712 2627 2363 2311 2026 1910
## 1960 1762 1815 2005 2089 2617 2828 2965 2891 2532 2363 2216 2026
## 1961 1804 1773 2015 2089 2627 2712 3007 2880 2490 2237 2205 1984
## 1962 1868 1815 2047 2142 2743 2775 3028 2965 2501 2501 2131 2015
## 1963 1910 1868 2121 2268 2690 2933 3218 3028 2659 2406 2258 2057
## 1964 1889 1984 2110 2311 2785 3039 3229 3070 2659 2543 2237 2142
## 1965 1962 1910 2216 2437 2817 3123 3345 3112 2659 2469 2332 2110
## 1966 1910 1941 2216 2342 2923 3229 3513 3355 2849 2680 2395 2205
## 1967 1994 1952 2290 2395 2965 3239 3608 3524 3018 2648 2363 2247
## 1968 1994 1941 2258 2332 3323 3608 3957 3672 3155 2933 2585 2384
## 1969 2057 2100 2458 2638 3292 3724 4652 4379 4231 3756 3429 3461
## 1970 3345 4220 4874 5064 5951 6774 7997 7523 7438 6879 6489 6288
## 1971 5919 6183 6594 6489 8040 9715 9714 9756 8595 7861 7753 8154
## 1972 7778 7402 8903 9742 11372 12741 13733 13691 12239 12502 11241 10829
## 1973 11569 10397 12493 11962 13974 14945 16805 16587 14225 14157 13016 12253
## 1974 11704 12275 13695 14082 16555 17339 17777 17592 16194 15336 14208 13116
## 1975 12354 12682 14141 14989 16159 18276 19157 18737 17109 17094 15418 14312
## 1976 13260 14990 15975 16770 19819 20983 22001 22337 20750 19969 17293 16498
## 1977 15117 16058 18137 18471 21398 23854 26025 25479 22804 19619 19627 18488
## 1978 17243 18284 20226 20903 23768 26323 28038 26776 22886 22813 22404 19795
## 1979 18839 18892 20823 22212 25076 26884 30611 30228 26762 25885 23328 21930
## 1980 21433 22369 24503 25905 30605 34984 37060 34502 31793 29275 28305 25248
## 1981 27730 27424 32684 31366 37459 41060 43558 42398 33827 34962 33480 32445
## 1982 30715 30400 31451 31306 40592 44133 47387 41310 37913 34355 34607 28729
## 1983 26138 30745 35018 34549 40980 42869 45022 40387 38180 38608 35308 30234
## 1984 28801 33034 35294 33181 40797 42355 46098 42430 41851 39331 37328 34514
## 1985 32494 33308 36805 34221 41020 44350 46173 44435 40943 39269 35901 32142
## 1986 31239 32261 34951 38109 43168 45547 49568 45387 41805 41281 36068 34879
## 1987 32791 34206 39128 40249 43519 46137 56709 52306 49397 45500 39857 37958
## 1988 35567 37696 42319 39137 47062 50610 54457 54435 48516 43225 42155 39995
## 1989 37541 37277 41778 41666 49616 57793 61884 62400 50820 51116 45731 42528
## 1990 40459 40295 44147 42697 52561 56572 56858 58363 45627 45622 41304 36016
## 1991 35592 35677 39864 41761 50380 49129 55066 55671 49058 44503 42145 38698
## 1992 38963 38690 39792 42545 50145 58164 59035 59408 55988 47321 42269 39606
## 1993 37059 37963 31043 41712 50366 56977 56807 54634 51367 48073 46251 43736
## 1994 39975 40478 46895 46147 55011 57799 62450 63896 57784 53231 50354 38410
## 1995 41600 41471 46287 49013 56624 61739 66600 60054
autoplot(gas, main = "Australia monthly gas production", ylab = "Gas production")
We can use window() function to extra data between Jan 1990 to Dec 1994.
gas.data <- window(gas, start = c(1990,1),end = c(1994,12))
autoplot(gas.data, main = "Australia monthly gas production", ylab = "Gas production")
fit <- auto.arima(gas.data)
fit
## Series: gas.data
## ARIMA(0,1,1)(1,1,0)[12]
##
## Coefficients:
## ma1 sar1
## -0.6818 -0.5791
## s.e. 0.1322 0.1385
##
## sigma^2 = 11940296: log likelihood = -451.38
## AIC=908.76 AICc=909.32 BIC=914.31
AIC and BIC of the fitted model.
AIC(fit)
## [1] 908.7579
BIC(fit)
## [1] 914.3083
Make 1-step ahead prediction for Jan 1996.
Jan1996_frc <- forecast(fit, h = 1)
Jan1996_frc
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## Jan 1995 40923.67 36495.3 45352.04 34151.07 47696.27
Make prediction for all months in 1996
frc_1996 <- forecast(fit, h =12)
frc_1996
## Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
## Jan 1995 40923.67 36495.30 45352.04 34151.07 47696.27
## Feb 1995 41658.90 37011.79 46306.00 34551.76 48766.03
## Mar 1995 40352.25 35496.24 45208.25 32925.63 47778.86
## Apr 1995 46216.00 41159.72 51272.27 38483.09 53948.90
## May 1995 54958.38 49709.47 60207.29 46930.86 62985.90
## Jun 1995 59960.34 54525.61 65395.06 51648.64 68272.03
## Jul 1995 61819.43 56205.03 67433.82 53232.95 70405.90
## Aug 1995 61169.61 55381.12 66958.10 52316.88 70022.34
## Sep 1995 56705.19 50747.70 62662.68 47593.99 65816.39
## Oct 1995 52881.30 46759.46 59003.13 43518.75 62243.84
## Nov 1995 50615.26 44333.38 56897.15 41007.95 60222.57
## Dec 1995 44131.73 37693.78 50569.68 34285.73 53977.73
We can also make a time series plot with point estimation and 95% prediction interval.
autoplot(frc_1996)