Raiz y Cointeg
Ana Karina Haique
2025-06-04
library(fpp3)## Registered S3 method overwritten by 'tsibble':
## method from
## as_tibble.grouped_df dplyr## ── Attaching packages ──────────────────────────────────────────── fpp3 1.0.1 ──## ✔ tibble 3.2.1 ✔ tsibble 1.1.6
## ✔ dplyr 1.1.4 ✔ tsibbledata 0.4.1
## ✔ tidyr 1.3.1 ✔ feasts 0.4.1
## ✔ lubridate 1.9.3 ✔ fable 0.4.1
## ✔ ggplot2 3.5.1## ── Conflicts ───────────────────────────────────────────────── fpp3_conflicts ──
## ✖ lubridate::date() masks base::date()
## ✖ dplyr::filter() masks stats::filter()
## ✖ tsibble::intersect() masks base::intersect()
## ✖ tsibble::interval() masks lubridate::interval()
## ✖ dplyr::lag() masks stats::lag()
## ✖ tsibble::setdiff() masks base::setdiff()
## ✖ tsibble::union() masks base::union()datos <- read.csv2("Table 21_1.csv", sep = ";", dec = ",")
datos<- datos |>
mutate(Tiempo = yearquarter(paste(Year, "Q", Quarter)))
datos_ts <- datos |>
select(-Year, -Quarter) |>
as_tsibble(index = Tiempo)Test de la Raiz Unitaria (Dickey-Fuller)
Versión simple: supuesto de residuos \(u_t\) no autocorrelados
L_gdp <- datos_ts|>
mutate(L_GDP = log(GDP)) |>
select(L_GDP)library(urca)
DF<- ur.df(L_gdp$L_GDP, type = "none", lags = 0)
#type = "none", "drift", "trend"
summary(DF)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression none
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.035015 -0.005217 -0.000404 0.005504 0.033038
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## z.lag.1 9.681e-04 7.489e-05 12.93 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009852 on 242 degrees of freedom
## Multiple R-squared: 0.4085, Adjusted R-squared: 0.406
## F-statistic: 167.1 on 1 and 242 DF, p-value: < 2.2e-16
##
##
## Value of test-statistic is: 12.927
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau1 -2.58 -1.95 -1.62\(\delta > 0\) 👎
DF<- ur.df(L_gdp$L_GDP, type = "drift", lags = 0)
#type = "none", "drift", "trend"
summary(DF)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression drift
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.036904 -0.004986 0.000110 0.005221 0.030495
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.022122 0.009088 2.434 0.0157 *
## z.lag.1 -0.001647 0.001077 -1.529 0.1275
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009753 on 241 degrees of freedom
## Multiple R-squared: 0.009612, Adjusted R-squared: 0.005503
## F-statistic: 2.339 on 1 and 241 DF, p-value: 0.1275
##
##
## Value of test-statistic is: -1.5294 88.2173
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau2 -3.46 -2.88 -2.57
## phi1 6.52 4.63 3.81\(\delta = -0.001647\)
Ho: \(\delta = 0\) raÃz unitaria \(=> Y_t\) no estacionaria
H1: \(\delta < 0\) raÃz unitaria \(=> Y_t\) estacionaria
\(\tau_o= -1.5294\) luego, No R Ho
Cómo se interpreta \(\phi\)?
| EstadÃstico | Hipótesis nula (\(H_0\)) |
|---|---|
| tau1, tau2, tau3 | \(\delta = 0\) → la serie tiene raÃz unitaria (no es estacionaria) |
| phi1 | El intercepto no es necesario en el modelo auxiliar |
| phi2 | La tendencia determinista no es necesaria en el modelo auxiliar |
| phi3 | Ni intercepto ni tendencia son necesarios en el modelo auxiliar |
DF<- ur.df(L_gdp$L_GDP, type = "trend", lags = 0)
#type = "none", "drift", "trend"
summary(DF)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression trend
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.037777 -0.004605 0.000245 0.005211 0.031461
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.2093680 0.1102607 1.899 0.0588 .
## z.lag.1 -0.0269282 0.0148755 -1.810 0.0715 .
## tt 0.0002100 0.0001232 1.704 0.0897 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009715 on 240 degrees of freedom
## Multiple R-squared: 0.02145, Adjusted R-squared: 0.0133
## F-statistic: 2.631 on 2 and 240 DF, p-value: 0.07412
##
##
## Value of test-statistic is: -1.8102 60.2439 2.6305
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2 6.22 4.75 4.07
## phi3 8.43 6.49 5.47DF Aumentado
ADF<- ur.df(L_gdp$L_GDP, type = "trend", lags = 4)
#agregamos 4 rezados porque se trata de trimestres
summary(ADF)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression trend
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.033396 -0.004510 -0.000113 0.004492 0.036007
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.2677251 0.1109526 2.413 0.0166 *
## z.lag.1 -0.0351573 0.0149967 -2.344 0.0199 *
## tt 0.0002800 0.0001241 2.256 0.0250 *
## z.diff.lag1 0.2989769 0.0646366 4.626 6.21e-06 ***
## z.diff.lag2 0.1450530 0.0672318 2.158 0.0320 *
## z.diff.lag3 -0.0620944 0.0674589 -0.920 0.3583
## z.diff.lag4 -0.0875919 0.0651833 -1.344 0.1803
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009094 on 232 degrees of freedom
## Multiple R-squared: 0.1617, Adjusted R-squared: 0.1401
## F-statistic: 7.461 on 6 and 232 DF, p-value: 2.552e-07
##
##
## Value of test-statistic is: -2.3443 15.7665 3.4543
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2 6.22 4.75 4.07
## phi3 8.43 6.49 5.47Transformaciones
Caso1: Procesos estacionarios en diferencia
Si una serie de tiempo tiene una raÃz unitaria, las primeras diferencias de tales series son estacionarias. La serie LPIB tiene una raiz unitaria.
L_gdp <- L_gdp |>
mutate(Dif1=difference(L_GDP)) |>
drop_na()
ADF<- ur.df(L_gdp$Dif1, type = "trend", lags = 4)
#type = "none", "drift", "trend"
summary(ADF)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression trend
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.031726 -0.004636 0.000130 0.005084 0.034803
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.232e-03 1.634e-03 5.038 9.5e-07 ***
## z.lag.1 -8.404e-01 1.087e-01 -7.729 3.3e-13 ***
## tt -1.067e-05 8.783e-06 -1.215 0.22549
## z.diff.lag1 1.203e-01 9.738e-02 1.236 0.21784
## z.diff.lag2 2.434e-01 8.729e-02 2.789 0.00573 **
## z.diff.lag3 1.743e-01 7.958e-02 2.190 0.02954 *
## z.diff.lag4 8.921e-02 6.533e-02 1.366 0.17338
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009181 on 231 degrees of freedom
## Multiple R-squared: 0.3664, Adjusted R-squared: 0.3499
## F-statistic: 22.26 on 6 and 231 DF, p-value: < 2.2e-16
##
##
## Value of test-statistic is: -7.7293 19.92 29.8722
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2 6.22 4.75 4.07
## phi3 8.43 6.49 5.47RHo, (en primeras diferencias es estacionaria)
L_gdp |> autoplot(Dif1)
Caso 2: Lo tratamos como un PET
Trataremos al L_GDP como si se tratara de un PET:
\[Y_t = \beta_1+\beta_2t+u_t\]
L_gdp <- datos_ts|>
mutate(L_GDP = log(GDP)) |>
select(L_GDP) |>
mutate(t = row_number())
L_gdp## # A tsibble: 244 x 3 [1Q]
## L_GDP Tiempo t
## <dbl> <qtr> <int>
## 1 7.36 1947 Q1 1
## 2 7.36 1947 Q2 2
## 3 7.36 1947 Q3 3
## 4 7.37 1947 Q4 4
## 5 7.39 1948 Q1 5
## 6 7.41 1948 Q2 6
## 7 7.41 1948 Q3 7
## 8 7.41 1948 Q4 8
## 9 7.40 1949 Q1 9
## 10 7.40 1949 Q2 10
## # ℹ 234 more rowsmodelo <- L_gdp |>
model(TSLM(L_GDP ~ t)) #es el equivalente al lm()
modelo |>
report() #es el equivalente al summary()## Series: L_GDP
## Model: TSLM
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.115056 -0.028572 -0.005401 0.028011 0.091774
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.412e+00 5.416e-03 1368.6 <2e-16 ***
## t 8.255e-03 3.833e-05 215.4 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.04217 on 242 degrees of freedom
## Multiple R-squared: 0.9948, Adjusted R-squared: 0.9948
## F-statistic: 4.639e+04 on 1 and 242 DF, p-value: < 2.22e-16Los residuos son ruido blanco?
uhat <- modelo |>
residuals()
autoplot(uhat)## Plot variable not specified, automatically selected `.vars = .resid`
Claramente los \(u_t\) no son estacionarios.
ADF<- ur.df(uhat$.resid, type = "none", lags = 4) # none son residuos
summary(ADF)##
## ###############################################
## # 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.032403 -0.004386 -0.000156 0.004714 0.036014
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## z.lag.1 -0.03517 0.01497 -2.350 0.0196 *
## z.diff.lag1 0.30415 0.06441 4.722 4.02e-06 ***
## z.diff.lag2 0.14892 0.06707 2.220 0.0273 *
## z.diff.lag3 -0.05874 0.06731 -0.873 0.3837
## z.diff.lag4 -0.08312 0.06499 -1.279 0.2022
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009082 on 234 degrees of freedom
## Multiple R-squared: 0.1567, Adjusted R-squared: 0.1387
## F-statistic: 8.695 on 5 and 234 DF, p-value: 1.408e-07
##
##
## Value of test-statistic is: -2.3496
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau1 -2.58 -1.95 -1.62No R Ho (No es estacionaria)
Cointegración
El caso del consumo vs el ingreso
EGdata <- datos_ts |>
mutate(LDPI = log(DPI), LPCE = log(PCE)) |>
select(LDPI, LPCE)- Primero, vemos si cada serie es de raÃz unitaria
# Ingreso Disponible
ldpi_adf <- ur.df(EGdata$LDPI, type = "trend", selectlags = "AIC")
summary(ldpi_adf)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression trend
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.033055 -0.005616 0.000189 0.005308 0.059535
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.1365188 0.0867179 1.574 0.117
## z.lag.1 -0.0177340 0.0123059 -1.441 0.151
## tt 0.0001341 0.0001067 1.257 0.210
## z.diff.lag -0.0583491 0.0634965 -0.919 0.359
##
## Residual standard error: 0.01008 on 238 degrees of freedom
## Multiple R-squared: 0.02879, Adjusted R-squared: 0.01655
## F-statistic: 2.352 on 3 and 238 DF, p-value: 0.07295
##
##
## Value of test-statistic is: -1.4411 41.2011 3.1349
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2 6.22 4.75 4.07
## phi3 8.43 6.49 5.47# Consumo privado
lpce_adf <- ur.df(EGdata$LPCE, type = "trend", selectlags = "AIC")
summary(lpce_adf)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression trend
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 + 1 + tt + z.diff.lag)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.039712 -0.003734 0.000254 0.004309 0.040272
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.1931977 0.1067023 1.811 0.0715 .
## z.lag.1 -0.0265705 0.0153994 -1.725 0.0857 .
## tt 0.0002277 0.0001351 1.686 0.0931 .
## z.diff.lag 0.0296930 0.0648432 0.458 0.6474
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.008254 on 238 degrees of freedom
## Multiple R-squared: 0.01439, Adjusted R-squared: 0.001971
## F-statistic: 1.159 on 3 and 238 DF, p-value: 0.3263
##
##
## Value of test-statistic is: -1.7254 41.5971 1.7008
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau3 -3.99 -3.43 -3.13
## phi2 6.22 4.75 4.07
## phi3 8.43 6.49 5.47- Segundo, regresion del Consumo sobre el Ingreso Disponible
mod_coint <- EGdata |>
model(TSLM(LPCE ~ LDPI))
mod_coint|>
report()## Series: LPCE
## Model: TSLM
##
## 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 ***
## LDPI 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.22e-16- Tercero, obtenemos los residuos para realizar la prueba de Engel-Granger
uhat <- mod_coint |>
residuals()
eg_test <- ur.df(uhat$.resid, type = "none", lags = 0)
summary(eg_test)##
## ###############################################
## # Augmented Dickey-Fuller Test Unit Root Test #
## ###############################################
##
## Test regression none
##
##
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.051728 -0.005155 -0.000042 0.004666 0.041086
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## z.lag.1 -0.07636 0.02507 -3.046 0.00258 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.01071 on 242 degrees of freedom
## Multiple R-squared: 0.03692, Adjusted R-squared: 0.03294
## F-statistic: 9.277 on 1 and 242 DF, p-value: 0.002578
##
##
## Value of test-statistic is: -3.0458
##
## Critical values for test statistics:
## 1pct 5pct 10pct
## tau1 -2.58 -1.95 -1.62p-valor
library(aTSA)##
## Attaching package: 'aTSA'## The following objects are masked from 'package:fabletools':
##
## estimate, forecast## The following object is masked from 'package:graphics':
##
## identifyattach(EGdata)
coint.test(y = LPCE, X = LDPI, d = 0, nlag = 1) ## Response: LPCE
## Input: LDPI
## Number of inputs: 1
## Model: y ~ X + 1
## -------------------------------
## Engle-Granger Cointegration Test
## alternative: cointegrated
##
## Type 1: no trend
## lag EG p.value
## 1.0000 -3.0458 0.0382
## -----
## Type 2: linear trend
## lag EG p.value
## 1.000 0.597 0.100
## -----
## Type 3: quadratic trend
## lag EG p.value
## 1.00 -3.12 0.10
## -----------
## Note: p.value = 0.01 means p.value <= 0.01
## : p.value = 0.10 means p.value >= 0.10#d=0 implica que esta trabajando con los niveles originales de las variablesRHO al 5% por lo tanto son cointegradas.
- Cuarto, análisis gráfico de los residuos
mod_coint |>
gg_tsresiduals()
Agrego una tendencia
EGdata <- EGdata|>
mutate(indice = row_number())
mod_coint <- EGdata |>
model(TSLM(LPCE ~ LDPI+indice))
mod_coint|>
report()## Series: LPCE
## Model: TSLM
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.048415 -0.009838 0.001095 0.011590 0.043480
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.8130531 0.1317646 21.35 <2e-16 ***
## LDPI 0.5843690 0.0186846 31.27 <2e-16 ***
## indice 0.0037110 0.0001618 22.94 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.01552 on 241 degrees of freedom
## Multiple R-squared: 0.9994, Adjusted R-squared: 0.9994
## F-statistic: 1.928e+05 on 2 and 241 DF, p-value: < 2.22e-16uhat <- mod_coint |>
residuals()
mod_coint |>
gg_tsresiduals()
Test EG:
attach(EGdata)## The following objects are masked from EGdata (pos = 3):
##
## LDPI, LPCE, Tiempocoint.test(y = LPCE,
X = cbind(LDPI, indice),
d = 0,
nlag = 1) ## 1 mas que adf lag## Response: LPCE
## Input: cbind(LDPI, indice)
## Number of inputs: 2
## Model: y ~ X + 1
## -------------------------------
## Engle-Granger Cointegration Test
## alternative: cointegrated
##
## Type 1: no trend
## lag EG p.value
## 1.00 -4.45 0.01
## -----
## Type 2: linear trend
## lag EG p.value
## 1.000 0.126 0.100
## -----
## Type 3: quadratic trend
## lag EG p.value
## 1.00 -1.03 0.10
## -----------
## Note: p.value = 0.01 means p.value <= 0.01
## : p.value = 0.10 means p.value >= 0.10🎯 R Ho al 5% Son cointegradas