Quiz II

#Cargamos los paquetes y los datos que se utilizarán en el trabajo

library(tidyverse); library(dynlm); library(dLagM); library(AER); library(xts)

## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.4.0      ✔ purrr   1.0.1 
## ✔ tibble  3.1.8      ✔ dplyr   1.0.10
## ✔ tidyr   1.3.0      ✔ stringr 1.5.0 
## ✔ readr   2.1.3      ✔ forcats 0.5.2 
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## Loading required package: zoo
## 
## 
## Attaching package: 'zoo'
## 
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## 
## Loading required package: nardl
## 
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo 
## 
## Loading required package: 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
## 
## 
## Loading required package: lmtest
## 
## Loading required package: sandwich
## 
## Loading required package: survival
## 
## 
## Attaching package: 'xts'
## 
## 
## The following objects are masked from 'package:dplyr':
## 
##     first, last

library(ecm); library(openxlsx); library(urca); library(fpp3); library(hunspell);

## ── Attaching packages ────────────────────────────────────────────── fpp3 0.5 ──
## ✔ lubridate   1.9.1     ✔ feasts      0.3.0
## ✔ tsibble     1.1.3     ✔ fable       0.3.2
## ✔ tsibbledata 0.4.1     ✔ fabletools  0.3.2
## ── Conflicts ───────────────────────────────────────────────── fpp3_conflicts ──
## ✖ lubridate::date()      masks base::date()
## ✖ dplyr::filter()        masks stats::filter()
## ✖ xts::first()           masks dplyr::first()
## ✖ fabletools::forecast() masks dLagM::forecast()
## ✖ tsibble::index()       masks zoo::index()
## ✖ tsibble::intersect()   masks base::intersect()
## ✖ tsibble::interval()    masks lubridate::interval()
## ✖ dplyr::lag()           masks stats::lag()
## ✖ xts::last()            masks dplyr::last()
## ✖ fabletools::MAPE()     masks dLagM::MAPE()
## ✖ fabletools::MASE()     masks dLagM::MASE()
## ✖ car::recode()          masks dplyr::recode()
## ✖ tsibble::setdiff()     masks base::setdiff()
## ✖ car::some()            masks purrr::some()
## ✖ tsibble::union()       masks base::union()

library(tidyverse); library(tseries); library(lmtest); library(forecast); library(stats);

## 
## Attaching package: 'forecast'
## 
## The following objects are masked from 'package:fabletools':
## 
##     accuracy, forecast
## 
## The following object is masked from 'package:dLagM':
## 
##     forecast

library(zoo); library(urca); library(googlesheets4); library(ggplot2); library(dplyr);
library(feasts); library(lubridate);library(ggpubr)

## 
## Attaching package: 'ggpubr'
## 
## The following object is masked from 'package:forecast':
## 
##     gghistogram

brent <- read_sheet("https://docs.google.com/spreadsheets/d/1cE0FlBU4LH7GaQCMiLRP_sS5kY7sjkSSbJ8w3KwM78k/edit#gid=0", sheet = 2) |>
  mutate(year = format(as.Date(date, '%Y-%m-%d'), "%d/%m/%Y"),
         UNRATE = as.numeric(as.character(UNRATE)),
         T10Y2YM = as.numeric(as.character(T10Y2YM)),
         FEDFUND = as.numeric(as.character(FEDFUND))) |>
  select(year, UNRATE, T10Y2YM, FEDFUND)

## ! Using an auto-discovered, cached token.
##   To suppress this message, modify your code or options to clearly consent to
##   the use of a cached token.
##   See gargle's "Non-interactive auth" vignette for more details:
##   <]8;;https://gargle.r-lib.org/articles/non-interactive-auth.htmlhttps://gargle.r-lib.org/articles/non-interactive-auth.html]8;;>
## ℹ The googlesheets4 package is using a cached token for ']8;;mailto:majoa1204@gmail.commajoa1204@gmail.com]8;;'.
## ✔ Reading from "AAPL".
## ✔ Range ''BRENT''.

## Warning in as.POSIXlt.POSIXct(x, tz = tz): unknown timezone '%Y-%m-%d'

Descripción de los datos

Para el presente trabajo se utilizarán 3 series de tiempo extraídas de un archivo de GoogleSheets específico. En la serie de tiempo se podrán encontrar las siguientes variables, todas con una incidencia mensual desde el año 1990 hasta el 2023.

• UNRATE (Unemployment Rate): Representa el número de desempleados como porcentaje de la fuerza laboral de Estados Unidos. Los datos de la fuerza laboral están restringidos a personas mayores de 16 años, que actualmente residen en 1 de los 50 estados

• T10Y2YM (10-Year Treasury Constant Maturity Minus 2-Year Treasury Constant Maturity): También conocido como el CMT Index, indica el valor teórico de un Tesoro de Estados Unidos que se basa en los valores recientes de los bonos del Tesos de EEUU subastados

• FEDFUND (Federal Funds Effective Rate): Es la tasa de interés a la que las instituciones de depósito (bancos, cajas de ahorro y cooperativas de credito) intercambian fondos federales (saldos mantenidos en los bancos de la Reserva Federal) entre sí durante la noche. Cuando una institución de depósito tiene saldos excedentes en su cuenta de reserva, presta a otros bancos que necesitan saldos mayores. En términos más simples, un banco con exceso de efectivo, lo que a menudo se conoce como liquidez, prestará a otro banco que necesite aumentar rápidamente la liquidez.

---

Gráfica inicial de la serie para analizar estacionariedad

plot.ts(brent[, c("UNRATE", "T10Y2YM", "FEDFUND")], main = "Figura 1.1: Análisis gráfico de las variables")

En general, se puede observar que las tres variables aparentan ser relativamente variables, siendo la Tasa de Interés de los Bonos la que resalta con mayor variabilidad a lo largo del tiempo. Por su parte, esta última variable parece ser una serie no estacionaria en media y varianza.

Por otro lado, únicamente desde el lado del análisis gráfico, se podría decir que la Tasa de Desempleo es una serie estacionaria en media pero no en varianza y que la Tasa Efectiva de la Reserva podría tratarse de otra serie no estacionaria en media y varianza.

---

Análisis de Correlogramas

ggarrange((ggAcf(brent$UNRATE, main = "AC de UNRATE", col="purple")),
          (ggPacf(brent$UNRATE, main = "PAC de UNRATE", col="purple")),
          (ggAcf(brent$FEDFUND, main = "AC de FEDFUND", col="purple")),
          (ggPacf(brent$FEDFUND, main = "PAC de FEDFUND", col="purple")),
          (ggAcf(brent$T10Y2YM, main = "AC de T10Y2YM", col="purple")),
          (ggPacf(brent$T10Y2YM, main = "PAC de T10Y2YM", col="purple"))
          + 
            rremove("x.text"), ncol = 2, nrow = 3)

## Warning in ggplot2::geom_segment(lineend = "butt", ...): Ignoring unknown parameters: `main`
## Ignoring unknown parameters: `main`
## Ignoring unknown parameters: `main`
## Ignoring unknown parameters: `main`
## Ignoring unknown parameters: `main`
## Ignoring unknown parameters: `main`

Podemos observar los siguientes resultados:

• La serie de la Tasa de Desempleo corresponde a un proceso autorregresivo de primer orden (AR(1))

• La serie de la Tasa de Interés de los Bonos corresponde a un proceso de tercer orden (AR(3)) con picos significativos en los rezagos 9 y 16

• La serie de la Tasa Efectiva de la Reserva corresponde a un proceso autorregresivo de tercer orden (AR(3))

---

Test de Dickey-Fuller

###Cargar funcion que interpreta resultados de dickey fuller
interp_urdf <- function(urdf, level="5pct") {
  if(class(urdf) != "ur.df") stop('parameter is not of class ur.df from urca package')
  if(!(level %in% c("1pct", "5pct", "10pct") ) ) stop('parameter level is not one of 1pct, 5pct, or 10pct')
  
  cat("========================================================================\n")
  cat( paste("At the", level, "level:\n") )
  if(urdf@model == "none") {
    cat("The model is of type none\n")
    tau1_crit = urdf@cval["tau1",level]
    tau1_teststat = urdf@teststat["statistic","tau1"]
    tau1_teststat_wi_crit = tau1_teststat > tau1_crit
    if(tau1_teststat_wi_crit) {
      cat("tau1: The null hypothesis is not rejected, unit root is present\n")
    } else {
      cat("tau1: The null hypothesis is rejected, unit root is not present\n")
    }
  } else if(urdf@model == "drift") {
    cat("The model is of type drift\n")
    tau2_crit = urdf@cval["tau2",level]
    tau2_teststat = urdf@teststat["statistic","tau2"]
    tau2_teststat_wi_crit = tau2_teststat > tau2_crit
    phi1_crit = urdf@cval["phi1",level]
    phi1_teststat = urdf@teststat["statistic","phi1"]
    phi1_teststat_wi_crit = phi1_teststat < phi1_crit
    if(tau2_teststat_wi_crit) {
      # Unit root present branch
      cat("tau2: The first null hypothesis is not rejected, unit root is present\n")
      if(phi1_teststat_wi_crit) {
        cat("phi1: The second null hypothesis is not rejected, unit root is present\n")
        cat("      and there is no drift.\n")
      } else {
        cat("phi1: The second null hypothesis is rejected, unit root is present\n")
        cat("      and there is drift.\n")
      }
    } else {
      # Unit root not present branch
      cat("tau2: The first null hypothesis is rejected, unit root is not present\n")
      if(phi1_teststat_wi_crit) {
        cat("phi1: The second null hypothesis is not rejected, unit root is present\n")
        cat("      and there is no drift.\n")
        warning("This is inconsistent with the first null hypothesis.")
      } else {
        cat("phi1: The second null hypothesis is rejected, unit root is not present\n")
        cat("      and there is drift.\n")
      }
    }
  } else if(urdf@model == "trend") {
    cat("The model is of type trend\n")
    tau3_crit = urdf@cval["tau3",level]
    tau3_teststat = urdf@teststat["statistic","tau3"]
    tau3_teststat_wi_crit = tau3_teststat > tau3_crit
    phi2_crit = urdf@cval["phi2",level]
    phi2_teststat = urdf@teststat["statistic","phi2"]
    phi2_teststat_wi_crit = phi2_teststat < phi2_crit
    phi3_crit = urdf@cval["phi3",level]
    phi3_teststat = urdf@teststat["statistic","phi3"]
    phi3_teststat_wi_crit = phi3_teststat < phi3_crit
    if(tau3_teststat_wi_crit) {
      # First null hypothesis is not rejected, Unit root present branch
      cat("tau3: The first null hypothesis is not rejected, unit root is present\n")
      if(phi3_teststat_wi_crit) {
        # Second null hypothesis is not rejected
        cat("phi3: The second null hypothesis is not rejected, unit root is present\n")
        cat("      and there is no trend\n")
        if(phi2_teststat_wi_crit) {
          # Third null hypothesis is not rejected
          # a0-drift = gamma = a2-trend = 0
          cat("phi2: The third null hypothesis is not rejected, unit root is present\n")
          cat("      there is no trend, and there is no drift\n")
        } else {
          # Third null hypothesis is rejected
          cat("phi2: The third null hypothesis is rejected, unit root is present\n")
          cat("      there is no trend, and there is drift\n")
        }
      }
      else {
        # Second null hypothesis is rejected
        cat("phi3: The second null hypothesis is rejected, unit root is present\n")
        cat("      and there is trend\n")
        if(phi2_teststat_wi_crit) {
          # Third null hypothesis is not rejected
          # a0-drift = gamma = a2-trend = 0
          cat("phi2: The third null hypothesis is not rejected, unit root is present\n")
          cat("      there is no trend, and there is no drift\n")
          warning("This is inconsistent with the second null hypothesis.")
        } else {
          # Third null hypothesis is rejected
          cat("phi2: The third null hypothesis is rejected, unit root is present\n")
          cat("      there is trend, and there may or may not be drift\n")
          warning("Presence of drift is inconclusive.")
        }
      }
    } else {
      # First null hypothesis is rejected, Unit root not present branch
      cat("tau3: The first null hypothesis is rejected, unit root is not present\n")
      if(phi3_teststat_wi_crit) {
        cat("phi3: The second null hypothesis is not rejected, unit root is present\n")
        cat("      and there is no trend\n")
        warning("This is inconsistent with the first null hypothesis.")
        if(phi2_teststat_wi_crit) {
          # Third null hypothesis is not rejected
          # a0-drift = gamma = a2-trend = 0
          cat("phi2: The third null hypothesis is not rejected, unit root is present\n")
          cat("      there is no trend, and there is no drift\n")
          warning("This is inconsistent with the first null hypothesis.")
        } else {
          # Third null hypothesis is rejected
          cat("phi2: The third null hypothesis is rejected, unit root is not present\n")
          cat("      there is no trend, and there is drift\n")
        }
      } else {
        cat("phi3: The second null hypothesis is rejected, unit root is not present\n")
        cat("      and there may or may not be trend\n")
        warning("Presence of trend is inconclusive.")
        if(phi2_teststat_wi_crit) {
          # Third null hypothesis is not rejected
          # a0-drift = gamma = a2-trend = 0
          cat("phi2: The third null hypothesis is not rejected, unit root is present\n")
          cat("      there is no trend, and there is no drift\n")
          warning("This is inconsistent with the first and second null hypothesis.")
        } else {
          # Third null hypothesis is rejected
          cat("phi2: The third null hypothesis is rejected, unit root is not present\n")
          cat("      there may or may not be trend, and there may or may not be drift\n")
          warning("Presence of trend and drift is inconclusive.")
        }
      }
    }
  } else warning('urdf model type is not one of none, drift, or trend')
  cat("========================================================================\n")
}

• DF de la variable UNRATE (Tasa de Desempleo)

#DF con Tendencia
DF.bUNRATE1 <- ur.df(brent$UNRATE, type = "trend", selectlags = "AIC")
summary(DF.bUNRATE1)

## 
## ############################################### 
## # 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 
## -1.6837 -0.1375 -0.0267  0.0901 10.2017 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.3472409  0.1147545   3.026 0.002642 ** 
## z.lag.1     -0.0553526  0.0166868  -3.317 0.000995 ***
## tt          -0.0001461  0.0002500  -0.585 0.559206    
## z.diff.lag   0.0530132  0.0506017   1.048 0.295444    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.566 on 391 degrees of freedom
## Multiple R-squared:  0.02845,    Adjusted R-squared:  0.02099 
## F-statistic: 3.816 on 3 and 391 DF,  p-value: 0.01021
## 
## 
## Value of test-statistic is: -3.3171 3.7427 5.6005 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.98 -3.42 -3.13
## phi2  6.15  4.71  4.05
## phi3  8.34  6.30  5.36

#Significativo al 10% pero no al 5 ni al 1


#DF con deriva
DF.bUNRATE2 <- ur.df(brent$UNRATE, type = "drift", selectlags = "AIC")
summary(DF.bUNRATE2)

## 
## ############################################### 
## # 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 
## -1.7104 -0.1410 -0.0137  0.0918 10.1779 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  0.31575    0.10124   3.119  0.00195 **
## z.lag.1     -0.05494    0.01666  -3.298  0.00106 **
## z.diff.lag   0.05342    0.05055   1.057  0.29133   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5655 on 392 degrees of freedom
## Multiple R-squared:  0.0276, Adjusted R-squared:  0.02264 
## F-statistic: 5.563 on 2 and 392 DF,  p-value: 0.004146
## 
## 
## Value of test-statistic is: -3.2981 5.4524 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.44 -2.87 -2.57
## phi1  6.47  4.61  3.79

#Estadístico mayor y significativo al 5% y al 10%


#DF sin tendencia y sin seriva
DF.bUNRATE3 <- ur.df(brent$UNRATE, type = "none", selectlags = "AIC")
summary(DF.bUNRATE3)

## 
## ############################################### 
## # 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 
## -2.0923 -0.0766  0.0239  0.1246 10.2980 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## z.lag.1    -0.005081   0.004734  -1.073    0.284
## z.diff.lag  0.027104   0.050396   0.538    0.591
## 
## Residual standard error: 0.5718 on 393 degrees of freedom
## Multiple R-squared:  0.003541,   Adjusted R-squared:  -0.00153 
## F-statistic: 0.6983 on 2 and 393 DF,  p-value: 0.4981
## 
## 
## Value of test-statistic is: -1.0733 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62

#Estadístico mayor para todos los niveles de significancia

Es necesario diferenciar la serie

#Diferenciar tasa de desempleo

brent <- brent |>
  mutate(dUNRATE = c(NA, diff(brent$UNRATE)))

DF.bUNRATE4 <- ur.df(brent[!is.na(brent$dUNRATE),]$dUNRATE, type = "trend",
                     selectlags = "AIC")
summary(DF.bUNRATE4)

## 
## ############################################### 
## # 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 
## -1.9889 -0.1090 -0.0107  0.0940 10.3014 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.0201445  0.0579014   0.348   0.7281    
## z.lag.1     -1.0806523  0.0703220 -15.367   <2e-16 ***
## tt          -0.0001264  0.0002532  -0.499   0.6180    
## z.diff.lag   0.1075201  0.0503402   2.136   0.0333 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5713 on 390 degrees of freedom
## Multiple R-squared:  0.4938, Adjusted R-squared:  0.4899 
## F-statistic: 126.8 on 3 and 390 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -15.3672 78.7176 118.0764 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.98 -3.42 -3.13
## phi2  6.15  4.71  4.05
## phi3  8.34  6.30  5.36

plot(DF.bUNRATE4)

####Interpretacion del test de Dicket Fuller para la serie diferenciada

intbUNRATE <- interp_urdf(urdf = ur.df(brent[!is.na(brent$dUNRATE),]$dUNRATE, type = "trend", selectlags = "AIC"))

## ========================================================================
## At the 5pct level:
## The model is of type trend
## tau3: The first null hypothesis is rejected, unit root is not present
## phi3: The second null hypothesis is rejected, unit root is not present
##       and there may or may not be trend

## Warning in interp_urdf(urdf = ur.df(brent[!is.na(brent$dUNRATE), ]$dUNRATE, :
## Presence of trend is inconclusive.

## phi2: The third null hypothesis is rejected, unit root is not present
##       there may or may not be trend, and there may or may not be drift

## Warning in interp_urdf(urdf = ur.df(brent[!is.na(brent$dUNRATE), ]$dUNRATE, :
## Presence of trend and drift is inconclusive.

## ========================================================================

Con la nueva serie diferenciada se rechazan las hipótesis nulas de tau3, phi3 y phi2, y es posible afirmar que no hay raíz unitaria presente; no obstante, para phi3 la presencia de tendencia es inconclusa y para phi2 la presencia tanto de tendencia como de deriva es inconclusa también.

• DF de la variable T10Y2YM (Tasa de Interés de los Bonos)

#DF con tendencia
DF.bT10Y2YM1 <- ur.df(brent$T10Y2YM, type = "trend", selectlags = "AIC")
summary(DF.bT10Y2YM1)

## 
## ############################################### 
## # 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.73643 -0.07320 -0.00426  0.06185  0.55638 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.021e-02  1.621e-02   1.863   0.0632 .  
## z.lag.1     -1.507e-02  7.822e-03  -1.926   0.0548 .  
## tt          -7.429e-05  5.990e-05  -1.240   0.2156    
## z.diff.lag   3.171e-01  4.815e-02   6.584 1.47e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1353 on 391 degrees of freedom
## Multiple R-squared:  0.1097, Adjusted R-squared:  0.1029 
## F-statistic: 16.06 on 3 and 391 DF,  p-value: 7.28e-10
## 
## 
## Value of test-statistic is: -1.9263 1.772 2.6384 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.98 -3.42 -3.13
## phi2  6.15  4.71  4.05
## phi3  8.34  6.30  5.36

#No se puede rechazar la hipótesis nula de raíz unitaria debido a que el estadístico es menor para todos los niveles de significancia

#DF con deriva
DF.bT10Y2YM2 <- ur.df(brent$T10Y2YM, type = "drift", selectlags = "AIC")
summary(DF.bT10Y2YM2)

## 
## ############################################### 
## # 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.73907 -0.07130 -0.00744  0.06082  0.55375 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.015495   0.011063   1.401   0.1621    
## z.lag.1     -0.015124   0.007827  -1.932   0.0541 .  
## z.diff.lag   0.322166   0.048009   6.711 6.81e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1353 on 392 degrees of freedom
## Multiple R-squared:  0.1062, Adjusted R-squared:  0.1016 
## F-statistic: 23.29 on 2 and 392 DF,  p-value: 2.776e-10
## 
## 
## Value of test-statistic is: -1.9322 1.8864 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.44 -2.87 -2.57
## phi1  6.47  4.61  3.79

#Tampoco se puede rechazar la hipótesis nula puesto que el estadístico del test es menor para todos los niveles de significancia

#DF sin tendencia y sin deriva
DF.bT10Y2YM3 <- ur.df(brent$T10Y2YM, type = "none", selectlags = "AIC")
summary(DF.bT10Y2YM3)

## 
## ############################################### 
## # 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.74289 -0.06127  0.00344  0.06807  0.55619 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -0.006484   0.004825  -1.344     0.18    
## z.diff.lag  0.316089   0.047871   6.603 1.31e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1355 on 393 degrees of freedom
## Multiple R-squared:  0.1019, Adjusted R-squared:  0.09733 
## F-statistic: 22.29 on 2 and 393 DF,  p-value: 6.741e-10
## 
## 
## Value of test-statistic is: -1.3441 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62

#Estadístico es mayor para todos los niveles de significancia

Es necesario diferenciar la serie de T10Y2YM

###Diferenciar T10Y2YM

brent <- brent |>
  mutate(dT10Y2YM = c(NA, diff(brent$T10Y2YM)))

DF.bT10Y2YM4 <- ur.df(brent[!is.na(brent$dT10Y2YM),]$dT10Y2YM, type = "trend",
                      selectlags = "AIC")
summary(DF.bT10Y2YM4)

## 
## ############################################### 
## # 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.69948 -0.07533 -0.00367  0.06140  0.54918 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.627e-02  1.375e-02   1.183   0.2373    
## z.lag.1     -7.606e-01  5.922e-02 -12.844   <2e-16 ***
## tt          -8.775e-05  6.023e-05  -1.457   0.1459    
## z.diff.lag   9.720e-02  5.032e-02   1.932   0.0541 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1352 on 390 degrees of freedom
## Multiple R-squared:  0.3536, Adjusted R-squared:  0.3487 
## F-statistic: 71.13 on 3 and 390 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -12.8439 54.9931 82.4887 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.98 -3.42 -3.13
## phi2  6.15  4.71  4.05
## phi3  8.34  6.30  5.36

plot(DF.bT10Y2YM4)

###Interpretacion del test de Dicket Fuller para la serie diferenciada de T10Y2YM

intbT10Y2YM <- interp_urdf(urdf = ur.df(brent[!is.na(brent$dT10Y2YM),]$dT10Y2YM, type = "trend", selectlags = "AIC"))

## ========================================================================
## At the 5pct level:
## The model is of type trend
## tau3: The first null hypothesis is rejected, unit root is not present
## phi3: The second null hypothesis is rejected, unit root is not present
##       and there may or may not be trend

## Warning in interp_urdf(urdf = ur.df(brent[!is.na(brent$dT10Y2YM), ]$dT10Y2YM, :
## Presence of trend is inconclusive.

## phi2: The third null hypothesis is rejected, unit root is not present
##       there may or may not be trend, and there may or may not be drift

## Warning in interp_urdf(urdf = ur.df(brent[!is.na(brent$dT10Y2YM), ]$dT10Y2YM, :
## Presence of trend and drift is inconclusive.

## ========================================================================

Se rechaza exitosamente la hipótesis nula y se concluye estacionariedad para la nueva serie diferenciada de T10Y2YM

• DF de la variable FEDFUND (Tasa Efectiva de la Reserva)

#DF con tendencia
DF.FEDFUND1 <- ur.df(brent$FEDFUND, type = "trend", selectlags = "AIC")
summary(DF.FEDFUND1)

## 
## ############################################### 
## # 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.95216 -0.03708 -0.00621  0.05071  0.52884 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.079e-02  3.108e-02   0.669   0.5040    
## z.lag.1     -7.929e-03  4.689e-03  -1.691   0.0917 .  
## tt          -1.385e-05  9.790e-05  -0.141   0.8875    
## z.diff.lag   6.305e-01  3.980e-02  15.841   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1441 on 391 degrees of freedom
## Multiple R-squared:  0.4084, Adjusted R-squared:  0.4039 
## F-statistic: 89.97 on 3 and 391 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -1.6908 2.005 2.9058 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.98 -3.42 -3.13
## phi2  6.15  4.71  4.05
## phi3  8.34  6.30  5.36

#El estadístico es menor para todos los niveles de significancia por lo que no se puede rechazar la hipótesis nula de raíz unitaria y no estacionariedad


#DF con deriva
DF.FEDFUND2 <- ur.df(brent$FEDFUND, type = "drift", selectlags = "AIC")
summary(DF.FEDFUND2)

## 
## ############################################### 
## # 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.95382 -0.03681 -0.00608  0.05127  0.52981 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.016676   0.011008   1.515   0.1306    
## z.lag.1     -0.007429   0.003083  -2.410   0.0164 *  
## z.diff.lag   0.629397   0.038967  16.152   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1439 on 392 degrees of freedom
## Multiple R-squared:  0.4084, Adjusted R-squared:  0.4054 
## F-statistic: 135.3 on 2 and 392 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -2.4096 3.005 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau2 -3.44 -2.87 -2.57
## phi1  6.47  4.61  3.79

#Al igual que para el DF con tendencia, no es posible rechazar la hipótesis nula de raíz unitaria


#DF sin tendencia y sin deriva
DF.FEDFUND3 <- ur.df(brent$FEDFUND, type = "none", selectlags = "AIC")
summary(DF.FEDFUND3)

## 
## ############################################### 
## # 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.94269 -0.02684  0.00660  0.05548  0.52976 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -0.003916   0.002035  -1.924   0.0551 .  
## z.diff.lag  0.629142   0.039031  16.119   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1442 on 393 degrees of freedom
## Multiple R-squared:  0.4066, Adjusted R-squared:  0.4036 
## F-statistic: 134.6 on 2 and 393 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -1.9242 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62

#Se rechaza la hipótesis nula al 10% pero no es posible hacer lo mismo para el 5% y el 1%

Es necesario diferenciar la serie de T10Y2YM

###Diferenciar FEDFUND

brent <- brent |>
  mutate(dFEDFUND = c(NA, diff(brent$FEDFUND)))

DF.FEDFUND4 <- ur.df(brent[!is.na(brent$dFEDFUND),]$dFEDFUND, type = "trend",
                     selectlags = "AIC")
summary(DF.FEDFUND4)

## 
## ############################################### 
## # 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.95894 -0.03236 -0.00046  0.04793  0.48853 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.227e-02  1.473e-02  -1.511   0.1315    
## z.lag.1     -3.156e-01  4.358e-02  -7.242 2.37e-12 ***
## tt           9.855e-05  6.430e-05   1.533   0.1262    
## z.diff.lag  -1.635e-01  5.009e-02  -3.263   0.0012 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1428 on 390 degrees of freedom
## Multiple R-squared:  0.2102, Adjusted R-squared:  0.2041 
## F-statistic: 34.59 on 3 and 390 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -7.2425 17.5082 26.2571 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau3 -3.98 -3.42 -3.13
## phi2  6.15  4.71  4.05
## phi3  8.34  6.30  5.36

plot(DF.FEDFUND4)

###Interpretacion del test de Dicket Fuller para la serie diferenciada de FEDFUND

intFEDFUND <- interp_urdf(urdf = ur.df(brent[!is.na(brent$dFEDFUND),]$dFEDFUND, type = "trend", selectlags = "AIC"))

## ========================================================================
## At the 5pct level:
## The model is of type trend
## tau3: The first null hypothesis is rejected, unit root is not present
## phi3: The second null hypothesis is rejected, unit root is not present
##       and there may or may not be trend

## Warning in interp_urdf(urdf = ur.df(brent[!is.na(brent$dFEDFUND), ]$dFEDFUND, :
## Presence of trend is inconclusive.

## phi2: The third null hypothesis is rejected, unit root is not present
##       there may or may not be trend, and there may or may not be drift

## Warning in interp_urdf(urdf = ur.df(brent[!is.na(brent$dFEDFUND), ]$dFEDFUND, :
## Presence of trend and drift is inconclusive.

## ========================================================================

Con el test de Dickey-Fuller para la serie diferenciada es posible rechazar la hipótesis nula para tau3, phi3 y phi2.

---

Ahora que las series están diferenciadas y son estacionarias, es posible trabajar con ellas y análizar sus rezagos a través del método de Koyck y proyectar con ADL y ARIMAX

---

Modelo de Rezago Distribuido

Un Modelo de Rezagos Distribuidos (ADL) son regresiones estándar de mínimos cuadrados que incluyen retrasos tanto de la variable dependiente como de las variables explicativas como regresores. Es un método para examinar las relaciones de cointegración entre variables.

Un modelo econométrico de este tipo es utilizado en series temporales en las que una o varias variables explicativas pueden tener efectos sobre la variable dependiente tras uno o más periodos.

lag1 <- dlm(formula = dFEDFUND ~ dT10Y2YM + dUNRATE, data = data.frame(brent),
                 q = 2)
summary(lag1)

## 
## Call:
## lm(formula = as.formula(model.formula), data = design)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.75061 -0.06571  0.00324  0.08102  0.59169 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.012143   0.007716  -1.574  0.11636    
## dT10Y2YM.t  -0.504561   0.057390  -8.792  < 2e-16 ***
## dT10Y2YM.1  -0.183752   0.060491  -3.038  0.00255 ** 
## dT10Y2YM.2  -0.318672   0.057418  -5.550 5.31e-08 ***
## dUNRATE.t   -0.064238   0.013647  -4.707 3.50e-06 ***
## dUNRATE.1   -0.005389   0.013588  -0.397  0.69188    
## dUNRATE.2   -0.006567   0.013627  -0.482  0.63013    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1531 on 387 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.3388, Adjusted R-squared:  0.3285 
## F-statistic: 33.05 on 6 and 387 DF,  p-value: < 2.2e-16
## 
## AIC and BIC values for the model:
##         AIC       BIC
## 1 -351.6746 -319.8638

Lo que buscamos al realizar una Regresión Lineal Dinámica es observar si las variables exógenas (Tasa de Interés de los Bonos y la Tasa de Desempleo) y los dos últimos rezagos de las variables exógenas influyen de alguna manera a la variable endógena (Tasa Efectiva de la Reserva). El summary del lag1 indica que la variable T10Y2YM es efecto importante e influye de manera significativa a la Tasa de Efectiva de la Reserva, dado que por su parte, se puede observar que el valor actual de la variable “Tasa de Interés de los Bonos” es significativa al 0.001, igual que el rezago 2 de la misma variable, mientras que el rezago 1 vendría a ser significante al 0.05.

En cuanto a la Tasa de Desempleo, se puede notar que esta variable no influye tanto como la anterior en nuestra variable endógena y que únicamente el valor actual de la serie influye significativamente a la Tasa Efectiva de la Reserva pues los valores de los rezagos 1 y 2 no son significativos para nada.

Por este motivo, se utilizarán las variables FEDFUND y T10Y2YM para el método de Koyck

---

Modelo de rezago distribuido: método de Koyck

También conocido como el método de distribución infinita de rezagos, consiste en introducir en el modelo original una restricción sobre los coeficientes de rezago.

#Koyck de x=T10Y2YM con y=FEDFUND

koyck1 <- koyckDlm(brent$dT10Y2YM, brent$dFEDFUND)
summary(koyck1, diagnostic = T)

## 
## Call:
## "Y ~ (Intercept) + Y.1 + X.t"
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.818363 -0.064612  0.004479  0.062796  0.420666 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.005184   0.006809  -0.761   0.4469    
## Y.1          0.543633   0.052221  10.410   <2e-16 ***
## X.t         -0.491048   0.205231  -2.393   0.0172 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1342 on 392 degrees of freedom
## Multiple R-Squared: 0.4858,  Adjusted R-squared: 0.4832 
## Wald test: 155.2 on 2 and 392 DF,  p-value: < 2.2e-16 
## 
## Diagnostic tests:
##                  df1 df2  statistic      p-value
## Weak instruments   1 392 23.4357538 1.861385e-06
## Wu-Hausman         1 391  0.1898209 6.633062e-01
## 
##                                alpha       beta       phi
## Geometric coefficients:  -0.01136002 -0.4910481 0.5436332

Se puede observar que aunque el intercepto no es significativo, Y.1 es significativo al 0.001 y que X.t es significativo al 0.05. ¿Qué se puede decir acerca de los estadísticos?

• Weak instruments = 23.4: Si el estadístico F de la “primera etapa” es menor a 10, esto indica la presencia de un instrumento débil. No es nuestro caso gracias a que el estadístico es de 23.4(>10) — En caso de que este estadístico sí indicara la presencia de un instrumento débil, se entendería que la variable instrumental no está suficientemente correlacionada con la variable endógena, en dicho caso, los estimadores del modelo lineal corregido son inconsistentes

• Wu-Hausman = 0.18: Dado que el p-valor es menor a 0.05, se tendría que rechazar la hipótesis nula de que el modelo de efectos fijos es apropiado y se tomaría la hipótesis alternativa de que el modelo de Mínimos Cuadrados Ordinarios (MCO) para datos agrupados es mejor. En otras palabras, bajo la hipótesis nula, los coeficientes de los modelos de efectos fijos y aleatorios son consistentes. Sin embargo, solo los coeficientes del modelo de efectos aleatorios son eficientes. Si no podemos rechazar la hipótesis nula utilizando la prueba de Wu-Hausman, significa que se debe preferir el modelo de efectos aleatorios.

koyck11 <- lm(dFEDFUND ~ lag(dFEDFUND) + dT10Y2YM, data = brent)
summary(koyck11)

## 
## Call:
## lm(formula = dFEDFUND ~ lag(dFEDFUND) + dT10Y2YM, data = brent)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.84082 -0.06013  0.00451  0.05982  0.43717 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   -0.004848   0.006737  -0.720    0.472    
## lag(dFEDFUND)  0.559476   0.037223  15.030  < 2e-16 ***
## dT10Y2YM      -0.404452   0.048548  -8.331 1.36e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1336 on 392 degrees of freedom
##   (2 observations deleted due to missingness)
## Multiple R-squared:  0.4899, Adjusted R-squared:  0.4873 
## F-statistic: 188.3 on 2 and 392 DF,  p-value: < 2.2e-16

Al correr el modeo de rezagos distribuidos de Koyck tomando como variable endógena la series diferencias Tasa Efectiva de la Reserva y como exógenas el rezago de la Tasa Efectiva y la serie diferenciada de la Tasa de Interés de los Bonos, podemos observar que ambas son significativas al 0.001 nivel de significancia.

Observamos también que el R-Cuadrado Multiple es del 48.99% lo que indica que el modelo de regresión explica aproximadamente el 50% o la mitad de los datos observados.

---

H de Durbin

La prueba “h” de Durbin rezada los residuales de MCO (como vimos en el test Wu-Hausman) en sus propios rezagos y la lista original de regresores. El coeficiente de la serie residual rezagada, en relación con su error estándar estimado, se distribuye ‘t’ bajo la autocorrelación nula de cero en el proceso de error.

h <- durbinH(koyck11, "lag(dFEDFUND)")
abs(h)

## [1] 3.961956

2*pnorm(-abs(h)) # Calculando el pvalor

## [1] 7.433819e-05

El test arroja una autocorrelación mayor a 2 lo cual indica autocorrelación negativa. Este tipo de autocorrelación implica que si un valor particular está por encima del promedio, es más probable que el siguiente valor (o el valor anterior) esté por debajo del promedio. Si un valor particular está por debajo del promedio, es probable que el siguiente valor esté por encima del promedio. Por otro lado, al ser el p-valor menor a 0.05, no se puede rechazar la hipótesis nula de no autocorrelación y se toma la hipótesis alternativa de que existe algún tipo de autocorrelación en los residuales que, en este caso, resulta ser negativa.