시계열 분석 벡터 자귀회귀

  실세계에서 서로 관련성이 있는 여러 시계열은 병렬로 존재할 수 있다.  이미 이러한 데이터를 정리하고 정렬하는 방법은 살펴 보았다. 이제 정리된 데이터를 최대한 활용하는 방법을 배워볼 차례이다.  그 중 여러 변수에 대한 AR(p) 모델을 생성해 볼 수 있다.  이러한 모델은 변수들에게 서로 영향을 주고 받는 상황을 제공한다는 제공한다는 장점을 가진다. 

즉,  모든 변수가 x  로 지정되며 특별히 y를 위한 변수는 없다는 것이다.  대신 적합은 모든 변수에 대해 균등하게 이루어진다.  한편 계열이 정상이 아닐 때,  전에 살펴본 모델처럼 차분이 적용될 수 있다는 것을 참고해야 한다. 

외생변수와 내생변수
변수가 서로 영향을 주는 모델을 사용할 때 통계에서는 이 변수들을 내생변수 endogenous라고 한다. 모델 내에서 보이는 현상을 통해 이 변수들의 값이 설명 될 수 있다는 의미이다. 외생변수 exogenous 는 이와 대조적인 개념으로, 모델 내에서는 변수들의 설명이 어렵다는 뜻이다. 즉, 가정을 통한 설명이 불가능 하다는 것이다. 따라서 이러한 변수의 값을 수용할 뿐, 값이 왜 그런지 의문을 제가하지 않는다.

그러면 모든 시계열은 자신 뿐만 아니라 다른 시계열도 추정하여 예측하므로 변수 한 개당 방정식 한 개를 갖게 된다.  시계열 세 개가 있다고 가정해보자,  시간 t에서 시계열의 값을 $y_{1,t}, y_{2,t}, y_{3,t}$ 라고 표기 하겠다 그러면 차수 2인 벡터자귀회귀(VAR) 방정식을 다음과 같이 작성 할 수 있습니다. 

VAR 모형 기본 방정식

행렬 곱셈

선형대수에 익숙하다면 이미 눈치 챘을 것이다. 행렬을 사용하면 방정식 세계의 관계를 훨씬 더 간단하게 표현 할 수 있다. 특히 AR과 매우 유사한 방식으로 VAR을 작성 할 수 있다. 상기 세 방정식은 다음과 같은 행렬 형식으로 표현 된다. $ y = ϕ_0 +ϕ_1 × y_{t-1} + ϕ_2 × y_{t_2}$

y 와 $ϕ_0$는 3 × I 행렬이며 ϕ 행렬의 모양은 3 × 3 입니다. 

이 식들을 통해 간단한 경우에도 모델의 파라미터 개수가 빠르게 증가하는 것을 알 수 있다.  예를 들면 지연 p 개와 변수 N 개가 있을 때 각 변수에 대한 예측 방정식은 1+p×N개의 값을 예측 한다.  예측해야 할 변수가 N 개 이므로 이는 $N + p × N^2 $ 총 변수로 변환된다.  즉 변수의 개수가 시계열의 개수에 비례하여 $O(N^2) $로 증가한다.  따라서 근거 없이 단순히 데이터가 더 있다고 해서 시계열을 추가 해서는 안된다.  어떤 관계를 정말로 기대할 수 있을때를 위한 방법으로 남겨 두어야 한다. 

VAR 모델을 가장 자주 사용하는 분야는 계량경제학 econometrics 이다.  그런데 VAR 모델은 모든 값이 서로 영향을 미친다는 가설을 뛰어 넘는,  그 이상의 구조가 없다는 비판을 받곤 한다.  모델의 적합도를 평가하기 어렵게 만들기 때문이다.  그래도 한 변수가 다른 변수를 야기 하는지의 여부를 검정하는 상황에서 VAR 모델은 여전히 유용하게 사용될 수 있다.  또한 변수들의 특정 관계를 주장하기 힘든 도메인 지식이 부족한 분석가가 여러 변수를 예측해야 하는 상황에서도 유용하다. 마지막으로 얼마나 근본적인 '원인'이 특정 값의 예측 변동을 초래하는지를 결정하는 데도  도움이 된다. 

  다음은 간단한 데모이다. UCI 수요 정보를 살펴보고, 두 번째 열로 자체 데이터 대신 은행 주문 banking order을 예측 하는 상황을 고려 한다. (변수를 균등하게 다루는 방식 때문에 해당열도 예측한다.)  그리고 교통 관제 부문 주문 order from the traffic control sector  변수의 사용도 고려 한다.  이름만 봐도 상당히 다른 것이기 때문에 재정 부문 fiscal sector의 과거 주문에 관련된 꽤 독립적인 정보를 제공할 가능성이 있다.  또한 각 열은 경제가 흘러가는 상황과 미래의 수요증가 또는 감소에 대한 근본적인 정보를 제공 한다. 

사용할 파라미터를 결정하기 위해서는 vars 패키지의 VARselect() 메서드를 사용한다. 

# 패키지 로딩 
  library(readxl)
  library(httr)
  library(forecast)

  library(vars)

# 데이터 셋 가져오기 
  url = "https://drive.google.com/u/0/uc?id=1txatfA9L88BFAfDKKXbqfGMSw6TBfwhQ&export=download"
  GET(url, write_disk(tf <- tempfile(fileext = ".xlsx")))
 demand <- read_excel(tf, 1L)

 VARselect(demand[, 11:12, with = FALSE], lag.max = 4, type = "const")  

## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      3      3      1      3 
## 
## $criteria
##                   1            2            3            4
## AIC(n) 3.975854e+01 3.967373e+01 3.957496e+01 3.968281e+01
## HQ(n)  3.984267e+01 3.981395e+01 3.977126e+01 3.993521e+01
## SC(n)  3.997554e+01 4.003540e+01 4.008130e+01 4.033382e+01
## FPE(n) 1.849280e+17 1.700189e+17 1.542863e+17 1.723729e+17

이 결과로 부터  이 함수는 선택 가능한 다양한 정보 기준을 제공한다는 사실을 알 수 있다.  또한 0이 아닌 평균을 수용하기 위해서 "const" 항의 사용을 명시 했다. drift, const 와 drift를 모두 사용하거나 두 개 모두 사용하지 않는 옵션이 있지만,  데이터에 가장 적합한 "const"를 선택 하였다.  다음은 약간 낮은 AIC를 가진 지연 3에 대해 VAR이 동작하는 방식을 살펴 보는 코드이다. 

 ※ 이 서적에서 다루는 데이터셋과  실제 코드는 조금 다르다.   fiscal sector orders  정보를 사용하였지만,  실제 코드에서는 'Banking orders (2)'를 사용하였고,  'Orders from the traffic controller sector' 는 'Banking order (3)'를 사용하였다.   이 쪽 서적은 특히 R 부분의 코드는 정리가 잘되지 않은 듯 하다.   -  그냥 코드를 바꾸고,  변수에 대한 내용은 그대로 기술 하겠다.

## R
  est.var <- VAR(demand[, 11:12, with = FALSE ], p =3, type = "const")
  est.var

## 
## VAR Estimation Results:
## ======================= 
## 
## Estimated coefficients for equation Banking.orders..2.: 
## ======================================================= 
## Call:
## Banking.orders..2. = Banking.orders..2..l1 + Banking.orders..3..l1 + Banking.orders..2..l2 + Banking.orders..3..l2 + Banking.orders..2..l3 + Banking.orders..3..l3 + const 
## 
## Banking.orders..2..l1 Banking.orders..3..l1 Banking.orders..2..l2 
##         -2.134351e-01         -1.070454e+00         -3.332587e-02 
## Banking.orders..3..l2 Banking.orders..2..l3 Banking.orders..3..l3 
##          6.778114e-01         -2.378322e-01          8.626345e-01 
##                 const 
##          1.055717e+05 
## 
## 
## Estimated coefficients for equation Banking.orders..3.: 
## ======================================================= 
## Call:
## Banking.orders..3. = Banking.orders..2..l1 + Banking.orders..3..l1 + Banking.orders..2..l2 + Banking.orders..3..l2 + Banking.orders..2..l3 + Banking.orders..3..l3 + const 
## 
## Banking.orders..2..l1 Banking.orders..3..l1 Banking.orders..2..l2 
##          9.985195e-02          4.954869e-01          3.819053e-02 
## Banking.orders..3..l2 Banking.orders..2..l3 Banking.orders..3..l3 
##          4.288644e-01         -4.107863e-03         -2.681979e-01 
##                 const 
##         -2.320958e+03

 par(mfrow = c(2, 1))
  plot(demand$`Banking orders (2)`, type = "l",
       ylab = "Fiscal sector orders")
  lines(fitted(est.var)[, 1], col =2, lty=2)

    plot(demand$`Banking orders (3)`, 
       type = "l", ylab = "Orders frome the traffic controller sector")
  lines(fitted(est.var)[, 2], col = 2, lty =2)

상단과 하단의 그래프는 각각 재정 부분의 주문과 교통관제 부문의 주문에 대한 것이다.  실선과 파선은 실제값과 예측값의 변화를 보여준다. 흥미롭게도 상단 그래프는 예측이 실제 데이터와 비교해서 다소 '느리게' 변화하여 전형적인 예측에 가까운 반면, 하단 그래프는 실제 변화가 발생하기 전에 미리 예측이 이루어진 것으로 보인다.   이것은 fisical sector orders가order from the traffic controller sector를 이끄는 상황이다. 

par(mfrow = c(2,1))
  acf(demand$`Banking orders (2)` -  fitted(est.var)[,1], main = "ACF of Prediction Errors\n Fiscal secoter orders")

## Warning in demand$`Banking orders (2)` - fitted(est.var)[, 1]: 두 객체의 길이가
## 서로 배수관계에 있지 않습니다

  acf(demand$`Banking orders (3)` -  fitted(est.var)[,2], main = "ACF of Prediction Error\n Orders frome the traffic controller sector")

## Warning in demand$`Banking orders (3)` - fitted(est.var)[, 2]: 두 객체의 길이가
## 서로 배수관계에 있지 않습니다

 

각 시계열 잔차에 대해 자기상관함수 그래프 이다.  두 계열 모두 모델로 완전한 설명이 어려울수 있는 지연 3의 오차에 대한 임계치를 넘는 자기 상관이 있다. 

ACF 는 우리가 원하는 정도로 잔차에 자기상관이 없음을 명확히 뒷받침하지 않는다.  이때 순차적 상관관계를 위해서 포트맨트 검정 portmanteau test을 적용해 본다.  이 검정은 vars 패키지의 serial.test() 메서드로 접근할 수 있다. 이 검정은 단변량 사례에서 살펴본 순차적 상관관계 검정과 유사하다. 

 serial.test(est.var, lags.pt = 8, type ="PT.asymptotic")

## 
##  Portmanteau Test (asymptotic)
## 
## data:  Residuals of VAR object est.var
## Chi-squared = 20.463, df = 20, p-value = 0.4293

P 값이 매우 크기 때문에 잔차에  연속적인 상관관계가 없다는 귀무가설을 기각 할 수 없다. 

즉  이 모델은 연속적인 상관관계가 없다는 것이다. 

ARMA 및  ARIMA까지 이르는 다양한 단변량 모델을 살펴 보았다.  이쯤 되면 VARIMA 모델이 있는지 궁금한 사람이 있을지 모른다.  실제로 VARIMA 모델이 존재하지만, 이미 꽤 복잡한 VAR의 성능이 상대적으로 괜찮기 때문에 VARIMA 모델을 많이 사용하지 않는다. 산업과 학계 모두에서 사용한 사례를 살펴보면 VARIMA 보다 VAR를 압도적으로 많이 사용한다는 사실을 알 수 있다. 

# 참고서적 :  실전 시계열 분석 - 한빛 미디어