library(astsa)
## Warning: package 'astsa' was built under R version 3.5.2
Forecast for quarterly U.S GDP
ts.plot(gdp)

diff2 <- diff(diff(log(gdp)))
ts.plot(diff2)

acf2(diff2)

## ACF PACF
## [1,] -0.40 -0.40
## [2,] 0.06 -0.11
## [3,] -0.10 -0.14
## [4,] -0.01 -0.12
## [5,] -0.13 -0.24
## [6,] 0.08 -0.12
## [7,] -0.04 -0.11
## [8,] -0.05 -0.20
## [9,] 0.06 -0.12
## [10,] 0.07 -0.01
## [11,] 0.06 0.08
## [12,] -0.10 -0.06
## [13,] -0.05 -0.15
## [14,] 0.06 0.00
## [15,] -0.12 -0.14
## [16,] 0.11 -0.05
## [17,] -0.04 -0.09
## [18,] 0.07 -0.02
## [19,] -0.03 -0.02
## [20,] 0.11 0.04
## [21,] -0.12 -0.07
## [22,] 0.06 0.01
## [23,] -0.10 -0.06
## [24,] 0.02 -0.07
## [25,] 0.04 0.01
## [26,] -0.04 -0.06
## [27,] 0.02 -0.07
sarima.for(gdp, 12,1,1,1,0,1,0,12)

## $pred
## Qtr1 Qtr2 Qtr3 Qtr4
## 2018 18732.30
## 2019 18824.83 18938.80 19031.71 19114.27
## 2020 19196.09 19329.96 19456.67 19560.18
## 2021 19660.73 19848.55 20008.59
##
## $se
## Qtr1 Qtr2 Qtr3 Qtr4
## 2018 91.04985
## 2019 152.39621 207.23757 256.47665 300.89612
## 2020 341.26148 378.24996 412.42725 444.25278
## 2021 474.09536 502.25023 528.95445
Forecase of crude oil, spot price in dollars per barrel
ts.plot(oil)

diff2 <- diff(diff(log(oil)))
sarima.for(oil, 12,1,1,1,0,1,0,12)

## $pred
## Time Series:
## Start = c(2010, 26)
## End = c(2010, 37)
## Frequency = 52
## [1] 69.79930 71.51837 71.12685 70.37116 72.71470 75.37260 74.04881
## [8] 72.61063 73.92975 69.95017 64.68997 58.85007
##
## $se
## Time Series:
## Start = c(2010, 26)
## End = c(2010, 37)
## Frequency = 52
## [1] 3.535961 5.409860 6.626322 7.716859 8.643360 9.491929 10.265211
## [8] 10.986615 11.662387 12.301590 12.908952 13.489100
Forecast of global mean land temperatures
?gtemp_land
ts.plot(gtemp_land)

diff <- diff(gtemp_land)
acf2(diff)

## ACF PACF
## [1,] -0.43 -0.43
## [2,] 0.02 -0.20
## [3,] -0.18 -0.33
## [4,] 0.18 -0.09
## [5,] -0.08 -0.11
## [6,] 0.00 -0.13
## [7,] -0.02 -0.09
## [8,] 0.12 0.06
## [9,] -0.17 -0.13
## [10,] 0.17 0.09
## [11,] -0.08 0.06
## [12,] 0.04 0.02
## [13,] -0.11 -0.03
## [14,] 0.07 -0.02
## [15,] 0.04 0.05
## [16,] -0.06 -0.05
## [17,] 0.02 0.01
## [18,] 0.08 0.09
## [19,] -0.13 -0.08
## [20,] 0.14 0.10
## [21,] -0.07 0.08
## [22,] 0.10 0.11
sarima.for(gtemp_land, 12,1,1,2,0,1,0,12)

## $pred
## Time Series:
## Start = 2018
## End = 2029
## Frequency = 1
## [1] 1.209286 1.473033 1.145609 1.235889 1.395919 1.245923 1.245923
## [8] 1.335923 1.355923 1.585923 1.765923 1.675923
##
## $se
## Time Series:
## Start = 2018
## End = 2029
## Frequency = 1
## [1] 0.2435807 0.2595829 0.2628782 0.2653672 0.2677569 0.2701172 0.2724563
## [8] 0.2747753 0.2770749 0.2793556 0.2816178 0.2838620
Forecast of global mean ocean temperatures
ts.plot(gtemp_ocean)

diff2 <- diff(gtemp_ocean)
ts.plot(diff2)

acf2(diff2)

## ACF PACF
## [1,] -0.04 -0.04
## [2,] -0.44 -0.44
## [3,] 0.01 -0.04
## [4,] 0.13 -0.08
## [5,] -0.07 -0.09
## [6,] -0.06 -0.05
## [7,] 0.03 -0.04
## [8,] 0.07 0.04
## [9,] -0.11 -0.12
## [10,] 0.03 0.07
## [11,] 0.04 -0.06
## [12,] -0.05 -0.03
## [13,] -0.06 -0.07
## [14,] 0.09 0.04
## [15,] -0.03 -0.09
## [16,] -0.05 -0.02
## [17,] -0.03 -0.09
## [18,] 0.18 0.16
## [19,] 0.05 0.05
## [20,] -0.10 0.06
## [21,] 0.05 0.14
## [22,] 0.00 -0.04
scard2 <- sarima(gtemp_ocean, 1,1,1,1,1,1,12)
## initial value -1.957890
## iter 2 value -2.118555
## iter 3 value -2.149903
## iter 4 value -2.157250
## iter 5 value -2.159665
## iter 6 value -2.161523
## iter 7 value -2.162280
## iter 8 value -2.164726
## iter 9 value -2.169686
## iter 10 value -2.179038
## iter 11 value -2.181247
## iter 12 value -2.182715
## iter 13 value -2.184555
## iter 14 value -2.185867
## iter 15 value -2.186245
## iter 16 value -2.186682
## iter 17 value -2.186839
## iter 18 value -2.186928
## iter 19 value -2.186938
## iter 20 value -2.186939
## iter 21 value -2.186939
## iter 21 value -2.186939
## iter 21 value -2.186939
## final value -2.186939
## converged
## initial value -2.199380
## iter 2 value -2.211360
## iter 3 value -2.234595
## iter 4 value -2.238227
## iter 5 value -2.238459
## iter 6 value -2.238502
## iter 7 value -2.238512
## iter 8 value -2.238515
## iter 9 value -2.238515
## iter 10 value -2.238516
## iter 11 value -2.238516
## iter 12 value -2.238516
## iter 13 value -2.238516
## iter 13 value -2.238516
## iter 13 value -2.238516
## final value -2.238516
## converged

scard2$ttable
## Estimate SE t.value p.value
## ar1 0.5366 0.1260 4.2587 0.0000
## ma1 -0.8171 0.0853 -9.5825 0.0000
## sar1 -0.0321 0.1127 -0.2844 0.7766
## sma1 -0.9541 0.4317 -2.2102 0.0290
sarima.for(gtemp_ocean, 12,1,1,1,0,1,0,12)

## $pred
## Time Series:
## Start = 2018
## End = 2029
## Frequency = 1
## [1] 0.6105146 0.5012268 0.4786098 0.6103564 0.6049575 0.4814258 0.5491155
## [8] 0.5876043 0.6766157 0.8059691 0.8255461 0.7652693
##
## $se
## Time Series:
## Start = 2018
## End = 2029
## Frequency = 1
## [1] 0.1333461 0.1620350 0.1746028 0.1809161 0.1844276 0.1865797 0.1880303
## [8] 0.1890990 0.1899499 0.1906716 0.1913143 0.1919075
Forecast of monthly totals of international airline passengers
ts.plot(log(AirPassengers))

diff2 <- diff(log(AirPassengers))
ts.plot(diff2)

sarima.for(AirPassengers, 12,3,1,1,0,1,0,12)

## $pred
## Jan Feb Mar Apr May Jun Jul
## 1961 443.9169 419.6465 448.7567 491.8838 503.7199 567.4318 654.9918
## Aug Sep Oct Nov Dec
## 1961 639.4503 541.8177 495.1153 424.3551 466.5488
##
## $se
## Jan Feb Mar Apr May Jun Jul
## 1961 11.33788 13.40446 15.36732 16.41562 17.19930 17.72060 18.10506
## Aug Sep Oct Nov Dec
## 1961 18.38579 18.60000 18.76618 18.89893 19.00738
Forecast for U.S monthly live births (adjusted) in thousands
ts.plot(birth)

diff2 <- diff(birth)
ts.plot(diff2)

sarima.for(birth, 12,1,1,1,0,1,0,12)

## $pred
## Jan Feb Mar Apr May Jun Jul
## 1979 261.5340 294.2369 268.8602 280.7508 273.7191 304.7098
## 1980 285.7061
## Aug Sep Oct Nov Dec
## 1979 320.7072 297.7064 290.7062 282.7061 289.7061
## 1980
##
## $se
## Jan Feb Mar Apr May Jun Jul
## 1979 8.351876 9.566134 10.218827 10.733095 11.197741 11.636695
## 1980 13.970426
## Aug Sep Oct Nov Dec
## 1979 12.057685 12.463903 12.857136 13.238649 13.609459
## 1980
Forecast for Johnson & Johnson quarterly earnings per share
ts.plot(jj)

diff2 <- diff(diff(log(jj)))
ts.plot(diff2)

acf2(diff2)

## ACF PACF
## [1,] -0.68 -0.68
## [2,] 0.33 -0.25
## [3,] -0.51 -0.78
## [4,] 0.72 -0.18
## [5,] -0.47 0.16
## [6,] 0.19 0.03
## [7,] -0.33 -0.07
## [8,] 0.53 0.02
## [9,] -0.36 -0.03
## [10,] 0.15 0.09
## [11,] -0.27 0.00
## [12,] 0.43 -0.05
## [13,] -0.28 0.06
## [14,] 0.09 0.01
## [15,] -0.19 0.01
## [16,] 0.31 -0.01
## [17,] -0.18 0.01
## [18,] 0.02 -0.02
## [19,] -0.12 -0.03
## [20,] 0.23 -0.02
sarima.for(jj, 12,1,1,0,0,1,1,12)

## $pred
## Qtr1 Qtr2 Qtr3 Qtr4
## 1981 18.09351 15.75556 17.75742 12.73033
## 1982 19.73216 17.20829 20.41156 14.27748
## 1983 21.58840 19.21297 21.32509 15.97085
##
## $se
## Qtr1 Qtr2 Qtr3 Qtr4
## 1981 0.7497062 0.7638479 0.9913778 1.0202682
## 1982 1.1617757 1.2005040 1.3028862 1.3470855
## 1983 1.4281305 1.4747817 1.5429391 1.5901719