9.2 R에서 ARIMA 오차를 고려하는 회귀

R 함수 Arima()에서 xreg 입력값을 사용하면 ARIMA 오차를 고려하여 회귀 모델을 맞춥니다. order 입력값으로 ARIMA 오차 모델의 차수를 정합니다. 차분을 정하면, 모델을 추정하기 전에 회귀 모델의 모든 변수에 대한 차분 계산을 수해합니다. 예를 들어, 다음과 같은 R 명령어는

fit <- Arima(y, xreg=x, order=c(1,1,0))

모델 \(y_t' = \beta_1 x'_t + \eta'_t\)을 맞춥니다. 여기에서 \(\eta'_t = \phi_1 \eta'_{t-1} + \varepsilon_t\)는 AR(1) 오차입니다. 이것은 다음과 같은 모델 \[ y_t = \beta_0 + \beta_1 x_t + \eta_t, \] 과 같습니다. 여기에서 \(\eta_t\)는 ARIMA(1,1,0) 오차입니다. 상수항이 차분 때문에 사라진다는 것에 주목하시길 바랍니다. 차분을 계산한 모델에서 상수를 넣기 위해서는 include.drift=TRUE라고 둡니다.

auto.arima() 함수도 xreg 입력값으로 회귀 항을 다룹니다. 사용하는 사람이 넣을 예측변수(predictor variable)를 반드시 정해야합니다만, auto.arima()는 오차에 대한 가장 좋은 ARIMA 모델을 고릅니다. 차분 계산이 필요한 경우에는, 최종 모델이 원래의 변수에 대한 식으로 표현될 것이지만, 추정 과정에서는 모든 변수에 대해 차분 계산을 수행합니다.

최종 모델에 대한 AICc를 계산하고, 이 값을 가장 좋은 예측변수(predictor variable)를 결정하는데 사용할 수 있습니다. 즉, 이 과정을 다룰 예측변수(predictor variable)의 모든 부분집합에 대해 반복해야만 하고, 최종적으로 AICc 값이 가장 낮은 모델이 선택될 것입니다.

예제: 미국 개인 소득과 소비

그림 9.1은 1970년부터 2016년 3분기까지 분기별 개인 소비 지출과 개인 가처분 소득의 변화를 나타냅니다. 소득의 변화에 따라 지출의 변화를 예측하려고 합니다. 소득의 변화는 곧바로 지출의 변화로 이어질 필요는 없습니다(즉, 직장을 잃은 경우, 새로운 상황에 적응하여 지출을 몇 달 동안 줄일 수도 있습니다). 하지만, 이 예제에서는 이렇게 복잡한 상황을 무시할 것이고 소득의 평균적인 변화의 효과가 소비 지출의 평균적인 변화에 곧바로 영향을 미치는 것을 측정하려고 합니다.

autoplot(uschange[,1:2], facets=TRUE) +
  xlab("연도") + ylab("") +
  ggtitle("미국 소비와 개인 소득의
    분기별 변화")
백분율 변화 미국에서 1970년부터 2016년 3분기까지의 개인 가처분 소득과 개인 소비 지출의 분기별 백분율 변화.

Figure 9.1: 백분율 변화 미국에서 1970년부터 2016년 3분기까지의 개인 가처분 소득과 개인 소비 지출의 분기별 백분율 변화.

(fit <- auto.arima(uschange[,"Consumption"],
  xreg=uschange[,"Income"]))
#> Series: uschange[, "Consumption"] 
#> Regression with ARIMA(1,0,2) errors 
#> 
#> Coefficients:
#>         ar1     ma1    ma2  intercept   xreg
#>       0.692  -0.576  0.198      0.599  0.203
#> s.e.  0.116   0.130  0.076      0.088  0.046
#> 
#> sigma^2 estimated as 0.322:  log likelihood=-156.9
#> AIC=325.9   AICc=326.4   BIC=345.3

데이터에 이미 분명하게 정상성(stationarity)이 나타납니다(원래의 수입과 지출 대신에 백분율 변화를 다루고 있기 때문에), 따라서 어떠한 차분값(difference)도 구할 필요가 없습니다. 맞춘 모델은 다음과 같습니다. \[\begin{align*} y_t &= 0.599 + 0.203 x_t + \eta_t, \\ \eta_t &= 0.692 \eta_{t-1} + \varepsilon_t -0.576 \varepsilon_{t-1} + 0.198 \varepsilon_{t-2},\\ \varepsilon_t &\sim \text{NID}(0,0.322). \end{align*}\]

residuals() 함수로 \(\eta_t\)\(\varepsilon_t\) 시계열 둘 다를 되돌릴 수 있습니다.

cbind("Regression Errors" = residuals(fit, type="regression"),
      "ARIMA errors" = residuals(fit, type="innovation")) %>%
  autoplot(facets=TRUE)
적합 모델에서 얻은 회귀 오차 ($\eta_t$)와 ARIMA 오차 ($\varepsilon_t$).

Figure 9.2: 적합 모델에서 얻은 회귀 오차 (\(\eta_t\))와 ARIMA 오차 (\(\varepsilon_t\)).

이것은 백색잡음(white noise) 시계열과 비슷해야 하는 ARIMA 오차입니다.

checkresiduals(fit)
잔차(즉, ARIMA 오차)가 백색잡음과 크게 다르지 않습니다.

Figure 9.3: 잔차(즉, ARIMA 오차)가 백색잡음과 크게 다르지 않습니다.

#> 
#> 	Ljung-Box test
#> 
#> data:  Residuals from Regression with ARIMA(1,0,2) errors
#> Q* = 5.89, df = 3, p-value = 0.12
#> 
#> Model df: 5.   Total lags used: 8