Comandos R para análises de séries temporais

9 - Co-integração

Os testes de co-integração desta seção serão realizados com o pacote urca do R.

Teste de Phillips & Ouliaris:

ca.po(X, demean = c("none", "constant", "trend"),lag = c("short", "long"))

X é um objeto do tipo data.frame contendo as duas séries as quais se tem interesse de testar a existência de cointegração.
Defina demean =constant para considerar as séries do modelo com constante (drift). Em outras palavras, utilize esta opção se as séries envolvidas possuírem média constante diferente de zero. A opção demean =trend considera o caso com séries com tendência e drift e o caso demean =none contempla a situação de séries ambas sem tendência e drift.
lag é a opção que define a ordem máxima do modelo VAR ajustado. Se lag='short' então a ordem máxima do VAR ajustado será \(\left[4\left(\frac{T}{100}\right)^{0.25}\right]\). Se lag='long' então teremos a ordem máxima igual a \(\left[12\left(\frac{T}{100}\right)^{0.25}\right]\), em que T indica o comprimento a série.
A forma padrão da função ca.po utiliza a estatística \(P_u\) (Veja Philips e Ouliaris, 1990) para realizar o teste. Desta forma, o teste rejeita H: As séries não são co-integradas para valores grandes da estatistica Pu.
Exemplo(Morettin, 2013):

\[ \begin{alignat}{1} X_{1t} & = & \beta_2 X_{2t} + \beta_3 X_{3t} + u_t \\ X_{2t} & = & X_{2,t-1} + v_t \\ X_{3t} & = & X_{3,t-1} +w_t \\ \underset{\sim}{\beta}^T & = & (1,-\beta_2,-\beta_3)=(1,-0.5,-0.5) \\ u_t & = & 0.75 u_{t-1}+\epsilon_t \nonumber \end{alignat} \]

# Simulando o sistema
beta=c(1,-0.5,-0.5)

u<-arima.sim(n = 250, list(ar = c(0.75),ma=NULL), sd = 0.5)
v<-rnorm(n=250, mean=0, sd=0.5)
w<-rnorm(n=250, mean=0, sd=0.5)

X2<-numeric()
X3<-numeric()
for(t in 1:250){
  if (t==1) {
    X2[t] = 0 + v[t] # Inicializando com X2[0]=0
    X3[t] = 0 + w[t] # Inicializando com X3[0]=0
  }
  else {
    X2[t] = X2[t-1] + v[t]
    X3[t] = X3[t-1] + w[t]
  }
}

X1= -beta[2]*X2-beta[3]*X3+u # Relação de cointegração
X<-cbind(X1,X2,X3) # Cria a matriz de dados

Gráfico das séries simuladas

library(urca)
cointest_X<-ca.po(X, demean="constant") # Teste de P&O para o caso com constante
summary(cointest_X)
## 
## ######################################## 
## # Phillips and Ouliaris Unit Root Test # 
## ######################################## 
## 
## Test of type Pu 
## detrending of series with constant only 
## 
## 
## Call:
## lm(formula = z[, 1] ~ z[, -1])
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.9464 -0.3185  0.0467  0.4428  1.5825 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept)  -0.4136     0.0841   -4.91  1.6e-06
## z[, -1]X2     0.4477     0.0115   38.91  < 2e-16
## z[, -1]X3     0.4821     0.0163   29.61  < 2e-16
## 
## Residual standard error: 0.63 on 247 degrees of freedom
## Multiple R-squared:  0.898,  Adjusted R-squared:  0.898 
## F-statistic: 1.09e+03 on 2 and 247 DF,  p-value: <2e-16
## 
## 
## Value of test-statistic is: 131.04 
## 
## Critical values of Pu are:
##                  10pct   5pct   1pct
## critical values 33.696 40.525 53.873

O valor da estatística é 131.04345 > 53.8731. Logo, rejeitamos H: Não há co-integração ao nível 1%. O vetor de co-integração é: (1, -0.44771, -0.48214) e a estimativa da constante é -0.41356.

É necessário ainda avaliar os resíduos do modelo de correção de erros. O gráfico dos resíduos pode ser obtido por:

plot(cointest_X)

Gráfico dos resíduos do modelo de correção de erros, considerando modelo *cointest_X* ajustado no exemplo

Se houver necessidade de se modelar estes resíduos, ou de se realizar análises adicionais, eles podem ser obtidos através da propriedade res, e armazenados em um outro objeto. Por exemplo:

resid_coint<-cointest_X@res # Armazenando os resíduos do modelo ajustado em 'resid_coint'
Estatísticas descritivas para os resíduos do MCE ajustado no exemplo
resid_coint
Min. :-1.9464
1st Qu.:-0.3185
Median : 0.0467
Mean : 0.0000
3rd Qu.: 0.4428
Max. : 1.5825

Teste de Johansen

Considerando o sistema multivariado dado por: \[ \underset{\sim}{X}_t=\Pi_1 \underset{\sim}{X}_{t-1}+\dots+ \Pi_k \underset{\sim}{X}_{t-k}+\mu+\Phi D_t +\epsilon_t , (t=1,\dots,T) \] O modelo vetorial de correção de erros é especificado por:
\[ \Delta\underset{\sim}{X}_t=\Gamma_1 \Delta\underset{\sim}{X}_{t-1}+\dots+ \Gamma_{k-1} \Delta\underset{\sim}{X}_{t-k+1}+\Gamma_{k}\underset{\sim}{X}_{t-k}+\Phi D_t +\epsilon_t , (t=1,\dots,T)\\ \Gamma_i=-(\underset{\sim}{I}-\Pi_1-\dots-\Pi_i),(i=1,\dots, k) \]
Em que os vetores \(\mu,D_t,\Phi,\epsilon_t,\Pi_i\) e \(\Gamma_i,i=1,\dots,k\) são especificados como em (Johansen,1988). A sintaxe para este teste, utilizando o pacote urca é:

joh_X<-ca.jo(X, type = c("eigen", "trace"), ecdet = c("none", "const", "trend"), K = 2)
summary(joh_X)

Novamente, X é do tipo data.frame. A opção type permite que se escolha qual estatística será utilizada. Utilize type="eigen" para a estatística do maior autovalor e type="trace" para a estatística do traço. A opção ecdet contempla os diversos casos possíveis de adição de constantes e tendências ao VECM. Veja (Zivot e Wang, 2006) para detalhes. K é a ordem do modelo VAR(k) a ser utilizada(Veja equação do sistema proposto acima) e é suposta igual a 2 se suprimida.
O comando summary permite visualizar os detalhes do ajuste, como a martriz de cargas, o(s) vetor(es) co-integrados e as estatísticas de teste. Novamente, rejeitamos a hipótese \(H:r\le r_0\) para valores grandes da estatística, seja ela a do traço ou a do máximo autovalor.