#COINTEGRACION Y CAUSALIDAD

#Prueba de Dickey y Fuller
#Definicion: El contraste de Dickey-Fuller es una prueba de raíz única que detecta estadísticamente la presencia de conducta tendencial estocástica en las series temporales de las variables mediante un contraste de hipótesis. 

#En otras palabras, el contraste de Dickey-Fuller nos permite saber si hay presencia significativa de tendencia en las series temporales de las variables mediante un contraste de hipótesis. 
#¿Donde lo aplicamos? El contraste Dickey-Fuller es comúnmente aplicado en econometría para comprobar la presencia de tendencia sobre las series temporales. La particularidad del contraste Dickey-Fuller es que es la herramienta más fácil de usar comparado con otros contrastes más complejos que también prueban la presencia de tendencia en los datos. 
#Hipotesis de la prueba(Nula y Alternativa)
#El punto de partida es el proceso estocástico de raíz unitaria se inicia con: 
#                    Yt = ρYt−1 + ut −1 ≤ ρ ≤ 1
#Restando en ambos lados de la ecuación Yt-1
#                   Yt - Yt-1 = pYt-1 - Yt-1 + ut-1
#                   = (ρ − 1)Yt−1 + ut
#La cual se expresa como: DeltaYt = δYt−1 + ut
#Donde δ = (ρ − 1) y Delta, es el operador de primeras diferencias.
#Por tanto: DeltaYt = (Yt − Yt−1) = ut (Econometria-DamodarN.Gujarati)
#Definicion de hipotesis nula; Si el parámetro alfa es igual a cero implica que la serie sigue un camino aleatorio, en consecuencia es no estacionaria; Si el parámetro alfa es distinto de cero y es negativo, entonces la serie es estacionaria
#En el contexto de la prueba Dickey-Fuller hay dos condiciones para que la serie de tiempo sea estacionaria:Rechazar la hipótesis nula y Que el estimador alfa sea negativo
#Es decir: Hipótesis nula: H0:δ = 0 (es decir, existe una raíz unitaria, la serie de tiempo es no estacionaria o tiene tendencia   estocástica).
#Es decir: Hipótesis alternativa: H1:δ < 0 (es decir, la serie de tiempo es estacionaria, posiblemente alrededor de una tendencia determinista).
#Esto se pude probar mediante un estadístico t de Student

#Consideraciones en la prueba Dickey-Fuller
# 1-Distribución del estadístico bajo la hipótesis nula, tiende a valores negativos
# 2-Los resultados de la estimación del parámetro α son afectados por la presencia de autocorrelación en los errores de la prueba
# 3-Los componentes de constante y tendencia deben ser incorporados en la especificación de la prueba

#EJEMPLO EN R(Estadistico de prueba)

#Cargamos librerias necesarias para series de tiempo 
##(Es opcional instalar todas estas librerias)##
library(lubridate)
## 
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union
library(tseries)
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.0     v dplyr   1.0.7
## v tidyr   1.1.4     v stringr 1.4.0
## v readr   2.0.2     v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x lubridate::as.difftime() masks base::as.difftime()
## x lubridate::date()        masks base::date()
## x dplyr::filter()          masks stats::filter()
## x lubridate::intersect()   masks base::intersect()
## x dplyr::lag()             masks stats::lag()
## x lubridate::setdiff()     masks base::setdiff()
## x lubridate::union()       masks base::union()
library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## The following object is masked from 'package:purrr':
## 
##     some
library(astsa)
library(foreign)
library(timsac)
library(vars)
## Loading required package: MASS
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
## Loading required package: strucchange
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Loading required package: sandwich
## 
## Attaching package: 'strucchange'
## The following object is masked from 'package:stringr':
## 
##     boundary
## Loading required package: urca
## Loading required package: lmtest
library(lmtest)
library(mFilter)
library(dynlm)
library(nlme)
## 
## Attaching package: 'nlme'
## The following object is masked from 'package:dplyr':
## 
##     collapse
library(broom)
library(kableExtra)
## 
## Attaching package: 'kableExtra'
## The following object is masked from 'package:dplyr':
## 
##     group_rows
library(knitr)
library(MASS)
library(parallel)
library(car)
library(mlogit)
## Loading required package: dfidx
## 
## Attaching package: 'dfidx'
## The following object is masked from 'package:MASS':
## 
##     select
## The following object is masked from 'package:stats':
## 
##     filter
library(dplyr)
library(tidyr)
library(forecast)
## 
## Attaching package: 'forecast'
## The following object is masked from 'package:nlme':
## 
##     getResponse
## The following object is masked from 'package:astsa':
## 
##     gas
library(fpp2)
## -- Attaching packages ---------------------------------------------- fpp2 2.4 --
## v fma       2.4     v expsmooth 2.3
## -- Conflicts ------------------------------------------------- fpp2_conflicts --
## x forecast::getResponse() masks nlme::getResponse()
## x car::some()             masks purrr::some()
## 
## Attaching package: 'fpp2'
## The following object is masked from 'package:astsa':
## 
##     oil
library(stats)
library(quantmod)
## Loading required package: xts
## 
## Attaching package: 'xts'
## The following objects are masked from 'package:dplyr':
## 
##     first, last
## Loading required package: TTR

#CARGAMOS BASE DATOS

library(readxl)
coint <- read_excel("C:/Users/naye/Downloads/coint.xls")
View(coint)

#LLAMAMOS EL ARCHIVO PARA QUE ESTE SEA LEIDO

#Leer Documento
attach(coint)
names(coint)
## [1] "DPI" "PCE"

#ANTES DE CONVERTIRLO A UN OBJETO DE SERIE DE TIEMPO SE PUEDEN CONVERTIR A LOGARITMO, DEPENDIENDO DE LOS VALORES ELEVADOS(NO SIEMPRE APLICA)

#Generar Logaritmos
lnPCE=log(PCE)
lnDPI=log(DPI)
#Generar variables de Series de Tiempo
DPI.ts=ts(lnDPI, start=c(1947,1),frequency = 4)
PCE.ts=ts(lnPCE, start=c(1947,1),frequency = 4)

#USAMOS UN DATAFRAME PARA LOS DATOS

#Hacer una tabla
datos1=cbind(DPI.ts,PCE.ts)
datos1
##           DPI.ts   PCE.ts
## 1947 Q1 6.999422 6.924809
## 1947 Q2 6.978027 6.941190
## 1947 Q3 7.005608 6.944569
## 1947 Q4 6.993658 6.944762
## 1948 Q1 7.009680 6.949473
## 1948 Q2 7.043422 6.960632
## 1948 Q3 7.063391 6.962338
## 1948 Q4 7.066382 6.970542
## 1949 Q1 7.045428 6.971762
## 1949 Q2 7.048734 6.987121
## 1949 Q3 7.054536 6.989335
## 1949 Q4 7.061077 7.003247
## 1950 Q1 7.133136 7.020101
## 1950 Q2 7.127212 7.035973
## 1950 Q3 7.142669 7.086152
## 1950 Q4 7.152582 7.056261
## 1951 Q1 7.152191 7.079521
## 1951 Q2 7.168195 7.051423
## 1951 Q3 7.174648 7.062963
## 1951 Q4 7.176637 7.068938
## 1952 Q1 7.176331 7.071658
## 1952 Q2 7.188337 7.090660
## 1952 Q3 7.207638 7.095064
## 1952 Q4 7.226936 7.129538
## 1953 Q1 7.240650 7.141562
## 1953 Q2 7.254531 7.147717
## 1953 Q3 7.250423 7.145354
## 1953 Q4 7.244085 7.138629
## 1954 Q1 7.249570 7.142274
## 1954 Q2 7.249499 7.154693
## 1954 Q3 7.260312 7.167886
## 1954 Q4 7.279733 7.188413
## 1955 Q1 7.294309 7.210449
## 1955 Q2 7.317810 7.229186
## 1955 Q3 7.336481 7.241438
## 1955 Q4 7.347751 7.253683
## 1956 Q1 7.357811 7.255238
## 1956 Q2 7.367645 7.258553
## 1956 Q3 7.371867 7.260804
## 1956 Q4 7.387647 7.274341
## 1957 Q1 7.389502 7.281179
## 1957 Q2 7.396028 7.282830
## 1957 Q3 7.400926 7.290975
## 1957 Q4 7.395230 7.291520
## 1958 Q1 7.385975 7.277455
## 1958 Q2 7.392463 7.285438
## 1958 Q3 7.412643 7.301283
## 1958 Q4 7.424821 7.313820
## 1959 Q1 7.431655 7.330340
## 1959 Q2 7.450661 7.347106
## 1959 Q3 7.449789 7.358321
## 1959 Q4 7.457725 7.359722
## 1960 Q1 7.469198 7.368718
## 1960 Q2 7.474091 7.380941
## 1960 Q3 7.474659 7.376821
## 1960 Q4 7.473751 7.377946
## 1961 Q1 7.483019 7.377884
## 1961 Q2 7.498094 7.392771
## 1961 Q3 7.511634 7.397623
## 1961 Q4 7.531231 7.417520
## 1962 Q1 7.540834 7.428155
## 1962 Q2 7.551765 7.440205
## 1962 Q3 7.557316 7.448334
## 1962 Q4 7.565483 7.462502
## 1963 Q1 7.573531 7.469141
## 1963 Q2 7.582891 7.478735
## 1963 Q3 7.594079 7.492203
## 1963 Q4 7.610655 7.500695
## 1964 Q1 7.630752 7.520343
## 1964 Q2 7.657661 7.537909
## 1964 Q3 7.670708 7.556219
## 1964 Q4 7.682206 7.559038
## 1965 Q1 7.693163 7.580853
## 1965 Q2 7.704091 7.591862
## 1965 Q3 7.731229 7.608871
## 1965 Q4 7.751174 7.636607
## 1966 Q1 7.760424 7.651215
## 1966 Q2 7.762894 7.653969
## 1966 Q3 7.775402 7.665332
## 1966 Q4 7.786801 7.669495
## 1967 Q1 7.801800 7.675360
## 1967 Q2 7.809744 7.688639
## 1967 Q3 7.819275 7.693710
## 1967 Q4 7.826483 7.699752
## 1968 Q1 7.843770 7.723253
## 1968 Q2 7.860301 7.738532
## 1968 Q3 7.862151 7.757137
## 1968 Q4 7.868522 7.761575
## 1969 Q1 7.869975 7.772753
## 1969 Q2 7.879859 7.779049
## 1969 Q3 7.899747 7.783641
## 1969 Q4 7.906952 7.791440
## 1970 Q1 7.911837 7.797456
## 1970 Q2 7.929271 7.802087
## 1970 Q3 7.942575 7.810799
## 1970 Q4 7.938945 7.807957
## 1971 Q1 7.959835 7.827002
## 1971 Q2 7.974051 7.836172
## 1971 Q3 7.978105 7.844123
## 1971 Q4 7.988475 7.860648
## 1972 Q1 7.994632 7.873826
## 1972 Q2 8.003530 7.892564
## 1972 Q3 8.024077 7.907799
## 1972 Q4 8.063063 7.930818
## 1973 Q1 8.074182 7.948738
## 1973 Q2 8.083421 7.948209
## 1973 Q3 8.088961 7.951771
## 1973 Q4 8.100950 7.948809
## 1974 Q1 8.085364 7.940157
## 1974 Q2 8.077106 7.944137
## 1974 Q3 8.080578 7.948597
## 1974 Q4 8.076298 7.934084
## 1975 Q1 8.072686 7.942575
## 1975 Q2 8.118088 7.958752
## 1975 Q3 8.104432 7.972880
## 1975 Q4 8.114325 7.983235
## 1976 Q1 8.128851 8.002928
## 1976 Q2 8.136665 8.011786
## 1976 Q3 8.145550 8.022208
## 1976 Q4 8.152774 8.035279
## 1977 Q1 8.154701 8.047062
## 1977 Q2 8.165477 8.052456
## 1977 Q3 8.180489 8.061928
## 1977 Q4 8.200316 8.076857
## 1978 Q1 8.206174 8.082495
## 1978 Q2 8.217789 8.103615
## 1978 Q3 8.226199 8.107961
## 1978 Q4 8.234379 8.116059
## 1979 Q1 8.245831 8.121272
## 1979 Q2 8.238854 8.120886
## 1979 Q3 8.244649 8.130736
## 1979 Q4 8.253618 8.133617
## 1980 Q1 8.260855 8.132031
## 1980 Q2 8.242756 8.109466
## 1980 Q3 8.252967 8.120291
## 1980 Q4 8.274051 8.133529
## 1981 Q1 8.270192 8.137776
## 1981 Q2 8.272571 8.138010
## 1981 Q3 8.294824 8.141976
## 1981 Q4 8.297245 8.134380
## 1982 Q1 8.297369 8.140957
## 1982 Q2 8.304470 8.144476
## 1982 Q3 8.308766 8.152140
## 1982 Q4 8.310464 8.170157
## 1983 Q1 8.318840 8.179901
## 1983 Q2 8.326686 8.199602
## 1983 Q3 8.340838 8.217196
## 1983 Q4 8.363226 8.232839
## 1984 Q1 8.386059 8.241413
## 1984 Q2 8.404472 8.255647
## 1984 Q3 8.420638 8.263358
## 1984 Q4 8.430305 8.276446
## 1985 Q1 8.429454 8.293099
## 1985 Q2 8.449642 8.302167
## 1985 Q3 8.442556 8.320959
## 1985 Q4 8.452761 8.323220
## 1986 Q1 8.464678 8.331514
## 1986 Q2 8.475079 8.342054
## 1986 Q3 8.479201 8.359252
## 1986 Q4 8.479159 8.365603
## 1987 Q1 8.487682 8.366045
## 1987 Q2 8.476850 8.379608
## 1987 Q3 8.494395 8.390791
## 1987 Q4 8.507971 8.393239
## 1988 Q1 8.520567 8.409742
## 1988 Q2 8.529379 8.416864
## 1988 Q3 8.537643 8.425188
## 1988 Q4 8.546597 8.436850
## 1989 Q1 8.559544 8.440571
## 1989 Q2 8.556241 8.445267
## 1989 Q3 8.561191 8.454743
## 1989 Q4 8.567867 8.459310
## 1990 Q1 8.576952 8.467394
## 1990 Q2 8.582737 8.470730
## 1990 Q3 8.583655 8.474828
## 1990 Q4 8.576707 8.467646
## 1991 Q1 8.577441 8.463391
## 1991 Q2 8.584310 8.472070
## 1991 Q3 8.586645 8.476392
## 1991 Q4 8.592189 8.475517
## 1992 Q1 8.607747 8.491875
## 1992 Q2 8.615154 8.497603
## 1992 Q3 8.619280 8.507506
## 1992 Q4 8.633945 8.519071
## 1993 Q1 8.614701 8.522638
## 1993 Q2 8.628770 8.531471
## 1993 Q3 8.630057 8.542334
## 1993 Q4 8.644213 8.551189
## 1994 Q1 8.635847 8.562224
## 1994 Q2 8.654012 8.569406
## 1994 Q3 8.660566 8.576537
## 1994 Q4 8.674351 8.586477
## 1995 Q1 8.680400 8.588062
## 1995 Q2 8.679601 8.596319
## 1995 Q3 8.684756 8.604984
## 1995 Q4 8.690020 8.611794
## 1996 Q1 8.701180 8.620580
## 1996 Q2 8.709432 8.631396
## 1996 Q3 8.717895 8.637763
## 1996 Q4 8.722987 8.646044
## 1997 Q1 8.732514 8.656799
## 1997 Q2 8.740641 8.661432
## 1997 Q3 8.751997 8.677729
## 1997 Q4 8.765084 8.688016
## 1998 Q1 8.786212 8.698981
## 1998 Q2 8.800656 8.714748
## 1998 Q3 8.811488 8.726757
## 1998 Q4 8.819222 8.740145
## 1999 Q1 8.826573 8.750097
## 1999 Q2 8.827923 8.765568
## 1999 Q3 8.832879 8.775966
## 1999 Q4 8.847302 8.788258
## 2000 Q1 8.869215 8.804070
## 2000 Q2 8.875916 8.810355
## 2000 Q3 8.888660 8.819961
## 2000 Q4 8.890080 8.828348
## 2001 Q1 8.893298 8.832456
## 2001 Q2 8.889046 8.834963
## 2001 Q3 8.916265 8.839349
## 2001 Q4 8.901911 8.856177
## 2002 Q1 8.929290 8.859676
## 2002 Q2 8.933954 8.865523
## 2002 Q3 8.930031 8.871112
## 2002 Q4 8.930534 8.874616
## 2003 Q1 8.934811 8.879737
## 2003 Q2 8.947117 8.888660
## 2003 Q3 8.962494 8.902850
## 2003 Q4 8.966650 8.908465
## 2004 Q1 8.975719 8.919333
## 2004 Q2 8.981569 8.925388
## 2004 Q3 8.988721 8.933994
## 2004 Q4 9.006852 8.944328
## 2005 Q1 8.998359 8.950196
## 2005 Q2 9.004656 8.958694
## 2005 Q3 9.001519 8.968625
## 2005 Q4 9.017484 8.971690
## 2006 Q1 9.029322 8.982423
## 2006 Q2 9.029849 8.988359
## 2006 Q3 9.034140 8.995140
## 2006 Q4 9.049079 9.004693
## 2007 Q1 9.062293 9.013802
## 2007 Q2 9.060343 9.017277
## 2007 Q3 9.070170 9.024276
## 2007 Q4 9.070526 9.029909
plot(datos1, main="Tendencia")

#Generar Modelo
modelo1=lm(PCE.ts ~DPI.ts)
summary(modelo1)
## 
## Call:
## lm(formula = PCE.ts ~ DPI.ts)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.062816 -0.019375 -0.001032  0.017186  0.078186 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.194233   0.023593  -8.233 1.14e-14 ***
## DPI.ts       1.011351   0.002902 348.542  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02764 on 242 degrees of freedom
## Multiple R-squared:  0.998,  Adjusted R-squared:  0.998 
## F-statistic: 1.215e+05 on 1 and 242 DF,  p-value: < 2.2e-16
residuales=modelo1$residuals
residualPlot(modelo1)

#CON LA PRUEBA CONOCEREMOS SI SON VALORES ESTACIONARIOS O NO ESTACIONARIOS

#Prueba de Dickey-Fuller
y=ur.df(residuales)
summary(y)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.050339 -0.004034 -0.000311  0.005048  0.042822 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -0.05978    0.02364  -2.529   0.0121 *  
## z.diff.lag -0.31809    0.05924  -5.369 1.87e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.009841 on 240 degrees of freedom
## Multiple R-squared:  0.1503, Adjusted R-squared:  0.1433 
## F-statistic: 21.23 on 2 and 240 DF,  p-value: 3.228e-09
## 
## 
## Value of test-statistic is: -2.5289 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62
y@teststat
##                tau1
## statistic -2.528947
y@cval
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62
#con intercepto
y2=ur.df(residuales, type="drift",selectlags="AIC")
summary(y2)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression drift 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + z.diff.lag)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.050207 -0.003904 -0.000181  0.005176  0.042952 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0001294  0.0006339  -0.204   0.8385    
## z.lag.1     -0.0598500  0.0236882  -2.527   0.0122 *  
## z.diff.lag  -0.3180310  0.0593632  -5.357 1.98e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.00986 on 239 degrees of freedom
## Multiple R-squared:  0.1504, Adjusted R-squared:  0.1433 
## F-statistic: 21.15 on 2 and 239 DF,  p-value: 3.478e-09
## 
## 
## Value of test-statistic is: -2.5266 3.2058 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.46 -2.88 -2.57
## phi1  6.52  4.63  3.81
#Sin tendencia y sin intercepto
y3=ur.df(residuales, type="none",selectlags="AIC")
summary(y3)
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.050339 -0.004034 -0.000311  0.005048  0.042822 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -0.05978    0.02364  -2.529   0.0121 *  
## z.diff.lag -0.31809    0.05924  -5.369 1.87e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.009841 on 240 degrees of freedom
## Multiple R-squared:  0.1503, Adjusted R-squared:  0.1433 
## F-statistic: 21.23 on 2 and 240 DF,  p-value: 3.228e-09
## 
## 
## Value of test-statistic is: -2.5289 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62
#Dickey Fuller Aumentada
adf.test(residuales)
## 
##  Augmented Dickey-Fuller Test
## 
## data:  residuales
## Dickey-Fuller = -1.2524, Lag order = 6, p-value = 0.8904
## alternative hypothesis: stationary

#Interpretacion de rechazo

#Al ubicar el valor critico en el grafico de t-student, se encuentra la zona de rechazar y no rechazar la hipotesis H0. Con este dato nosotros interpretamos si H0 es raiz unitaria (No es estacionario); Y por el otro lado H1 nos dice que no hay raiz unitaria(Estacionaria)
#En cada caso, las hipótesis son: Hipótesis nula: H0:δ = 0 (es decir, existe una raíz unitaria, la serie de tiempo es no estacionaria o tiene tendencia estocástica); Hipótesis alternativa: H1:δ < 0 (es decir, la serie de tiempo es estacionaria, posiblemente alrededor de una tendencia determinista).