library(readxl)
library(tidyverse) #manipulación y visualización de datos
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(forecast) #análisis y pronóstico de series de tiempo
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
library(quantmod) #recuperación y manipulación de datos financieros
## Cargando paquete requerido: xts
## Cargando paquete requerido: zoo
## 
## Adjuntando el paquete: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## 
## ######################### Warning from 'xts' package ##########################
## #                                                                             #
## # The dplyr lag() function breaks how base R's lag() function is supposed to  #
## # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or       #
## # source() into this session won't work correctly.                            #
## #                                                                             #
## # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
## # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop           #
## # dplyr from breaking base R's lag() function.                                #
## #                                                                             #
## # Code in packages is not affected. It's protected by R's namespace mechanism #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning.  #
## #                                                                             #
## ###############################################################################
## 
## Adjuntando el paquete: 'xts'
## 
## The following objects are masked from 'package:dplyr':
## 
##     first, last
## 
## Cargando paquete requerido: TTR
library(AER) #métodos econométricos aplicados
## Cargando paquete requerido: car
## Cargando paquete requerido: carData
## 
## Adjuntando el paquete: 'car'
## 
## The following object is masked from 'package:dplyr':
## 
##     recode
## 
## The following object is masked from 'package:purrr':
## 
##     some
## 
## Cargando paquete requerido: lmtest
## Cargando paquete requerido: sandwich
## Cargando paquete requerido: survival
library(MASS) #análisis estadístico aplicado
## 
## Adjuntando el paquete: 'MASS'
## 
## The following object is masked from 'package:dplyr':
## 
##     select
library(tseries) #análisis de series de tiempo
library(stats) #funciones estadísticas básicas
library(car) #diagnósticos y pruebas de regresión lineal
library(lmtest) #pruebas adicionales para modelos de regresión lineal
library(urca) #pruebas de raíz unitaria y cointegración
library(aTSA) #análisis alternativo de series de tiempo
## 
## Adjuntando el paquete: 'aTSA'
## 
## The following objects are masked from 'package:tseries':
## 
##     adf.test, kpss.test, pp.test
## 
## The following object is masked from 'package:forecast':
## 
##     forecast
## 
## The following object is masked from 'package:graphics':
## 
##     identify
library(urca)
library(vars)
## Cargando paquete requerido: strucchange
## 
## Adjuntando el paquete: 'strucchange'
## 
## The following object is masked from 'package:stringr':
## 
##     boundary
## 
## 
## Adjuntando el paquete: 'vars'
## 
## The following object is masked from 'package:aTSA':
## 
##     arch.test
library(rugarch)
## Cargando paquete requerido: parallel
## 
## Adjuntando el paquete: 'rugarch'
## 
## The following object is masked from 'package:purrr':
## 
##     reduce
## 
## The following object is masked from 'package:stats':
## 
##     sigma
Crudo <- read_excel("Crudo_Datos.xlsx")
Acero <- read_excel("Acero_Datos.xlsx")

summary(Crudo)
##      Fecha                         Precio      
##  Min.   :1996-01-03 00:00:00   Min.   : -2.37  
##  1st Qu.:2003-02-06 18:00:00   1st Qu.: 25.32  
##  Median :2010-03-15 12:00:00   Median : 49.53  
##  Mean   :2010-03-14 21:35:46   Mean   : 52.54  
##  3rd Qu.:2017-04-19 06:00:00   3rd Qu.: 72.81  
##  Max.   :2024-05-24 00:00:00   Max.   :132.71  
##                                NA's   :554
summary(Acero)
##      Date                Open            High            Low       
##  Length:1361        Min.   : 4.80   Min.   : 5.34   Min.   : 4.54  
##  Class :character   1st Qu.:13.81   1st Qu.:14.07   1st Qu.:13.54  
##  Mode  :character   Median :22.28   Median :22.84   Median :21.74  
##                     Mean   :21.95   Mean   :22.41   Mean   :21.51  
##                     3rd Qu.:26.39   3rd Qu.:26.91   3rd Qu.:25.79  
##                     Max.   :49.77   Max.   :50.20   Max.   :49.24  
##      Close         Adj Close          Volume         
##  Min.   : 4.90   Min.   : 4.775   Min.   :   787700  
##  1st Qu.:13.79   1st Qu.:13.345   1st Qu.:  7287100  
##  Median :22.33   Median :21.821   Median : 11265200  
##  Mean   :21.95   Mean   :21.620   Mean   : 13127632  
##  3rd Qu.:26.36   3rd Qu.:25.923   3rd Qu.: 16903400  
##  Max.   :49.59   Max.   :49.472   Max.   :112841400
Crudo <- data.frame(Fecha = as.Date(Crudo$Fecha), Precio = Crudo$Precio)
Acero <- data.frame(Fecha = as.Date(Acero$Date), Precio = Acero$Close)
observaciones_adicionales <- setdiff(rownames(Crudo), rownames(Acero))
Crudo <- Crudo[!(rownames(Crudo) %in% observaciones_adicionales), ]
PETROLEOSP <- ts(Crudo$Precio, start = c(2019, 01), frequency = 265)
ACEROSP <- ts(Acero$Precio, start = c(2019, 01), frequency = 265)

plot(ACEROSP, main = "Precio histórico Acero", xlab = "Fecha", ylab = "Precio")

plot(PETROLEOSP, main = "Precio histórico Petróleo", xlab = "Fecha", ylab = "Precio")

head(Crudo)
##        Fecha Precio
## 1 1996-01-03  17.40
## 2 1996-01-04  17.41
## 3 1996-01-05  17.70
## 4 1996-01-08  17.54
## 5 1996-01-09  17.41
## 6 1996-01-10  17.21
head(Acero)
##        Fecha Precio
## 1 2019-01-02  18.51
## 2 2019-01-03  18.48
## 3 2019-01-04  20.34
## 4 2019-01-07  20.45
## 5 2019-01-08  20.70
## 6 2019-01-09  20.84
Acero$Fecha <- as.Date(Acero$Fecha, format="01/01/2019")
Crudo$Fecha <- as.Date(Crudo$Fecha, format="01/01/2019")

petroleo_data_recortado <- Crudo %>% filter(Fecha >= as.Date("01/01/2019"))
acero_data_recortado <- Acero %>% filter(Fecha >= as.Date("01/01/2019"))

# Calcular el promedio de los precios
promedio_acero <- mean(Acero$Precio, na.rm=TRUE)
promedio_petroleo <- mean(Crudo$Precio, na.rm=TRUE)

# Reemplazar valores faltantes con el promedio
Acero$Precio[is.na(Acero$Precio)] <- promedio_acero
Crudo$Precio[is.na(Crudo$Precio)] <- promedio_petroleo

# Si se desean identificar y reemplazar valores atípicos, se puede usar el siguiente código:
# Definir una función para identificar valores atípicos
es_atipico <- function(x) {
  abs(x - mean(x, na.rm=TRUE)) > 3 * sd(x, na.rm=TRUE)
}
# Convertir "N/A" a NA
Acero$Precio[Acero$Precio == "N/A"] <- NA
Crudo$Precio[Crudo$Precio == "N/A"] <- NA

# Asegurarse de que las columnas de precios sean numéricas
Acero$Precio <- as.numeric(as.character(Acero$Precio))
Crudo$Precio <- as.numeric(as.character(Crudo$Precio))

# Calcular el promedio de los precios
promedio_acero <- mean(Acero$Precio, na.rm=TRUE)
promedio_petroleo <- mean(Crudo$Precio, na.rm=TRUE)

# Reemplazar valores atípicos con el promedio
Acero$Precio[es_atipico(Acero$Precio)] <- promedio_acero
Crudo$Precio[es_atipico(Crudo$Precio)] <- promedio_petroleo
# Graficar la serie temporal original de ACEROSP
plot(ACEROSP, type="l", main="Acero 2019 - 2024")

# Graficar la función de autocorrelación de ACEROSP
acf(ACEROSP)

# Realizar la prueba de Dickey-Fuller aumentada en ACEROSP
adf.test(ACEROSP)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag     ADF p.value
## [1,]   0  0.1173   0.678
## [2,]   1  0.0535   0.659
## [3,]   2  0.0272   0.652
## [4,]   3  0.0156   0.648
## [5,]   4 -0.0114   0.641
## [6,]   5 -0.0458   0.631
## [7,]   6 -0.0190   0.638
## [8,]   7 -0.0451   0.631
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -1.10   0.667
## [2,]   1 -1.19   0.635
## [3,]   2 -1.14   0.653
## [4,]   3 -1.15   0.647
## [5,]   4 -1.18   0.638
## [6,]   5 -1.23   0.622
## [7,]   6 -1.18   0.638
## [8,]   7 -1.22   0.624
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -2.59   0.327
## [2,]   1 -2.70   0.282
## [3,]   2 -2.77   0.251
## [4,]   3 -2.80   0.240
## [5,]   4 -2.85   0.216
## [6,]   5 -2.92   0.187
## [7,]   6 -2.89   0.203
## [8,]   7 -2.94   0.180
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
# Determinar el número de diferencias necesarias
ndiffs(ACEROSP)
## [1] 1
# Calcular la primera diferencia de ACEROSP
d3 <- diff(ACEROSP)

# Graficar la primera diferencia de ACEROSP
plot(d3, type="l", main="Primera Diferencia del Acero 2019 - 2024")

# Graficar la función de autocorrelación de la primera diferencia
acf(d3)

# Realizar la prueba de Dickey-Fuller aumentada en la primera diferencia
adf.test(d3)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -17.9    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -17.9    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -18.0    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
# Ajustar un modelo ARIMA a la primera diferencia
MRLSd3 <- auto.arima(d3)
summary(MRLSd3)
## Series: d3 
## ARIMA(0,0,0) with zero mean 
## 
## sigma^2 = 0.6386:  log likelihood = -1624.76
## AIC=3251.52   AICc=3251.52   BIC=3256.74
## 
## Training set error measures:
##                      ME      RMSE       MAE MPE MAPE      MASE       ACF1
## Training set 0.01368382 0.7991062 0.5189191 100  100 0.6368498 0.04857202
# Graficar la primera diferencia con la tendencia ajustada por el modelo ARIMA
plot(d3, type="l", main="Primera Diferencia del Acero 2019 - 2024")
lines(fitted(MRLSd3), col="red")

# Calcular los residuos del modelo ARIMA
d3_t <- residuals(MRLSd3)

# Graficar los residuos del modelo ARIMA
plot(d3_t, type="l", main="Primera Diferencia del Acero 2019 - 2024 sin tendencia")

# Graficar la función de autocorrelación de los residuos
acf(d3_t)

# Realizar la prueba de Dickey-Fuller aumentada en los residuos
adf.test(d3_t)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -17.9    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -17.9    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -18.0    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01

Pruebas de Estacionariedad

Augmented Dickey-Fuller (ADF)

library(tseries)

adf_acero <- adf.test(Acero$Precio)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag     ADF p.value
## [1,]   0  0.1173   0.678
## [2,]   1  0.0535   0.659
## [3,]   2  0.0272   0.652
## [4,]   3  0.0156   0.648
## [5,]   4 -0.0114   0.641
## [6,]   5 -0.0458   0.631
## [7,]   6 -0.0190   0.638
## [8,]   7 -0.0451   0.631
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -1.10   0.667
## [2,]   1 -1.19   0.635
## [3,]   2 -1.14   0.653
## [4,]   3 -1.15   0.647
## [5,]   4 -1.18   0.638
## [6,]   5 -1.23   0.622
## [7,]   6 -1.18   0.638
## [8,]   7 -1.22   0.624
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -2.59   0.327
## [2,]   1 -2.70   0.282
## [3,]   2 -2.77   0.251
## [4,]   3 -2.80   0.240
## [5,]   4 -2.85   0.216
## [6,]   5 -2.92   0.187
## [7,]   6 -2.89   0.203
## [8,]   7 -2.94   0.180
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
print(adf_acero)
## $type1
##      lag         ADF   p.value
## [1,]   0  0.11727193 0.6776460
## [2,]   1  0.05350881 0.6592981
## [3,]   2  0.02718279 0.6517227
## [4,]   3  0.01563846 0.6484008
## [5,]   4 -0.01141459 0.6406162
## [6,]   5 -0.04580448 0.6307205
## [7,]   6 -0.01897033 0.6384421
## [8,]   7 -0.04506030 0.6309346
## 
## $type2
##      lag       ADF   p.value
## [1,]   0 -1.097359 0.6673065
## [2,]   1 -1.188444 0.6350640
## [3,]   2 -1.137659 0.6530410
## [4,]   3 -1.153650 0.6473806
## [5,]   4 -1.179917 0.6380825
## [6,]   5 -1.225000 0.6221239
## [7,]   6 -1.180995 0.6377011
## [8,]   7 -1.219835 0.6239523
## 
## $type3
##      lag       ADF   p.value
## [1,]   0 -2.591050 0.3269262
## [2,]   1 -2.697910 0.2819326
## [3,]   2 -2.772046 0.2507175
## [4,]   3 -2.797702 0.2399148
## [5,]   4 -2.854646 0.2159386
## [6,]   5 -2.923674 0.1868740
## [7,]   6 -2.886523 0.2025165
## [8,]   7 -2.939848 0.1800638
adf_petroleo <- adf.test(Crudo$Precio)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag    ADF p.value
## [1,]   0 -2.092  0.0374
## [2,]   1 -1.420  0.1714
## [3,]   2 -1.094  0.2879
## [4,]   3 -0.830  0.3821
## [5,]   4 -0.703  0.4276
## [6,]   5 -0.618  0.4579
## [7,]   6 -0.564  0.4773
## [8,]   7 -0.517  0.4939
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -7.48  0.0100
## [2,]   1 -5.11  0.0100
## [3,]   2 -3.94  0.0100
## [4,]   3 -2.99  0.0381
## [5,]   4 -2.55  0.1095
## [6,]   5 -2.26  0.2246
## [7,]   6 -2.12  0.2808
## [8,]   7 -2.02  0.3210
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -7.88  0.0100
## [2,]   1 -5.39  0.0100
## [3,]   2 -4.16  0.0100
## [4,]   3 -3.17  0.0939
## [5,]   4 -2.69  0.2835
## [6,]   5 -2.39  0.4114
## [7,]   6 -2.24  0.4761
## [8,]   7 -2.12  0.5241
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
print(adf_petroleo)
## $type1
##      lag        ADF    p.value
## [1,]   0 -2.0915193 0.03736435
## [2,]   1 -1.4201150 0.17139301
## [3,]   2 -1.0940423 0.28785653
## [4,]   3 -0.8301153 0.38212343
## [5,]   4 -0.7028176 0.42759041
## [6,]   5 -0.6179412 0.45790577
## [7,]   6 -0.5637323 0.47726759
## [8,]   7 -0.5170827 0.49392944
## 
## $type2
##      lag       ADF    p.value
## [1,]   0 -7.476875 0.01000000
## [2,]   1 -5.110905 0.01000000
## [3,]   2 -3.937512 0.01000000
## [4,]   3 -2.994031 0.03806567
## [5,]   4 -2.546272 0.10949138
## [6,]   5 -2.258592 0.22456321
## [7,]   6 -2.118040 0.28078392
## [8,]   7 -2.017493 0.32100274
## 
## $type3
##      lag       ADF    p.value
## [1,]   0 -7.876712 0.01000000
## [2,]   1 -5.389902 0.01000000
## [3,]   2 -4.158146 0.01000000
## [4,]   3 -3.165340 0.09390502
## [5,]   4 -2.694266 0.28346681
## [6,]   5 -2.390471 0.41138067
## [7,]   6 -2.236736 0.47611119
## [8,]   7 -2.123425 0.52407663

Augmented Dickey-Fuller (ADF) diferencia

adf_diff_acero <- adf.test(diff(Acero$Precio))
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -17.9    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -17.9    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -35.1    0.01
## [2,]   1 -25.9    0.01
## [3,]   2 -21.0    0.01
## [4,]   3 -18.0    0.01
## [5,]   4 -15.7    0.01
## [6,]   5 -14.8    0.01
## [7,]   6 -13.4    0.01
## [8,]   7 -13.3    0.01
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
print(adf_diff_acero)
## $type1
##      lag       ADF p.value
## [1,]   0 -35.09090    0.01
## [2,]   1 -25.87921    0.01
## [3,]   2 -20.97032    0.01
## [4,]   3 -17.93927    0.01
## [5,]   4 -15.72716    0.01
## [6,]   5 -14.78377    0.01
## [7,]   6 -13.42327    0.01
## [8,]   7 -13.31492    0.01
## 
## $type2
##      lag       ADF p.value
## [1,]   0 -35.08784    0.01
## [2,]   1 -25.87837    0.01
## [3,]   2 -20.97193    0.01
## [4,]   3 -17.94220    0.01
## [5,]   4 -15.73107    0.01
## [6,]   5 -14.78883    0.01
## [7,]   6 -13.42912    0.01
## [8,]   7 -13.32266    0.01
## 
## $type3
##      lag       ADF p.value
## [1,]   0 -35.08660    0.01
## [2,]   1 -25.88582    0.01
## [3,]   2 -20.98268    0.01
## [4,]   3 -17.95613    0.01
## [5,]   4 -15.74725    0.01
## [6,]   5 -14.80772    0.01
## [7,]   6 -13.44952    0.01
## [8,]   7 -13.34528    0.01
adf_diff_petroleo <- adf.test(diff(Crudo$Precio))
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -54.4    0.01
## [2,]   1 -39.8    0.01
## [3,]   2 -35.6    0.01
## [4,]   3 -30.1    0.01
## [5,]   4 -26.3    0.01
## [6,]   5 -22.8    0.01
## [7,]   6 -20.3    0.01
## [8,]   7 -19.3    0.01
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -54.4    0.01
## [2,]   1 -39.8    0.01
## [3,]   2 -35.6    0.01
## [4,]   3 -30.1    0.01
## [5,]   4 -26.3    0.01
## [6,]   5 -22.8    0.01
## [7,]   6 -20.3    0.01
## [8,]   7 -19.2    0.01
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -54.3    0.01
## [2,]   1 -39.7    0.01
## [3,]   2 -35.6    0.01
## [4,]   3 -30.1    0.01
## [5,]   4 -26.3    0.01
## [6,]   5 -22.8    0.01
## [7,]   6 -20.3    0.01
## [8,]   7 -19.2    0.01
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
print(adf_diff_petroleo)
## $type1
##      lag       ADF p.value
## [1,]   0 -54.38621    0.01
## [2,]   1 -39.76721    0.01
## [3,]   2 -35.59793    0.01
## [4,]   3 -30.10996    0.01
## [5,]   4 -26.27168    0.01
## [6,]   5 -22.78835    0.01
## [7,]   6 -20.32788    0.01
## [8,]   7 -19.25354    0.01
## 
## $type2
##      lag       ADF p.value
## [1,]   0 -54.36619    0.01
## [2,]   1 -39.75255    0.01
## [3,]   2 -35.58479    0.01
## [4,]   3 -30.09883    0.01
## [5,]   4 -26.26194    0.01
## [6,]   5 -22.77990    0.01
## [7,]   6 -20.32037    0.01
## [8,]   7 -19.24648    0.01
## 
## $type3
##      lag       ADF p.value
## [1,]   0 -54.34615    0.01
## [2,]   1 -39.73788    0.01
## [3,]   2 -35.57164    0.01
## [4,]   3 -30.08768    0.01
## [5,]   4 -26.25220    0.01
## [6,]   5 -22.77137    0.01
## [7,]   6 -20.31260    0.01
## [8,]   7 -19.23897    0.01

Phillips-Perron (PP)

acf(Crudo)

acf(Acero)

PP.test(PETROLEOSP)
## 
##  Phillips-Perron Unit Root Test
## 
## data:  PETROLEOSP
## Dickey-Fuller = NA, Truncation lag parameter = 7, p-value = NA
PP.test(ACEROSP)
## 
##  Phillips-Perron Unit Root Test
## 
## data:  ACEROSP
## Dickey-Fuller = -2.7131, Truncation lag parameter = 7, p-value = 0.2764
pp.test(PETROLEOSP)
## Phillips-Perron Unit Root Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##  lag  Z_rho p.value
##    7 -0.471   0.587
## ----- 
##  Type 2: with drift no trend 
##  lag Z_rho p.value
##    7 -5.03   0.461
## ----- 
##  Type 3: with drift and trend 
##  lag Z_rho p.value
##    7 -5.65   0.759
## --------------- 
## Note: p-value = 0.01 means p.value <= 0.01
pp.test(ACEROSP)
## Phillips-Perron Unit Root Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##  lag  Z_rho p.value
##    7 0.0682   0.707
## ----- 
##  Type 2: with drift no trend 
##  lag Z_rho p.value
##    7 -3.69   0.576
## ----- 
##  Type 3: with drift and trend 
##  lag Z_rho p.value
##    7 -13.3   0.316
## --------------- 
## Note: p-value = 0.01 means p.value <= 0.01
PPt_Petroleo = ur.pp(PETROLEOSP)
PPt_Acero = ur.pp(ACEROSP)

Phillips-Perron (PP) diferencia

PP_diff_crudo <- pp.test(diff(Crudo$Precio))
## Phillips-Perron Unit Root Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##  lag Z_rho p.value
##    7 -1451    0.01
## ----- 
##  Type 2: with drift no trend 
##  lag Z_rho p.value
##    7 -1451    0.01
## ----- 
##  Type 3: with drift and trend 
##  lag Z_rho p.value
##    7 -1451    0.01
## --------------- 
## Note: p-value = 0.01 means p.value <= 0.01
print(PP_diff_crudo)
##        lag     Z_rho p.value
## type 1   7 -1451.369    0.01
## type 2   7 -1451.369    0.01
## type 3   7 -1451.369    0.01
PP_diff_acero <- pp.test(diff(Acero$Precio))
## Phillips-Perron Unit Root Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##  lag Z_rho p.value
##    7 -1301    0.01
## ----- 
##  Type 2: with drift no trend 
##  lag Z_rho p.value
##    7 -1300    0.01
## ----- 
##  Type 3: with drift and trend 
##  lag Z_rho p.value
##    7 -1299    0.01
## --------------- 
## Note: p-value = 0.01 means p.value <= 0.01
print(PP_diff_acero)
##        lag     Z_rho p.value
## type 1   7 -1300.649    0.01
## type 2   7 -1300.303    0.01
## type 3   7 -1299.465    0.01

Pruebas de Cointegración

Engle-Granger

crudo_tsm <- ts(Crudo$Precio, start = c(2019, 01), frequency = 12)
acero_tsm <- ts(Acero$Precio, start = c(2019, 01), frequency = 12)
ts_data <- cbind(crudo_tsm, acero_tsm)
eg_test <- ca.po(ts_data)
summary(eg_test)
## 
## ######################################## 
## # Phillips and Ouliaris Unit Root Test # 
## ######################################## 
## 
## Test of type Pu 
## detrending of series none 
## 
## 
## Call:
## lm(formula = z[, 1] ~ z[, -1] - 1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -18.085  -4.469   3.507   9.468  15.923 
## 
## Coefficients:
##         Estimate Std. Error t value Pr(>|t|)    
## z[, -1] 0.669582   0.009167   73.04   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 8.163 on 1360 degrees of freedom
## Multiple R-squared:  0.7969, Adjusted R-squared:  0.7967 
## F-statistic:  5335 on 1 and 1360 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: 25.8562 
## 
## Critical values of Pu are:
##                   10pct    5pct    1pct
## critical values 20.3933 25.9711 38.3413

Engle-Granger

# Cargar los datos
Crudo <- read_excel("Crudo_Datos.xlsx")
Acero <- read_excel("Acero_Datos.xlsx")

# Asegurarse de que las fechas están en el formato correcto y convertir a data.frame
Crudo <- data.frame(Fecha = as.Date(Crudo$Fecha), Precio = Crudo$Precio)
Acero <- data.frame(Fecha = as.Date(Acero$Date), Precio = Acero$Close)

# Fusionar los datos por fecha común
merged_data <- merge(Crudo, Acero, by.x = "Fecha", by.y = "Fecha", all = FALSE)

# Renombrar columnas para mayor claridad
colnames(merged_data) <- c("Fecha", "Precio_Crudo", "Precio_Acero")

# Ajuste del modelo de cointegración de Engle-Granger
eg_model <- lm(Precio_Crudo ~ Precio_Acero, data = merged_data)
residuals_eg <- residuals(eg_model)

# Prueba ADF en los residuos
adf_test_residuals <- adf.test(residuals_eg)
## Augmented Dickey-Fuller Test 
## alternative: stationary 
##  
## Type 1: no drift no trend 
##      lag   ADF p.value
## [1,]   0 -2.52  0.0126
## [2,]   1 -2.53  0.0119
## [3,]   2 -2.37  0.0191
## [4,]   3 -2.35  0.0198
## [5,]   4 -2.43  0.0163
## [6,]   5 -2.36  0.0192
## [7,]   6 -2.26  0.0237
## [8,]   7 -2.31  0.0216
## Type 2: with drift no trend 
##      lag   ADF p.value
## [1,]   0 -2.52   0.120
## [2,]   1 -2.53   0.115
## [3,]   2 -2.37   0.182
## [4,]   3 -2.35   0.187
## [5,]   4 -2.43   0.155
## [6,]   5 -2.36   0.183
## [7,]   6 -2.26   0.224
## [8,]   7 -2.31   0.204
## Type 3: with drift and trend 
##      lag   ADF p.value
## [1,]   0 -2.47   0.379
## [2,]   1 -2.49   0.371
## [3,]   2 -2.31   0.445
## [4,]   3 -2.30   0.450
## [5,]   4 -2.38   0.414
## [6,]   5 -2.32   0.441
## [7,]   6 -2.21   0.486
## [8,]   7 -2.26   0.465
## ---- 
## Note: in fact, p.value = 0.01 means p.value <= 0.01
print(adf_test_residuals)
## $type1
##      lag       ADF    p.value
## [1,]   0 -2.520130 0.01256584
## [2,]   1 -2.534568 0.01194708
## [3,]   2 -2.366549 0.01914789
## [4,]   3 -2.352280 0.01975945
## [5,]   4 -2.434053 0.01625489
## [6,]   5 -2.364212 0.01924806
## [7,]   6 -2.261050 0.02366927
## [8,]   7 -2.309791 0.02158040
## 
## $type2
##      lag       ADF   p.value
## [1,]   0 -2.519248 0.1203007
## [2,]   1 -2.533704 0.1145186
## [3,]   2 -2.365822 0.1816713
## [4,]   3 -2.351571 0.1873716
## [5,]   4 -2.433339 0.1546645
## [6,]   5 -2.363438 0.1826247
## [7,]   6 -2.260310 0.2238758
## [8,]   7 -2.309079 0.2043682
## 
## $type3
##      lag       ADF   p.value
## [1,]   0 -2.466764 0.3792574
## [2,]   1 -2.486924 0.3707688
## [3,]   2 -2.310829 0.4449141
## [4,]   3 -2.299239 0.4497941
## [5,]   4 -2.384379 0.4139456
## [6,]   5 -2.319595 0.4412230
## [7,]   6 -2.214121 0.4856332
## [8,]   7 -2.262828 0.4651251

Los datos muestran que los precios del acero y del petróleo han experimentado fluctuaciones significativas desde 2019. La prueba de Dickey-Fuller indicó que las series no eran estacionarias, por lo que se diferenciaron para realizar un análisis más robusto. Los modelos ARIMA ajustados a las series diferenciadas mostraron que ambas series tienen componentes de autoregresión y media móvil que explican sus comportamientos temporales. La prueba de cointegración de Engle-Granger sugiere que existe una relación de equilibrio a largo plazo entre los precios del acero y del petróleo, lo cual es crucial para estrategias de cobertura y planificación financiera a largo plazo.

Modelo de Corrección de Errores

Modelos ARIMA y Proyecciones

# Cargar los datos
Crudo <- read_excel("Crudo_Datos.xlsx")
Acero <- read_excel("Acero_Datos.xlsx")

# Asegurarse de que las fechas están en el formato correcto y convertir a data.frame
Crudo <- data.frame(Fecha = as.Date(Crudo$Fecha), Precio = Crudo$Precio)
Acero <- data.frame(Fecha = as.Date(Acero$Date), Precio = Acero$Close)

# Fusionar los datos por fecha común
merged_data <- merge(Crudo, Acero, by = "Fecha", all = FALSE)

# Renombrar columnas para mayor claridad
colnames(merged_data) <- c("Fecha", "Precio_Crudo", "Precio_Acero")

# Convertir a series temporales
ts_crudo <- ts(merged_data$Precio_Crudo, start = c(2019, 1), frequency = 12)
ts_acero <- ts(merged_data$Precio_Acero, start = c(2019, 1), frequency = 12)

# Ajustar el modelo ARIMA para los datos originales
auto_arima_crudo <- auto.arima(ts_crudo)
auto_arima_acero <- auto.arima(ts_acero)

Modelo AR(1)

# Ajustar modelo AR(1) para el precio del Crudo
ar1_crudo <- arima(Crudo$Precio, order = c(1, 0, 0))
summary(ar1_crudo)
## 
## Call:
## arima(x = Crudo$Precio, order = c(1, 0, 0))
## 
## Coefficients:
##          ar1  intercept
##       0.9990    52.5149
## s.e.  0.0005    13.5587
## 
## sigma^2 estimated as 1.711:  log likelihood = -11737.7,  aic = 23481.4
## 
## Training set error measures:
##                       ME     RMSE       MAE        MPE     MAPE      MASE
## Training set 0.005386238 1.307863 0.8049095 0.05659783 1.853651 0.9887456
##                    ACF1
## Training set 0.04380896
# Ajustar modelo AR(1) para el precio del acero
ar1_acero <- arima(Acero$Precio, order = c(1, 0, 0))
summary(ar1_acero)
## 
## Call:
## arima(x = Acero$Precio, order = c(1, 0, 0))
## 
## Coefficients:
##          ar1  intercept
##       0.9969    23.4728
## s.e.  0.0019     5.8038
## 
## sigma^2 estimated as 0.6376:  log likelihood = -1627.46,  aic = 3260.91
## 
## Training set error measures:
##                       ME      RMSE       MAE        MPE     MAPE     MASE
## Training set 0.008662804 0.7984932 0.5193006 -0.1426747 2.671228 1.000735
##                    ACF1
## Training set 0.05072827

Modelo MA(1)

# Modelo MA(1) para el precio del crudo
ma1_crudo <- arima(Crudo$Precio, order = c(0, 0, 1))
summary(ma1_crudo)
## 
## Call:
## arima(x = Crudo$Precio, order = c(0, 0, 1))
## 
## Coefficients:
##          ma1  intercept
##       0.9736    52.2306
## s.e.  0.0032     0.3571
## 
## sigma^2 estimated as 237.1:  log likelihood = -28877.7,  aic = 57761.4
## 
## Training set error measures:
##                      ME     RMSE      MAE       MPE   MAPE     MASE      ACF1
## Training set 0.09445574 15.39706 12.88013 -27.02645 42.272 15.82186 0.9051288
# Modelo MA(1) para el precio del acero
ma1_acero <- arima(Acero$Precio, order = c(0, 0, 1))
summary(ma1_acero)
## 
## Call:
## arima(x = Acero$Precio, order = c(0, 0, 1))
## 
## Coefficients:
##          ma1  intercept
##       0.9608    21.9548
## s.e.  0.0056     0.2791
## 
## sigma^2 estimated as 27.6:  log likelihood = -4190.29,  aic = 8386.58
## 
## Training set error measures:
##                         ME     RMSE      MAE       MPE     MAPE     MASE
## Training set -0.0005655064 5.253794 4.092059 -15.19846 27.18266 7.885735
##                   ACF1
## Training set 0.9035304

Modelo ARMA(1,1)

#Modelo ARMA(1,1) para el precio del crudo
arma11_crudo <- arima(Crudo$Precio, order = c(1, 0, 1))
summary(arma11_crudo)
## 
## Call:
## arima(x = Crudo$Precio, order = c(1, 0, 1))
## 
## Coefficients:
##       ar1     ma1  intercept
##         1  0.0310    52.5227
## s.e.    0  0.0125  5065.5854
## 
## sigma^2 estimated as 1.707:  log likelihood = -11734.87,  aic = 23477.74
## 
## Training set error measures:
##                       ME     RMSE      MAE      MPE     MAPE      MASE
## Training set 0.005590886 1.306399 0.803101 0.113994 1.848134 0.9865241
##                    ACF1
## Training set 0.01392774
#Modelo ARMA(1,1) para el precio del acero
arma11_acero <- arima(Acero$Precio, order = c(1, 0, 1))
## Warning in arima(Acero$Precio, order = c(1, 0, 1)): possible convergence
## problem: optim gave code = 1
summary(arma11_acero)
## 
## Call:
## arima(x = Acero$Precio, order = c(1, 0, 1))
## 
## Coefficients:
##          ar1     ma1  intercept
##       0.9967  0.0521    24.3862
## s.e.  0.0021  0.0276     5.9028
## 
## sigma^2 estimated as 0.6359:  log likelihood = -1625.68,  aic = 3259.35
## 
## Training set error measures:
##                     ME      RMSE       MAE        MPE     MAPE     MASE
## Training set 0.0049734 0.7974404 0.5189512 -0.1612654 2.673319 1.000062
##                       ACF1
## Training set -0.0002199933

Modelos de Volatilidad

Modelo GARCH

# Librerías necesarias


# Carga de datos
Crudo <- read_excel("Crudo_Datos.xlsx")
Acero <- read_excel("Acero_Datos.xlsx")

# Conversión a data.frames con fechas en formato Date
Crudo <- data.frame(Fecha = as.Date(Crudo$Fecha), Precio = Crudo$Precio)
Acero <- data.frame(Fecha = as.Date(Acero$Date), Precio = Acero$Close)

# Eliminación de filas con valores NA en ambas columnas
Crudo <- Crudo %>% drop_na()
Acero <- Acero %>% drop_na()

# Verificación de varianza
if (var(Crudo$Precio) == 0) {
  stop("La serie de precios del crudo es constante.")
}

if (var(Acero$Precio) == 0) {
  stop("La serie de precios del acero es constante.")
}

# Definición de las especificaciones del modelo GARCH
spec_garch <- ugarchspec(
  variance.model = list(model = "sGARCH", garchOrder = c(1, 1)),
  mean.model = list(armaOrder = c(1, 1), include.mean = TRUE),
  distribution.model = "norm"
)

# Ajuste del modelo GARCH para la serie de precios del crudo
garch_crudo <- ugarchfit(spec = spec_garch, data = Crudo$Precio)
summary(garch_crudo)
##    Length     Class      Mode 
##         1 uGARCHfit        S4
# Ajuste del modelo GARCH para la serie de precios del acero
garch_acero <- ugarchfit(spec = spec_garch, data = Acero$Precio)
summary(garch_acero)
##    Length     Class      Mode 
##         1 uGARCHfit        S4

Interpretación de los Resultados

1. Graficar las series temporales de precios históricos

Se han graficado los precios históricos del acero y del petróleo desde enero de 2019 hasta mayo de 2024. Estos gráficos permiten observar la tendencia general y las fluctuaciones en los precios a lo largo del tiempo. Por ejemplo, los precios del acero muestran una tendencia ascendente significativa con picos y caídas notables, especialmente alrededor de 2022, mientras que los precios del petróleo también muestran una alta volatilidad con un máximo en 2022 seguido de una tendencia descendente.

2. Análisis de autocorrelación y prueba de Dickey-Fuller aumentada

Autocorrelation Function (ACF): - La función de autocorrelación (ACF) para el acero y el petróleo muestra una alta correlación en los primeros rezagos, lo que sugiere que los valores de la serie están altamente correlacionados con sus propios valores en el pasado. Esto indica la presencia de una estructura de dependencia temporal significativa.

Augmented Dickey-Fuller Test (ADF): - La prueba ADF se utiliza para verificar la estacionariedad de las series temporales. Los resultados indican que ambas series, en su nivel original, no son estacionarias (es decir, tienen una tendencia subyacente). Tras la diferenciación, las series se vuelven estacionarias, lo cual es crucial para construir modelos predictivos fiables.

3. Cálculo de las diferencias y ajuste de un modelo ARIMA

Diferencias: - La primera diferencia de las series temporales se calculó para eliminar la tendencia y convertir las series en estacionarias. Esto se observa en los gráficos de las primeras diferencias, que muestran variaciones alrededor de una media constante.

Modelo ARIMA: - Se ajustaron modelos ARIMA a las series diferenciadas. Los modelos ARIMA combinan componentes de autoregresión (AR), integración (I) y media móvil (MA) para modelar la estructura temporal de las series.

Resultados del Modelo ARIMA: - Los modelos ARIMA ajustados mostraron un buen desempeño en la captura de la estructura temporal de las series de precios del acero y del petróleo. Las líneas ajustadas en los gráficos de las diferencias indican que los modelos capturan adecuadamente la tendencia subyacente.

4. Análisis de los modelos y proyecciones

Ajuste del Modelo ARIMA: - Los modelos ARIMA fueron evaluados en términos de sus coeficientes y errores estándar. Los residuos del modelo ARIMA fueron analizados para asegurar que no contienen patrones significativos, lo cual sugiere que el modelo ha capturado bien la estructura de la serie temporal.

Proyecciones: - A partir de los modelos ARIMA ajustados, se pueden realizar proyecciones futuras de los precios del acero y del petróleo. Esto es útil para la toma de decisiones estratégicas y la gestión del riesgo en los mercados financieros.

Conclusión

El análisis de series temporales realizado muestra que los precios del acero y del petróleo han experimentado una alta volatilidad en el período estudiado, con tendencias y fluctuaciones significativas. La aplicación de modelos ARIMA ha permitido capturar la estructura temporal de las series y realizar proyecciones futuras. Estos modelos y análisis son esenciales para la toma de decisiones informadas en contextos de inversión y comercio de commodities, así como para la planificación estratégica y la gestión de riesgos.