Actividad Repaso

Librerías y carga de datos

Instalar paquetes y llamar librerías

library(readxl)
library(dplyr)
library(tidyr)
library(ggplot2)
library(corrplot)
library(stringr)
library(plm)
library(psych)
library(ggthemes)
library(caret)
library(glmnet)
library(car)
library(regclass)
library(lmtest)
library(mctest)
library(caret)
library(lubridate)
library(forecast)
library(tseries)
library(urca)
library(reshape2)

Carga de datasets

exports <- read_excel("C:/Users/karla/Desktop/CONCENTRACION/Modulo_progra/inegi_exports_dataset.xlsx", 
    sheet = "exports")
ts_exports <- read_excel("C:/Users/karla/Desktop/CONCENTRACION/Modulo_progra/inegi_exports_dataset.xlsx", 
    sheet = "ts_exports")
data <- read_excel("C:/Users/karla/Desktop/CONCENTRACION/Modulo_progra/inegi_exports_dataset.xlsx", 
    sheet = "data")
fdi <- read_excel("C:/Users/karla/Desktop/CONCENTRACION/Modulo_progra/inegi_exports_dataset.xlsx", 
    sheet = "fdi")

Limpieza y tratamiento inicial de datos

Diseño y composicion del dataset tipo panel

exports <- exports %>%
  mutate(state = str_trim(state))

fdi <- fdi %>%
  mutate(state = str_trim(state))

data <- data %>%
  mutate(state = str_trim(state))

# EXPORTS
exports_long <- exports %>%
  pivot_longer(
    cols = starts_with("real_exports_"),
    names_to = "year",
    values_to = "real_exports"
  ) %>%
  mutate(
    year = str_extract(year, "\\d{4}") %>% as.numeric()
  ) %>%
  select(state, region, year, real_exports)

# FDI
fdi_long <- fdi %>%
  pivot_longer(
    cols = starts_with("fdi_"),
    names_to = "year",
    values_to = "fdi"
  ) %>%
  mutate(
    year = str_extract(year, "\\d{4}") %>% as.numeric()
  ) %>%
  select(state, region, year, fdi)

# Verificar duplicados en exports
dup_exports <- exports_long %>%
  count(state, year) %>%
  filter(n > 1)

# Verificar duplicados en fdi
dup_fdi <- fdi_long %>%
  count(state, year) %>%
  filter(n > 1)

# Verificar duplicados en base principal
dup_data <- data %>%
  count(state, year) %>%
  filter(n > 1)

# Si cualquiera tiene filas, hay problema estructural
print(dup_exports)

## # A tibble: 0 × 3
## # ℹ 3 variables: state <chr>, year <dbl>, n <int>

print(dup_fdi)

## # A tibble: 0 × 3
## # ℹ 3 variables: state <chr>, year <dbl>, n <int>

print(dup_data)

## # A tibble: 0 × 3
## # ℹ 3 variables: state <chr>, year <dbl>, n <int>

ts_exports_annual <- ts_exports %>%
  mutate(
    year = substr(date, 1, 4) %>% as.numeric()
  ) %>%
  group_by(year) %>%
  summarise(
    ts_exports = sum(ts_exports, na.rm = TRUE),
    .groups = "drop"
  )

panel_data <- data %>%
  left_join(exports_long, by = c("state", "year")) %>%
  left_join(fdi_long %>% select(-region), 
            by = c("state", "year")) %>%
  left_join(ts_exports_annual, by = "year")

Variables de estudio:

Para propósitos de este análisis nos enfocaremos en las siguientes variables clave:

• real_exports: total de exportaciones no procedentes del petroleo o sus derivados a valor del dinero del 2018.
• fdi: Inversión extranjera directa con un valor del dinero relativo a 2018.
• year: año en el que sucedió el evento.
• real_public_investment: Inversión pública por la población del estado.
• state: Estado de la república.
• border_distance: distrancia entre la capital del estado y la frontera con Estados Unidos más cercana.

Estadística Descriptiva

psych::describe(panel_data)

##                           vars   n         mean           sd      median
## state*                       1 288        16.50         9.25       16.50
## year                         2 288      2019.00         2.59     2019.00
## pop_density                  3 288       308.82      1063.93       65.02
## gdp_per_capita_2018          4 288      1986.91      1062.05     1788.31
## lq_primary                   5 288         1.06         0.88        0.84
## lq_secondary                 6 288         1.00         0.38        0.97
## lq_tertiary                  7 288         1.00         0.09        1.00
## average_daily_salary         8 288       332.49        46.96      324.12
## real_public_investment_pc    9 288       627.38       446.39      506.42
## border_economic_activity    10 288        -1.77         0.87       -1.99
## crime_rate                  11 288        27.44        24.26       18.34
## college_education           12 256         0.25         0.05        0.25
## exchange_rate               13 288        19.20         1.28       19.47
## border_distance             14 288       704.92       274.21      751.64
## inpc                        15 288       108.19        14.01      105.93
## region*                     16 256         4.09         1.51        4.00
## real_exports                17 256 228822083.78 274331648.02 97535120.79
## fdi                         18 256     18535.61     26141.08    10478.15
## ts_exports                  19 288   7822529.33    606507.99  7997670.57
##                                trimmed          mad        min          max
## state*                           16.50        11.86       1.00 3.200000e+01
## year                           2019.00         2.97    2015.00 2.023000e+03
## pop_density                      96.05        67.14       9.60 6.233410e+03
## gdp_per_capita_2018            1831.36       728.79     603.71 8.216560e+03
## lq_primary                        0.92         0.64       0.01 4.630000e+00
## lq_secondary                      0.97         0.45       0.37 1.990000e+00
## lq_tertiary                       1.00         0.10       0.79 1.180000e+00
## average_daily_salary            328.56        43.29     249.97 5.057300e+02
## real_public_investment_pc       570.97       384.18       4.14 2.423360e+03
## border_economic_activity         -1.92         0.26      -2.77 2.530000e+00
## crime_rate                       23.32        14.07       1.98 1.169500e+02
## college_education                 0.25         0.05       0.15 4.400000e-01
## exchange_rate                    19.32         0.95      16.89 2.052000e+01
## border_distance                 720.74       223.45       8.83 1.252660e+03
## inpc                            107.59        16.87      89.05 1.323700e+02
## region*                           4.16         1.48       1.00 6.000000e+00
## real_exports              178098052.86 136535980.71  268782.64 1.163117e+09
## fdi                           13034.20     11111.04   -7412.68 1.712389e+05
## ts_exports                  7860733.19    558612.44 6535587.76 8.792925e+06
##                                  range  skew kurtosis          se
## state*                    3.100000e+01  0.00    -1.21        0.55
## year                      8.000000e+00  0.00    -1.24        0.15
## pop_density               6.223810e+03  5.20    25.70       62.69
## gdp_per_capita_2018       7.612850e+03  2.49     9.25       62.58
## lq_primary                4.620000e+00  1.46     1.96        0.05
## lq_secondary              1.620000e+00  0.42    -0.58        0.02
## lq_tertiary               4.000000e-01 -0.13    -0.57        0.01
## average_daily_salary      2.557600e+02  0.83     0.55        2.77
## real_public_investment_pc 2.419220e+03  1.32     2.08       26.30
## border_economic_activity  5.300000e+00  3.60    14.55        0.05
## crime_rate                1.149700e+02  1.50     1.69        1.43
## college_education         2.900000e-01  0.50     0.37        0.00
## exchange_rate             3.630000e+00 -0.86    -0.71        0.08
## border_distance           1.243830e+03 -0.54     0.04       16.16
## inpc                      4.332000e+01  0.35    -1.06        0.83
## region*                   5.000000e+00 -0.32    -1.14        0.09
## real_exports              1.162848e+09  1.38     0.99 17145728.00
## fdi                       1.786516e+05  3.26    12.47     1633.82
## ts_exports                2.257337e+06 -0.56     0.08    35738.83

panel_data$region <- as.factor(panel_data$region)
panel_data$state <- as.factor(panel_data$state)
summary(panel_data)

##                  state          year       pop_density      
##  Aguascalientes     :  9   Min.   :2015   Min.   :   9.601  
##  Baja California    :  9   1st Qu.:2017   1st Qu.:  42.067  
##  Baja California Sur:  9   Median :2019   Median :  65.017  
##  Campeche           :  9   Mean   :2019   Mean   : 308.816  
##  Chiapas            :  9   3rd Qu.:2021   3rd Qu.: 161.291  
##  Chihuahua          :  9   Max.   :2023   Max.   :6233.409  
##  (Other)            :234                                    
##  gdp_per_capita_2018   lq_primary       lq_secondary     lq_tertiary    
##  Min.   : 603.7      Min.   :0.01493   Min.   :0.3685   Min.   :0.7867  
##  1st Qu.:1307.2      1st Qu.:0.43199   1st Qu.:0.6470   1st Qu.:0.9360  
##  Median :1788.3      Median :0.83619   Median :0.9693   Median :0.9993  
##  Mean   :1986.9      Mean   :1.05887   Mean   :0.9969   Mean   :0.9973  
##  3rd Qu.:2314.7      3rd Qu.:1.32423   3rd Qu.:1.2569   3rd Qu.:1.0669  
##  Max.   :8216.6      Max.   :4.63480   Max.   :1.9911   Max.   :1.1847  
##                                                                         
##  average_daily_salary real_public_investment_pc border_economic_activity
##  Min.   :250.0        Min.   :   4.139          Min.   :-2.772          
##  1st Qu.:297.5        1st Qu.: 289.064          1st Qu.:-2.137          
##  Median :324.1        Median : 506.424          Median :-1.985          
##  Mean   :332.5        Mean   : 627.380          Mean   :-1.765          
##  3rd Qu.:360.2        3rd Qu.: 879.463          3rd Qu.:-1.791          
##  Max.   :505.7        Max.   :2423.360          Max.   : 2.530          
##                                                                         
##    crime_rate      college_education exchange_rate   border_distance  
##  Min.   :  1.984   Min.   :0.1453    Min.   :16.89   Min.   :   8.83  
##  1st Qu.: 10.270   1st Qu.:0.2077    1st Qu.:19.11   1st Qu.: 613.26  
##  Median : 18.335   Median :0.2460    Median :19.47   Median : 751.64  
##  Mean   : 27.443   Mean   :0.2482    Mean   :19.20   Mean   : 704.92  
##  3rd Qu.: 37.212   3rd Qu.:0.2798    3rd Qu.:20.11   3rd Qu.: 875.76  
##  Max.   :116.950   Max.   :0.4376    Max.   :20.52   Max.   :1252.66  
##                    NA's   :32                                         
##       inpc                       region    real_exports            fdi        
##  Min.   : 89.05   CdMx              : 8   Min.   :2.688e+05   Min.   : -7413  
##  1st Qu.: 98.27   Centro_Sur_Oriente:48   1st Qu.:2.152e+07   1st Qu.:  4048  
##  Median :105.93   Noreste           :32   Median :9.754e+07   Median : 10478  
##  Mean   :108.19   Noroeste          :48   Mean   :2.288e+08   Mean   : 18536  
##  3rd Qu.:117.31   Occidente_Bajio   :64   3rd Qu.:3.502e+08   3rd Qu.: 22572  
##  Max.   :132.37   Sur               :56   Max.   :1.163e+09   Max.   :171239  
##                   NA's              :32   NA's   :32          NA's   :32      
##    ts_exports     
##  Min.   :6535588  
##  1st Qu.:7588156  
##  Median :7997671  
##  Mean   :7822529  
##  3rd Qu.:8098039  
##  Max.   :8792925  
## 

Resultados

El análisis descriptivo arrojó los siguientes resultados:

• real_exports:

• Media: 228822083.78
  

• Mediana: 97535120.79

• Máximo: 1163116946.59
  

• Mínimo: 268782.64

• Desviación estándar: 274331648.02

• fdi:

• Media: 18535.61
  

• Mediana: 10478.15
  

• Máximo: 171238.91
  

• Mínimo: -7412.68

• Desviación estándar: 26141.08

• year:

• Media: 2019.50
  

• Mediana: 2019.50

• Máximo: 2023.00
  

• Mínimo: 2016.00

• Desviación estándar: 2.30

• real_public_investment_pc:

• Media: 599.27
  

• Mediana: 483.94
  

• Máximo: 2423.36
  

• Mínimo: 4.14

• Desviación estándar: 423.74

• border_distance:

• Media: 704.92
  

• Mediana: 751.64
  

• Máximo: 1252.66
  

• Mínimo: 8.83

• Desviación estándar: 274.34

Interpretaciones

• real_exports: La media es mucho mayor que la mediana, lo que indica que pocos estados exportan montos muy altos y elevan el promedio.

• fdi: Existe mucha variación entre estados y algunos concentran grandes montos de inversión, mientras que otros incluso presentan valores negativos.

• year: Los datos están distribuidos de forma uniforme entre 2016 y 2023, lo que confirma un panel balanceado en el tiempo.

• real_public_investment_pc: Hay diferencias importantes entre estados en inversión pública por habitante, con algunos valores considerablemente más altos que el resto.

• border_distance: La distancia a la frontera varía bastante entre estados, aunque la mayoría se concentra en distancias intermedias o altas.

Variable de mayor dispersión: real_exports.

panel_data <- panel_data[-(1:32), , drop = FALSE]
str(panel_data)

## tibble [256 × 19] (S3: tbl_df/tbl/data.frame)
##  $ state                    : Factor w/ 32 levels "Aguascalientes",..: 1 2 3 4 5 6 7 8 9 10 ...
##  $ year                     : num [1:256] 2016 2016 2016 2016 2016 ...
##  $ pop_density              : num [1:256] 235.01 50.85 9.79 14.92 70.76 ...
##  $ gdp_per_capita_2018      : num [1:256] 2329 2281 2254 7731 713 ...
##  $ lq_primary               : num [1:256] 0.403 0.845 0.599 0.486 0.731 ...
##  $ lq_secondary             : num [1:256] 1.191 1.508 0.651 0.709 0.965 ...
##  $ lq_tertiary              : num [1:256] 0.992 0.866 1.097 1.091 1.024 ...
##  $ average_daily_salary     : num [1:256] 314 333 318 424 307 ...
##  $ real_public_investment_pc: num [1:256] 1027 275 1360 887 1022 ...
##  $ border_economic_activity : num [1:256] -1.85 2.42 -2.16 -2.15 -2.36 ...
##  $ crime_rate               : num [1:256] 3.71 31.68 32.89 10.84 10.64 ...
##  $ college_education        : num [1:256] 0.249 0.258 0.287 0.229 0.146 ...
##  $ exchange_rate            : num [1:256] 20.5 20.5 20.5 20.5 20.5 ...
##  $ border_distance          : num [1:256] 625.59 8.83 800.32 978.33 1111.82 ...
##  $ inpc                     : num [1:256] 92 92 92 92 92 ...
##  $ region                   : Factor w/ 6 levels "CdMx","Centro_Sur_Oriente",..: 5 4 4 6 6 4 1 3 5 4 ...
##  $ real_exports             : num [1:256] 1.75e+08 8.03e+08 5.30e+06 2.12e+08 1.31e+07 ...
##  $ fdi                      : num [1:256] 13225 33796 11420 3162 3091 ...
##  $ ts_exports               : num [1:256] 7399443 7399443 7399443 7399443 7399443 ...

Visualizacion de los Datos

Histogramas

variables <- c("year","real_public_investment_pc","fdi","real_exports","border_distance","inpc")

panel_data %>%
  pivot_longer(cols = all_of(variables), names_to = "variable", values_to = "valor") %>%
  ggplot(aes(x = valor)) +
  geom_histogram(bins = 20) +
  facet_wrap(~ variable, scales = "free") +
  labs(title = "Distribucion de variables clave",
       x = "Valor de exportaciones", y = "Frecuencia") +
  theme_fivethirtyeight()

## Interpretación de histogramas

En términos generales, la mayoría de las variables presentan distribuciones asimétricas (sesgo positivo).
fdi, real_exports y real_public_investment_pc muestran una clara asimetría a la derecha: la mayor parte de las observaciones se concentran en valores bajos o medios, pero existen algunos valores muy altos que alargan la cola derecha.
inpc y year son relativamente más simétricas y estables.

Respecto a los outliers, las colas largas en los histogramas, evidencian valores atípicos principalmente en:
* fdi * real_exports * real_public_investment

Boxplots

panel_data %>%
  pivot_longer(cols = all_of(variables), names_to = "variable", values_to = "valor") %>%
  ggplot(aes(x = region, y = valor)) +
  geom_boxplot(outlier.alpha = 0.6) +
  facet_wrap(~ variable, scales = "free_y") +
  labs(title = "Dispersion de variables clave por region (boxplots)",
       x = "Region", y = "Valor de exportaciones") +
  theme_fivethirtyeight() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Interpretación de los boxplots

Observando los boxplots:

• La región Norte muestra la mayor dispersión, especialmente en real_exports y fdi, ya que tiene cajas más amplias y valores más extremos.
• La región Sur presenta la menor dispersión, con cajas más compactas y menor variabilidad en exportaciones.

¿Es relevante la región?

Sí. La región parece ser una variable relevante para analizar el desempeño exportador de México, ya que existen diferencias claras en niveles y dispersión de exportaciones entre regiones, particularmente entre el Norte y el Sur.

Matriz de Correlacion

corr_mat <- panel_data %>%
  select(all_of(variables)) %>%
  cor(use = "complete.obs")

corr_df <- melt(corr_mat)

ggplot(corr_df, aes(Var1, Var2, fill = value)) +
  geom_tile(color = "white") +
  geom_text(aes(label = round(value, 2)), size = 3) +
  scale_fill_gradient2(low = "#3B4CC0",
                       mid = "white",
                       high = "#B40426",
                       midpoint = 0,
                       limits = c(-1,1),
                       name = "Correlacion") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        panel.grid = element_blank()) +
  labs(title = "Matriz de correlacion: Variables Clave",
       x = NULL, y = NULL)

Interpretación de la Matriz de Correlacion

Las variables con mayor correlación son:

• year e inpc (0.98): correlación positiva muy fuerte.
• real_exports y border_distance (-0.75): correlación negativa fuerte.

Las variables con menor correlación (cercanas a cero) son:

• real_exports y year (0.03)
• real_exports e inpc (0.03)
• real_public_investment_pc y border_distance (0.05)
• fdi y year (-0.06)

Gráficos de dispersion

scatter_data <- panel_data %>%
  select(real_exports, all_of(variables)) %>%
  pivot_longer(
    cols = -real_exports,
    names_to = "variable",
    values_to = "value"
  )

ggplot(scatter_data, aes(x = value, y = real_exports)) +
  geom_point(alpha = 0.5, color = "#2C3E50") +
  geom_smooth(method = "lm", se = FALSE, color = "#E74C3C") +
  facet_wrap(~ variable, scales = "free_x")

Interpretación de gráficos de dispersión

A simple vista y similar a lo observado en la matriz de correlación, la única variable con la que parece tener una correlación los suficientement fuerte o significativa es con border_distance.

Preguntas clave

1. ¿Qué tan relevante es la distancia a la frontera con Estados Unidos en el crecimiento futuro de las exportaciones estatales?
2. ¿En qué medida la inversión extranjera directa (fdi) impulsa el desempeño exportador de los estados en los próximos años?
   
3. ¿La estabilidad de precios (inpc) influye directamente en la evolución de las exportaciones?

Hipótesis

Hipótesis 1: Los estados con menor distancia a la frontera con Estados Unidos presentan mayores niveles y tasas de crecimiento en sus exportaciones en comparación con estados más alejados.

Hipótesis 2: Un mayor flujo de inversión extranjera directa (fdi) se asocia positivamente con un mejor desempeño exportador estatal en los próximos años.

Hipótesis 3: Una mayor estabilidad de precios (menor variabilidad en el inpc) favorece el crecimiento sostenido de las exportaciones estatales.

Modelos de regresion (MLRM, MLRM + log(), LASSO)

Modelo base de regresion lineal múltiple (1 explicativa + 1 control)

modelo_base <- lm(real_exports ~ fdi + inpc, data = panel_data)
rmse_base <- sqrt(mean(residuals(modelo_base)^2, na.rm = TRUE))
summary(modelo_base)

## 
## Call:
## lm(formula = real_exports ~ fdi + inpc, data = panel_data)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -637143530 -164696862 -107592203  100346822  894220473 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 105723561  143428770   0.737    0.462    
## fdi              2958        635   4.658 5.16e-06 ***
## inpc           617390    1275907   0.484    0.629    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 264300000 on 253 degrees of freedom
## Multiple R-squared:  0.07904,    Adjusted R-squared:  0.07176 
## F-statistic: 10.86 on 2 and 253 DF,  p-value: 2.995e-05

rmse_base

## [1] 262751920

El modelo de regresión lineal múltiple estima las exportaciones reales en función de fdi e inpc. Los resultados indican que fdi tiene un efecto positivo y estadísticamente significativo, lo que sugiere que un aumento en la inversión extranjera directa se asocia con un incremento en las exportaciones, manteniendo constante el nivel de precios.

En contraste, inpc no resulta estadísticamente significativo, por lo que no se puede afirmar que tenga un efecto sistemático sobre las exportaciones en este modelo.

Sin embargo, el R2 ajustado es bajo, lo que indica que el modelo explica aproximadamente el 7% de la variabilidad en las exportaciones, por lo que existen otros factores relevantes no incluidos. El error estándar residual es de 264,300,000, reflejando una dispersión considerable en los datos.

MLRM con 2–4 controles (explicativa + controles)

modelo_full <- lm(real_exports ~ state + fdi + real_public_investment_pc + border_distance + inpc + year, data = panel_data)
rmse_full <- sqrt(mean(residuals(modelo_full)^2, na.rm = TRUE))
summary(modelo_full)

## 
## Call:
## lm(formula = real_exports ~ state + fdi + real_public_investment_pc + 
##     border_distance + inpc + year, data = panel_data)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -113742023   -8899649     874857    9831552  158504436 
## 
## Coefficients: (1 not defined because of singularities)
##                             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               -2.250e+09  1.088e+10  -0.207 0.836330    
## stateBaja California       5.825e+08  1.663e+07  35.020  < 2e-16 ***
## stateBaja California Sur  -1.787e+08  1.582e+07 -11.295  < 2e-16 ***
## stateCampeche              7.367e+07  1.620e+07   4.546 9.01e-06 ***
## stateChiapas              -1.633e+08  1.622e+07 -10.068  < 2e-16 ***
## stateChihuahua             8.182e+08  1.675e+07  48.847  < 2e-16 ***
## stateCiudad de Mexico     -1.715e+08  3.291e+07  -5.213 4.27e-07 ***
## stateCoahuila              6.612e+08  1.618e+07  40.878  < 2e-16 ***
## stateColima               -1.665e+08  1.626e+07 -10.238  < 2e-16 ***
## stateDurango              -1.378e+08  1.595e+07  -8.642 1.14e-15 ***
## stateGuanajuato            2.808e+08  1.626e+07  17.272  < 2e-16 ***
## stateGuerrero             -1.623e+08  1.633e+07  -9.936  < 2e-16 ***
## stateHidalgo              -1.385e+08  1.628e+07  -8.511 2.70e-15 ***
## stateJalisco               1.986e+08  1.660e+07  11.964  < 2e-16 ***
## stateMexico                1.410e+08  1.736e+07   8.120 3.32e-14 ***
## stateMichoacan            -1.026e+08  1.637e+07  -6.268 1.90e-09 ***
## stateMorelos              -1.197e+08  1.617e+07  -7.403 2.80e-12 ***
## stateNayarit              -1.756e+08  1.610e+07 -10.901  < 2e-16 ***
## stateNuevo Leon            5.296e+08  1.972e+07  26.858  < 2e-16 ***
## stateOaxaca               -1.616e+08  1.639e+07  -9.860  < 2e-16 ***
## statePuebla                9.394e+07  1.612e+07   5.827 1.99e-08 ***
## stateQueretaro             4.504e+07  1.583e+07   2.845 0.004863 ** 
## stateQuintana Roo         -1.779e+08  1.620e+07 -10.980  < 2e-16 ***
## stateSan Luis Potosi       7.629e+07  1.583e+07   4.819 2.69e-06 ***
## stateSinaloa              -1.450e+08  1.585e+07  -9.150  < 2e-16 ***
## stateSonora                1.732e+08  1.611e+07  10.751  < 2e-16 ***
## stateTabasco              -5.498e+07  1.613e+07  -3.409 0.000775 ***
## stateTamaulipas            3.282e+08  1.590e+07  20.640  < 2e-16 ***
## stateTlaxcala             -1.530e+08  1.601e+07  -9.556  < 2e-16 ***
## stateVeracruz             -6.786e+07  1.632e+07  -4.157 4.62e-05 ***
## stateYucatan              -1.541e+08  1.670e+07  -9.226  < 2e-16 ***
## stateZacatecas            -1.207e+08  1.608e+07  -7.506 1.50e-12 ***
## fdi                        2.727e+02  2.397e+02   1.138 0.256485    
## real_public_investment_pc  1.004e+04  6.582e+03   1.525 0.128773    
## border_distance                   NA         NA      NA       NA    
## inpc                       6.167e+04  9.482e+05   0.065 0.948199    
## year                       1.196e+06  5.437e+06   0.220 0.826073    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 31620000 on 220 degrees of freedom
## Multiple R-squared:  0.9885, Adjusted R-squared:  0.9867 
## F-statistic: 542.1 on 35 and 220 DF,  p-value: < 2.2e-16

rmse_full

## [1] 29312002

Al incorporar efectos fijos por estado, el ajuste del modelo mejora sustancialmente, alcanzando un R2 ajustado de 0.9867, lo que implica que aproximadamente el 98.7% de la variación en real_exports es explicada por el modelo. Esto sugiere que las diferencias estructurales entre estados capturan gran parte de la variabilidad en las exportaciones. El error estándar residual es de 31,620,000, indicando una reducción importante en la dispersión respecto al modelo base.

Sin embargo, las variables macroeconómicas pierden significancia estadística. En particular, fdi no resulta significativa (p = 0.256), al igual que real_public_investment_pc (p = 0.129), inpc (p = 0.948) y year (p = 0.826). Esto indica que, una vez controlando por efectos específicos de cada estado, estas variables no presentan un efecto sistemático adicional sobre las exportaciones.

En síntesis, el alto poder explicativo del modelo proviene principalmente de los efectos fijos estatales, más que de las variables económicas incluidas.

MLRM con escala logarítmica con 2–4 controles (explicativa + controles)

modelo_log <- lm(real_exports ~ state + log(1+fdi) + log(real_public_investment_pc) + log(border_distance) + log(inpc) + year, data = panel_data)
rmse_log <- sqrt(mean(residuals(modelo_log)^2, na.rm = TRUE))
summary(modelo_log)

## 
## Call:
## lm(formula = real_exports ~ state + log(1 + fdi) + log(real_public_investment_pc) + 
##     log(border_distance) + log(inpc) + year, data = panel_data)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -114197179   -9293080     107110    9653188  157965861 
## 
## Coefficients: (1 not defined because of singularities)
##                                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                     1.262e+10  1.344e+10   0.939 0.348540    
## stateBaja California            5.841e+08  1.668e+07  35.013  < 2e-16 ***
## stateBaja California Sur       -1.782e+08  1.608e+07 -11.081  < 2e-16 ***
## stateCampeche                   7.199e+07  1.643e+07   4.381 1.84e-05 ***
## stateChiapas                   -1.692e+08  1.657e+07 -10.211  < 2e-16 ***
## stateChihuahua                  8.200e+08  1.677e+07  48.884  < 2e-16 ***
## stateCiudad de Mexico          -1.357e+08  1.723e+07  -7.875 1.64e-13 ***
## stateCoahuila                   6.651e+08  1.614e+07  41.199  < 2e-16 ***
## stateColima                    -1.719e+08  1.751e+07  -9.820  < 2e-16 ***
## stateDurango                   -1.406e+08  1.610e+07  -8.730 7.00e-16 ***
## stateGuanajuato                 2.819e+08  1.629e+07  17.306  < 2e-16 ***
## stateGuerrero                  -1.671e+08  1.694e+07  -9.867  < 2e-16 ***
## stateHidalgo                   -1.434e+08  1.650e+07  -8.692 9.00e-16 ***
## stateJalisco                    2.033e+08  1.634e+07  12.442  < 2e-16 ***
## stateMexico                     1.506e+08  1.638e+07   9.195  < 2e-16 ***
## stateMichoacan                 -1.077e+08  1.654e+07  -6.509 5.20e-10 ***
## stateMorelos                   -1.243e+08  1.641e+07  -7.576 1.03e-12 ***
## stateNayarit                   -1.805e+08  1.626e+07 -11.104  < 2e-16 ***
## stateNuevo Leon                 5.420e+08  1.760e+07  30.799  < 2e-16 ***
## stateOaxaca                    -1.681e+08  1.721e+07  -9.766  < 2e-16 ***
## statePuebla                     9.283e+07  1.637e+07   5.671 4.52e-08 ***
## stateQueretaro                  4.681e+07  1.610e+07   2.908 0.004014 ** 
## stateQuintana Roo              -1.824e+08  1.626e+07 -11.215  < 2e-16 ***
## stateSan Luis Potosi            7.831e+07  1.604e+07   4.882 2.04e-06 ***
## stateSinaloa                   -1.460e+08  1.601e+07  -9.122  < 2e-16 ***
## stateSonora                     1.707e+08  1.630e+07  10.473  < 2e-16 ***
## stateTabasco                   -6.001e+07  1.626e+07  -3.690 0.000284 ***
## stateTamaulipas                 3.303e+08  1.608e+07  20.541  < 2e-16 ***
## stateTlaxcala                  -1.571e+08  1.623e+07  -9.684  < 2e-16 ***
## stateVeracruz                  -7.066e+07  1.662e+07  -4.251 3.17e-05 ***
## stateYucatan                   -1.589e+08  1.791e+07  -8.871 2.78e-16 ***
## stateZacatecas                 -1.265e+08  1.675e+07  -7.554 1.18e-12 ***
## log(1 + fdi)                   -7.205e+05  2.481e+06  -0.290 0.771779    
## log(real_public_investment_pc)  3.059e+06  3.473e+06   0.881 0.379443    
## log(border_distance)                   NA         NA      NA       NA    
## log(inpc)                       1.491e+08  1.370e+08   1.088 0.277886    
## year                           -6.512e+06  6.967e+06  -0.935 0.350976    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 32010000 on 216 degrees of freedom
##   (4 observations deleted due to missingness)
## Multiple R-squared:  0.9884, Adjusted R-squared:  0.9865 
## F-statistic: 524.4 on 35 and 216 DF,  p-value: < 2.2e-16

rmse_log

## [1] 29632525

En el modelo con variables en escala logarítmica, el ajuste se mantiene muy alto (R2 ajustado = 0.9865), lo que indica que casi toda la variación en real_exports sigue siendo explicada principalmente por los efectos por estado.

Sin embargo, ninguna de las variables transformadas resulta significativa (log(fdi), log(inversión pública), log(inpc) y year), por lo que no se evidencia un efecto estadísticamente relevante adicional una vez controlando por las diferencias estructurales entre estados.

LASSO

set.seed(123)
idx <- createDataPartition(panel_data$real_exports, p = 0.75, list = FALSE)
train_data <- panel_data[idx, ]
test_data  <- panel_data[-idx, ]

# --- TRAIN: x e y alineados (mismas filas) ---
mf_train <- model.frame(real_exports ~ ., data = train_data, na.action = na.omit)
y_train  <- model.response(mf_train)
x_train  <- model.matrix(real_exports ~ ., data = mf_train)[, -1]

# --- TEST: x e y alineados (mismas filas) ---
mf_test <- model.frame(real_exports ~ ., data = test_data, na.action = na.omit)
y_test  <- model.response(mf_test)
x_test  <- model.matrix(real_exports ~ ., data = mf_test)[, -1]

set.seed(123)
cv_lasso <- cv.glmnet(x_train, y_train, alpha = 1)
lasso_model <- glmnet(x_train, y_train, alpha = 1, lambda = cv_lasso$lambda.min)

pred_lasso <- as.vector(predict(lasso_model, newx = x_test))

RMSE_lasso <- RMSE(pred_lasso, y_test)
R2_lasso   <- R2(pred_lasso, y_test)

RMSE_lasso

## [1] 32077848

R2_lasso

## [1] 0.9858817

coef(lasso_model)

## 53 x 1 sparse Matrix of class "dgCMatrix"
##                                      s0
## (Intercept)                2.434832e+10
## stateBaja California       1.646669e+08
## stateBaja California Sur  -9.417562e+07
## stateCampeche              1.740530e+08
## stateChiapas               .           
## stateChihuahua             7.184962e+08
## stateCiudad de Mexico     -1.519623e+08
## stateCoahuila              5.336884e+08
## stateColima               -1.150422e+07
## stateDurango              -8.313408e+07
## stateGuanajuato            3.148007e+08
## stateGuerrero             -1.339583e+07
## stateHidalgo              -1.034491e+08
## stateJalisco               2.590956e+08
## stateMexico                2.084174e+08
## stateMichoacan            -3.059761e+07
## stateMorelos              -4.885840e+07
## stateNayarit              -5.774609e+07
## stateNuevo Leon            3.214024e+08
## stateOaxaca               -7.960520e+06
## statePuebla                1.724382e+08
## stateQueretaro             2.445335e+07
## stateQuintana Roo          6.838772e+07
## stateSan Luis Potosi       .           
## stateSinaloa              -3.404641e+07
## stateSonora                5.810535e+07
## stateTabasco               7.354367e+07
## stateTamaulipas            1.734331e+08
## stateTlaxcala             -1.107964e+08
## stateVeracruz             -3.546283e+07
## stateYucatan              -4.593071e+07
## stateZacatecas            -1.023260e+08
## year                      -1.221207e+07
## pop_density                .           
## gdp_per_capita_2018       -5.817087e+03
## lq_primary                -7.187881e+06
## lq_secondary               2.631034e+07
## lq_tertiary                .           
## average_daily_salary       5.493408e+05
## real_public_investment_pc  4.211319e+03
## border_economic_activity   6.272532e+07
## crime_rate                -1.670094e+05
## college_education         -1.463661e+07
## exchange_rate              2.700079e+06
## border_distance           -2.245104e+05
## inpc                       .           
## regionCentro_Sur_Oriente   .           
## regionNoreste              2.964960e+07
## regionNoroeste             .           
## regionOccidente_Bajio      .           
## regionSur                  .           
## fdi                        1.080623e+02
## ts_exports                 6.138800e+01

Como modelo alternativo, el LASSO presenta un buen desempeño predictivo en el conjunto de prueba, con RMSE = 32,077,848 y R2 = 0.9859, lo que indica una alta capacidad para predecir real_exports.

Además, se observa su principal ventaja metodológica: realiza selección automática de variables, dejando varios coeficientes en cero, lo que reduce complejidad y posible sobreajuste.

Entre los coeficientes retenidos predominan los efectos por estado, reforzando la idea de que las diferencias estructurales entre entidades explican gran parte del nivel de exportaciones. En cuanto a variables numéricas, fdi se mantiene con efecto positivo, mientras que border_distance aparece con efecto negativo. Otras variables como inpc son eliminadas del modelo, lo que sugiere que no aportan poder predictivo adicional una vez capturados los efectos estructurales y demás controles.

Tabla Comparativa de los Modelos

resultados <- data.frame(
  Regression_Model = c("a) Modelo Base", "b) Modelo Full", "c) Modelo log", "d) LASSO"),
  Adjusted_R2      = c(summary(modelo_base)$adj.r.squared,
                       summary(modelo_full)$adj.r.squared,
                       summary(modelo_log)$adj.r.squared,
                       R2_lasso),
  RMSE             = c(rmse_base, rmse_full, rmse_log, RMSE_lasso)
)

resultados

##   Regression_Model Adjusted_R2      RMSE
## 1   a) Modelo Base  0.07176185 262751920
## 2   b) Modelo Full  0.98671516  29312002
## 3    c) Modelo log  0.98648286  29632525
## 4         d) LASSO  0.98588168  32077848

El Modelo Base tiene muy bajo ajuste (R2 = 0.07) y alto error, por lo que explica poco las exportaciones. Al incluir controles y efectos por estado (Modelo Full y Log), el ajuste aumenta drásticamente (R2= 0.99) y el RMSE baja considerablemente, mostrando que las diferencias estructurales entre estados explican casi toda la variación.

El LASSO mantiene un desempeño alto y competitivo, aunque con un RMSE ligeramente mayor, aportando además selección automática de variables.

Diagnosticos

Multicolinealidad

VIF(modelo_base)

##     fdi    inpc 
## 1.00578 1.00578

Los VIF de fdi e inpc son aproximadamente 1, lo que indica ausencia de multicolinealidad relevante en el modelo base. En términos prácticos, estas variables no presentan una correlación lineal fuerte entre sí que distorsione las estimaciones de los coeficientes.

Heterocedasticidad

bptest(modelo_base)

## 
##  studentized Breusch-Pagan test
## 
## data:  modelo_base
## BP = 45.532, df = 2, p-value = 1.297e-10

bptest(modelo_full)

## 
##  studentized Breusch-Pagan test
## 
## data:  modelo_full
## BP = 88.661, df = 35, p-value = 1.511e-06

bptest(modelo_log)

## 
##  studentized Breusch-Pagan test
## 
## data:  modelo_log
## BP = 87.2, df = 35, p-value = 2.407e-06

Heterocedasticidad (modelo_base): El test de Breusch–Pagan arroja BP = 45.532, df = 2 y p-value = 1.297e-10. Dado que el p-value es mucho menor a 0.05, se rechaza la hipótesis de homocedasticidad, confirmando la presencia de heterocedasticidad; es decir, la varianza de los errores no es constante.

Heterocedasticidad (modelo_full): Para este modelo, el test reporta BP = 88.661, df = 35 y p-value = 1.511e-06. Nuevamente, al ser p < 0.05, se concluye que persiste la heterocedasticidad incluso al incorporar efectos por estado y más controles.

Heterocedasticidad (modelo_log): En el modelo con transformación logarítmica, el test muestra BP = 87.2, df = 35 y p-value = 2.407e-06. El resultado sigue siendo significativo, por lo que la transformación logarítmica no corrige el problema; la varianza del error continúa dependiendo de las variables explicativas.

Normalidad de Residuales

hist(modelo_base$residuals,xlab="Residuales estimados de la regresion", main="Distribucion de Residuales en Modelo Base", col="lightgreen", border="white")

hist(modelo_full$residuals,xlab="Residuales estimados de la regresion", main="Distribucion de Residuales en Modelo Full", col="lightgreen", border="white")

hist(modelo_log$residuals,xlab="Residuales estimados de la regresion", main="Distribucion de Residuales en Modelo logaritmico", col="lightgreen", border="white")

Modelo Base (normalidad de residuales): El histograma muestra una distribución poco simétrica, con una cola derecha pronunciada y presencia de valores extremos. Esto sugiere que los residuales no siguen una distribución normal y que el modelo base deja una parte importante de la variabilidad sin explicar, consistente con su bajo R² y alto RMSE.

Modelo Full (normalidad de residuales): Los residuales se concentran mucho más alrededor de cero y la forma se aproxima a una distribución más simétrica. Aunque persisten algunas colas y outliers, la normalidad mejora notablemente respecto al modelo base, en línea con su mejor ajuste.

Modelo Log (normalidad de residuales): El patrón es similar al modelo full: los residuales están centrados en cero y muestran mayor simetría que en el modelo base. No obstante, todavía se observan colas y valores extremos. La transformación logarítmica no modifica sustancialmente la normalidad en comparación con el modelo full.

library(tseries)

# Jarque-Bera test
jb_test_base <- jarque.bera.test(modelo_base$residuals)
jb_test_base

## 
##  Jarque Bera Test
## 
## data:  modelo_base$residuals
## X-squared = 74.541, df = 2, p-value < 2.2e-16

jb_test_full <- jarque.bera.test(modelo_full$residuals)
jb_test_full

## 
##  Jarque Bera Test
## 
## data:  modelo_full$residuals
## X-squared = 304.1, df = 2, p-value < 2.2e-16

jb_test_log <- jarque.bera.test(modelo_log$residuals)
jb_test_log

## 
##  Jarque Bera Test
## 
## data:  modelo_log$residuals
## X-squared = 265.67, df = 2, p-value < 2.2e-16

En los tres modelos, el test Jarque–Bera rechaza claramente la hipótesis de normalidad de los residuales. Para el modelo_base, se obtiene X² = 74.541 con p < 2.2e-16; en el modelo_full, X² = 304.1 con p < 2.2e-16; y en el modelo_log, X2 = 265.67 con p < 2.2e-16.

En todos los casos, el p-value extremadamente bajo indica que los residuales no siguen una distribución normal. Esto implica que, aunque visualmente los histogramas del modelo full y log parezcan más simétricos que el del modelo base, estadísticamente la normalidad sigue siendo rechazada, posiblemente debido a colas pesadas o presencia de outliers.

AIC(modelo_base, modelo_full, modelo_log)

##             df       AIC
## modelo_base  4 10660.498
## modelo_full 37  9603.572
## modelo_log  37  9460.154

Aquí los resultados son más claros. Dado que un AIC más bajo es preferible, el orden de los modelos queda así:

modelo_log (9460.154) < modelo_full (9603.572) < modelo_base (10660.498).

Esto sugiere que, considerando simultáneamente el ajuste y la complejidad del modelo, el modelo_log ofrece el mejor balance, seguido por el modelo_full, mientras que el modelo_base queda claramente rezagado en términos de desempeño relativo.

Análisis de series del tiempo

Limpieza de datos

ts_df <- ts_exports %>%
  mutate(
    year  = as.integer(substr(date, 1, 4)),
    month = as.integer(substr(date, 6, 7)),
    date_m = as.Date(paste0(year, "-", month, "-01"))
  ) %>%
  arrange(date_m)

Visualización de series del tiempo

p_ts <- ggplot(ts_df, aes(x = date_m, y = ts_exports)) +
  geom_line() +
  theme_minimal() +
  labs(
    title = "Exportaciones a traves del tiempo",
    x = "Date",
    y = "ts_exports"
  )

p_ts

La serie de tiempo de exportaciones muestra una tendencia claramente creciente de largo plazo, especialmente a partir de mediados de los años 2000. Esto sugiere un crecimiento estructural sostenido en el nivel de exportaciones a lo largo del tiempo.

Se observan también fluctuaciones cíclicas y episodios de alta volatilidad, particularmente alrededor de 2008–2009 y 2020, donde se presentan caídas abruptas seguidas de recuperaciones, lo cual podría asociarse a choques económicos externos.

Además, la variabilidad parece aumentar conforme crece el nivel de la serie, lo que sugiere posible heterocedasticidad y no estacionariedad en nivel, indicando que podría ser necesario diferenciar la serie o trabajar en logaritmos para análisis más formales de modelado temporal.

Descomposicion de series del tiempo

ts_y <- ts(
  ts_df$ts_exports,
  start = c(min(ts_df$year), min(ts_df$month[ts_df$year == min(ts_df$year)])),
  frequency = 12
)

decomp_classic <- decompose(ts_y, type = "additive")
plot(decomp_classic)

La descomposición aditiva muestra que la serie de exportaciones está compuesta por tres elementos principales: tendencia, componente estacional y componente aleatorio.

La tendencia presenta un crecimiento sostenido de largo plazo, especialmente a partir de los años 2000, lo que confirma que el aumento observado en la serie no es temporal sino estructural.

El componente estacional es bastante marcado y repetitivo a lo largo del tiempo, lo que indica un patrón mensual estable en las exportaciones (probablemente relacionado con ciclos productivos o comerciales).

Por último, el componente aleatorio muestra fluctuaciones relativamente moderadas en la mayor parte del período, aunque destaca un choque muy fuerte alrededor de 2020, evidenciando un evento atípico que afectó significativamente la serie.

En conjunto, la serie presenta tendencia creciente y estacionalidad clara, con algunos shocks puntuales que alteran temporalmente su comportamiento.

Test de estacionariedad y autocorrelacion serial

adf_level  <- tseries::adf.test(ts_y)
adf_level

## 
##  Augmented Dickey-Fuller Test
## 
## data:  ts_y
## Dickey-Fuller = -2.7908, Lag order = 7, p-value = 0.2433
## alternative hypothesis: stationary

Acf(ts_y, main = "ACF: ts_exports")

Pacf(ts_y, main = "PACF: ts_exports")

El test ADF (p = 0.2433) indica que la serie no es estacionaria en nivel, ya que no se rechaza la presencia de raíz unitaria. Esto es consistente con la ACF, que muestra autocorrelaciones altas y decrecimiento lento, reflejando fuerte persistencia y tendencia.

La PACF presenta un pico importante en el rezago 1 y señales alrededor del lag 12, lo que sugiere componente autorregresivo de corto plazo y estacionalidad anual. En conjunto, la serie requiere diferenciación regular y posiblemente estacional para poder modelarse adecuadamente.

ts_y_d1 <- diff(ts_y, differences = 1)

modelo_arma <- Arima(
  ts_y_d1,
  order = c(1,0,1),
  include.constant = FALSE
)

summary(modelo_arma)

## Series: ts_y_d1 
## ARIMA(1,0,1) with zero mean 
## 
## Coefficients:
##          ar1      ma1
##       0.1397  -0.7685
## s.e.  0.0648   0.0388
## 
## sigma^2 = 1.969e+09:  log likelihood = -4931.87
## AIC=9869.73   AICc=9869.79   BIC=9881.76
## 
## Training set error measures:
##                    ME     RMSE      MAE      MPE     MAPE      MASE
## Training set 6095.921 44259.69 28757.12 11.34027 211.5769 0.8993409
##                      ACF1
## Training set -0.005080321

checkresiduals(modelo_arma)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,0,1) with zero mean
## Q* = 195.82, df = 22, p-value < 2.2e-16
## 
## Model df: 2.   Total lags used: 24

Después de aplicar una diferenciación para corregir la no estacionariedad, se estimó un ARIMA(1,0,1) que logra capturar parte de la dinámica de corto plazo de la serie. Los coeficientes son significativos, lo que indica que existe dependencia temporal entre periodos consecutivos.

Sin embargo, el test de Ljung-Box resulta significativo, lo que implica que todavía hay autocorrelación en los residuales. En otras palabras, el modelo no logra explicar completamente la estructura temporal de la serie.

En conclusión, aunque la diferenciación mejora el ajuste respecto a trabajar en niveles, el modelo aún puede perfeccionarse, posiblemente incorporando un componente estacional (SARIMA) para capturar mejor los patrones mensuales.

modelo_arima <- Arima(
  ts_y,
  order = c(1,1,1),
  include.constant = TRUE
)

summary(modelo_arima)

## Series: ts_y 
## ARIMA(1,1,1) with drift 
## 
## Coefficients:
##          ar1      ma1      drift
##       0.1647  -0.8148  1652.2861
## s.e.  0.0649   0.0382   485.8929
## 
## sigma^2 = 1.928e+09:  log likelihood = -4927.24
## AIC=9862.47   AICc=9862.57   BIC=9878.51
## 
## Training set error measures:
##                    ME     RMSE      MAE        MPE     MAPE      MASE
## Training set 1.093527 43697.63 28231.01 -0.9572928 7.338796 0.6850747
##                    ACF1
## Training set 0.01418405

checkresiduals(modelo_arima)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,1,1) with drift
## Q* = 185.54, df = 22, p-value < 2.2e-16
## 
## Model df: 2.   Total lags used: 24

El ARIMA(1,1,1) con drift muestra que los coeficientes AR(1) y MA(1) son estadísticamente significativos, lo que confirma que existe dependencia de corto plazo en la serie. El término drift positivo (1652) indica una tendencia promedio creciente en el tiempo, coherente con el comportamiento observado en la serie original.

En términos de ajuste, el modelo mejora ligeramente respecto al anterior (AIC = 9862.47, menor que antes), lo que sugiere una mejor especificación.

Sin embargo, el test de Ljung-Box sigue siendo significativo (p < 2.2e-16), lo que indica que aún existe autocorrelación en los residuales. Es decir, aunque el modelo captura parte de la dinámica y la tendencia, no elimina completamente la dependencia temporal, por lo que podría requerirse un modelo con componente estacional (SARIMA).

modelo_sarima <- Arima(
  ts_y,
  order = c(1,1,1),
  seasonal = c(1,1,1),
  include.constant = TRUE
)

summary(modelo_sarima)

## Series: ts_y 
## ARIMA(1,1,1)(1,1,1)[12] 
## 
## Coefficients:
##          ar1      ma1    sar1     sma1
##       0.1122  -0.7140  0.0104  -0.8390
## s.e.  0.0757   0.0523  0.0585   0.0306
## 
## sigma^2 = 1.369e+09:  log likelihood = -4720.77
## AIC=9451.55   AICc=9451.7   BIC=9471.44
## 
## Training set error measures:
##                     ME     RMSE      MAE        MPE    MAPE      MASE
## Training set -97.56791 36215.02 22759.48 -0.4223074 5.68124 0.5522986
##                    ACF1
## Training set 0.01070795

checkresiduals(modelo_sarima)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,1,1)(1,1,1)[12]
## Q* = 150.49, df = 20, p-value < 2.2e-16
## 
## Model df: 4.   Total lags used: 24

El SARIMA incorpora tanto la tendencia como la estacionalidad anual detectada previamente, lo que le permite capturar mejor la dinámica completa de la serie. Presenta un AIC más bajo y menor RMSE en comparación con los modelos ARIMA anteriores, indicando un mejor balance entre ajuste y complejidad.

Además, los coeficientes estacionales reflejan que existe una dependencia sistemática cada 12 periodos, coherente con el patrón mensual observado. Sin embargo, el test de Ljung-Box sigue siendo significativo, lo que sugiere que aún queda cierta autocorrelación en los residuales. En conjunto, es el modelo más adecuado hasta ahora, aunque todavía podría perfeccionarse.

library(forecast)

fc_5 <- forecast(modelo_sarima, h = (5*12))

fc_5

##          Point Forecast    Lo 80     Hi 80    Lo 95     Hi 95
## Jan 2026       626786.3 579376.6  674196.0 554279.4  699293.1
## Feb 2026       701484.9 650453.8  752516.0 623439.6  779530.2
## Mar 2026       770569.9 717184.5  823955.3 688924.0  852215.8
## Apr 2026       714768.1 659228.2  770308.0 629827.1  799709.0
## May 2026       727821.7 670218.5  785424.9 639725.2  815918.3
## Jun 2026       738255.1 678661.1  797849.1 647114.0  829396.3
## Jul 2026       731630.7 670110.4  793150.9 637543.6  825717.8
## Aug 2026       751347.6 687959.6  814735.6 654404.1  848291.1
## Sep 2026       749738.1 684535.9  814940.3 650019.9  849456.3
## Oct 2026       802541.5 735574.2  869508.8 700123.8  904959.2
## Nov 2026       753991.8 685304.7  822678.9 648943.9  859039.6
## Dec 2026       762075.9 691711.0  832440.7 654462.1  869689.6
## Jan 2027       643327.4 569173.6  717481.2 529919.0  756735.8
## Feb 2027       718024.4 641595.3  794453.5 601136.1  834912.6
## Mar 2027       786930.0 708418.7  865441.3 666857.4  907002.7
## Apr 2027       731014.2 650488.3  811540.1 607860.4  854168.0
## May 2027       744186.0 661696.1  826675.9 618028.6  870343.4
## Jun 2027       755060.9 670652.9  839468.8 625970.0  884151.7
## Jul 2027       748059.9 661776.4  834343.3 616100.8  880019.0
## Aug 2027       768174.3 680055.3  856293.2 633407.9  902940.6
## Sep 2027       766603.5 676686.5  856520.6 629087.3  904119.8
## Oct 2027       818583.8 726904.0  910263.7 678371.6  958796.1
## Nov 2027       770978.1 677568.7  864387.5 628120.7  913835.5
## Dec 2027       778723.7 683616.2  873831.2 633269.3  924178.1
## Jan 2028       659963.4 561468.0  758458.8 509327.7  810599.2
## Feb 2028       734660.4 633958.2  835362.5 580649.8  888671.0
## Mar 2028       803564.2 700807.6  906320.7 646411.6  960716.7
## Apr 2028       747647.2 642888.0  852406.4 587431.8  907862.6
## May 2028       760820.2 654097.1  867543.3 597601.3  924039.0
## Jun 2028       771699.6 663048.3  880350.9 605531.7  937867.5
## Jul 2028       764694.7 654148.8  875240.6 595629.3  933760.1
## Aug 2028       784813.2 672404.6  897221.8 612899.1  956727.3
## Sep 2028       783242.9 669002.0  897483.8 608526.5  957959.3
## Oct 2028       835214.7 719170.4  951258.9 657740.3 1012689.1
## Nov 2028       787618.7 669798.7  905438.8 607428.5  967809.0
## Dec 2028       795360.8 675791.4  914930.3 612495.1  978226.5
## Jan 2029       676600.4 553775.0  799425.8 488755.2  864445.6
## Feb 2029       751297.4 626253.3  876341.5 560059.0  942535.8
## Mar 2029       820201.1 693072.1  947330.2 625774.0 1014628.3
## Apr 2029       764284.1 635114.1  893454.2 566735.6  961832.7
## May 2029       777457.2 646279.0  908635.3 576837.5  978076.8
## Jun 2029       788336.6 655180.8  921492.4 584692.4  991980.9
## Jul 2029       781331.7 646227.2  916436.2 574707.2  987956.2
## Aug 2029       801450.3 664424.7  938475.8 591887.8 1011012.7
## Sep 2029       799879.9 660960.0  938799.9 587420.2 1012339.7
## Oct 2029       851851.6 711062.7  992640.5 636533.5 1067169.7
## Nov 2029       804255.8 661622.4  946889.2 586116.9 1022394.7
## Dec 2029       811997.8 667543.5  956452.1 591074.1 1032921.6
## Jan 2030       693237.4 545582.1  840892.7 467418.1  919056.7
## Feb 2030       767934.4 618018.5  917850.3 538657.8  997210.9
## Mar 2030       836838.1 684783.2  988893.1 604290.2 1069386.1
## Apr 2030       780921.1 626766.5  935075.8 545161.9 1016680.4
## May 2030       794094.2 637869.0  950319.3 555168.4 1033019.9
## Jun 2030       804973.6 646705.3  963242.0 562923.0 1047024.3
## Jul 2030       797968.7 637683.1  958254.3 552833.0 1043104.4
## Aug 2030       818087.3 655809.5  980365.0 569904.9 1066269.6
## Sep 2030       816516.9 652271.2  980762.6 565324.8 1067709.0
## Oct 2030       868488.6 702298.3 1034679.0 614322.4 1122654.8
## Nov 2030       820892.8 652780.3  989005.3 563786.9 1077998.7
## Dec 2030       828634.8 658621.9  998647.8 568622.5 1088647.2

autoplot(fc_5) +
  theme_minimal() +
  labs(
    title = "SARIMA Forecast: Proximos 5 periodos",
    x = "Time",
    y = "ts_exports"
  )

El pronóstico generado con el SARIMA proyecta las exportaciones para los próximos 5 años. La trayectoria esperada mantiene una tendencia creciente con patrón estacional, coherente con el comportamiento histórico de la serie.

Los valores puntuales muestran que las exportaciones seguirían aumentando gradualmente, manteniendo fluctuaciones mensuales.

Las bandas azules representan los intervalos de confianza al 80% y 95%, que se amplían conforme avanza el horizonte temporal. Esto refleja el aumento de la incertidumbre en el largo plazo, algo típico en modelos de series de tiempo.

En resumen, el modelo proyecta crecimiento sostenido con estacionalidad, aunque con mayor variabilidad e incertidumbre conforme se extiende el pronóstico.