Caso 2. Extracción de señales y Pronóstico ARIMA para el sector de vehículos

La variable “ventas de vehículos” es un indicador clave del dinamismo económico, tanto a nivel local como nacional, en particular, su análisis desagregado en Cali, el departamento del Valle del Cauca y Colombia en general, permite comprender patrones de consumo, capacidad adquisitiva, ciclos económicos y tendencias del mercado automotor en diferentes escalas territoriales.

En el caso de Cali, como principal núcleo urbano del suroccidente colombiano, las ventas de vehículos reflejan el comportamiento del sector terciario y la evolución del ingreso disponible. A nivel departamental (Valle), permiten identificar el impacto de factores regionales como la infraestructura vial, las dinámicas agroindustriales y la conectividad intermunicipal. Finalmente, en el ámbito nacional, estas ventas se ven afectadas por variables macroeconómicas como la tasa de interés, el dólar, la inflación y la confianza del consumidor.

Además, este tipo de variable es especialmente útil por su frecuencia mensual y alta disponibilidad estadística, lo que favorece el uso de modelos econométricos y análisis de series temporales, como estos modelos ARIMA o SARIMA, para evaluar tendencias, prever escenarios y orientar decisiones tanto del sector público como del privado.

logo1

#Cargar librerías necesarias
library(readxl)  # Para leer archivos Excel
library(tseries)  # Para pruebas de estacionariedad
library(forecast)  # Para modelado ARIMA y pronósticos
library(ggplot2)  # Para visualización de datos
library(plotly)  # Para gráficos interactivos
library(timetk)   #timetk simplifica y acelera el análisis exploratorio, visualización, y preparación de datos temporales para modelado. Es ideal para quienes trabajan con series temporales en un flujo de trabajo "tidy" y buscan integrar análisis visuales, detección de patrones y forecasting en un solo paquete.

Cargar base de datos

library(readxl)
BASE_CASO_2_AJUSTADA_VEHICULOS <- read_excel("C:/Users/Jazmin/OneDrive/ESPECIALIZACION EN FINANZAS/I SEMESTRE/ANALITICA DE NEGOCIOS/CASO 2/BASE CASO 2 AJUSTADA_VEHICULOS.xlsx", 
    col_types = c("date", "numeric", "numeric", 
        "numeric"))
View(BASE_CASO_2_AJUSTADA_VEHICULOS)

PASO INDISPENSABLE: Declarar la (s) variable (s) como serie (s) temporal (es):

VEH_COLOMBIA Variable 1

# Convertir/declarar =VEH_COLOMBIA en serie de tiempo mensual
VEH_COLOMBIA1_ts <- ts(BASE_CASO_2_AJUSTADA_VEHICULOS$VEH_COLOMBIA, start = c(2014, 1), frequency = 12)

VEH_VALLE Variable 2

# Convertir/declarar VEH_VALLE en serie de tiempo  mensual
VEH_VALLE2_ts <- ts(BASE_CASO_2_AJUSTADA_VEHICULOS$VEH_VALLE, start = c(2014, 1), frequency = 12)

VEH_CALI Variable 3

# Convertir/declarar VEH_CALI en serie de tiempo mensual
VEH_CALI3_ts <- ts(BASE_CASO_2_AJUSTADA_VEHICULOS$VEH_CALI, start = c(2014, 1), frequency = 12)

Extracción de señales

Gráfico inicial de la variable 1 en niveles -Original

library(ggplot2)
library(plotly)

# Convertir la serie temporal a un vector numérico para lograr graficar con ggplot2
BASE_CASO_2_AJUSTADA_VEHICULOS$VEH_COLOMBIA1_ts <- as.numeric(VEH_COLOMBIA1_ts)

# Crear el gráfico
grafico_serie <- ggplot(BASE_CASO_2_AJUSTADA_VEHICULOS, aes(x = seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = nrow(BASE_CASO_2_AJUSTADA_VEHICULOS)), 
                                      y = VEH_COLOMBIA1_ts)) +
  geom_line(color = "grey", linewidth = 0.4) +  # Cambiado 'size' por 'linewidth'
  geom_point(color = "black", size = 0.1) +
  ggtitle("VEH_COLOMBIA Variable 1: Serie original") +
  xlab("Tiempo") +
  ylab("VARIABLE VEH_COLOMBIA") +
  theme_minimal()

ggplotly(grafico_serie)

Extracción señales variable 1 VH_COLOMBIA

# Cargar librerías necesarias
library(ggplot2)
library(plotly)

# Descomposición de la serie temporal
stl_decomp_VEH_COLOMBIA1_ts <- stl(VEH_COLOMBIA1_ts, s.window = "periodic")

# Convertir la descomposición a un data frame para graficar con ggplot2
stl_df_var1 <- data.frame(
  Time = rep(time(VEH_COLOMBIA1_ts), 12),  # Tiempo repetido para cada componente (son 12 componentes)
  Value = c(stl_decomp_VEH_COLOMBIA1_ts$time.series[, "seasonal"], 
            stl_decomp_VEH_COLOMBIA1_ts$time.series[, "trend"], 
            stl_decomp_VEH_COLOMBIA1_ts$time.series[, "remainder"], 
            VEH_COLOMBIA1_ts),
  Component = rep(c("Estacional", "Tendencia", "Residuo", "Serie Original"), each = length(VEH_COLOMBIA1_ts))
)

# Crear gráfico con ggplot2
 stl_df_var1 <- ggplot(stl_df_var1, aes(x = Time, y = Value, color = Component)) +
  geom_line() +
  facet_wrap(~Component, scales = "free_y", ncol = 1) + 
  theme_minimal() +
  labs(title = "Descomposición temporal de la variable VEH_COLOMBIA",
       x = "Tiempo",
       y = "Valor")

# Convertir a gráfico interactivo con plotly
ggplotly(stl_df_var1)

interpretación varibale escogida

La variable escogida VEH_COLOMBIA representa la cantidad de vehículos nuevos registrados mensualmente en Colombia, es un dato clave para el comportamiento del sector automotriz, reflejando la demanda, el consumo y el dinamismo económico.

El gráfico de la serie original muestra un crecimiento sostenido entre 2014 y 2019, una caída abrupta en 2020 por la pandemia del COVID 19 y recuperación progresiva posterior, este patrón refleja la vulnerabilidad del sector ante choques externos (como COVID 19), pero también su capacidad de recuperación.

Descomposición STL

• Estacionalidad: Se observa un comportamiento cíclico con picos regulares, probablemente por dinámicas de consumo, por ejemplo como las promociones de fin de año y aumento en el consumo de temporada de fin de año.

• Residuo: Las variaciones aleatorias confirman eventos no explicados por la tendencia o la estacionalidad, como políticas públicas o condiciones crediticias.

• Tendencia: Es claramente sostenido con interrupciones en 2020. Esto indica una mejora estructural en el sector, aunque sujeta a ciclos económicos.

Relación con la empresa (ejemplo: Concesionario de vehículos):

Al plantear hipotéticamente un concesionario de vehículos, entender estas fluctuaciones le permite:

• Planificar stock.

• Ajustar campañas promocionales según estacionalidad.

• Prepararse para desaceleraciones anticipadas por shocks.

Factores clave clave: • La señal de estacionalidad identificada podría explotarse comercialmente. • La recuperación post 2020 sugiere un entorno de inversión atractivo.

Extracción señales variable 2 VEH_VALLE

# Cargar librerías necesarias
library(ggplot2)
library(plotly)

# Descomposición de la serie temporal
stl_decomp_VEH_VALLE2_ts <- stl(VEH_VALLE2_ts, s.window = "periodic")

# Convertir la descomposición a un data frame para graficar con ggplot2
stl_df_var2 <- data.frame(
  Time = rep(time(VEH_VALLE2_ts), 12),  # Tiempo repetido para cada componente (son 12 componentes)
  Value = c(stl_decomp_VEH_VALLE2_ts$time.series[, "seasonal"], 
            stl_decomp_VEH_VALLE2_ts$time.series[, "trend"], 
            stl_decomp_VEH_VALLE2_ts$time.series[, "remainder"], 
            VEH_VALLE2_ts),
  Component = rep(c("Estacional", "Tendencia", "Residuo", "Serie Original"), each = length(VEH_VALLE2_ts))
)

# Crear gráfico con ggplot2
 stl_df_var2 <- ggplot(stl_df_var2, aes(x = Time, y = Value, color = Component)) +
  geom_line() +
  facet_wrap(~Component, scales = "free_y", ncol = 1) + 
  theme_minimal() +
  labs(title = "Descomposición temporal de la variable VEH_VALLE",
       x = "Tiempo",
       y = "Valor")

# Convertir a gráfico interactivo con plotly
ggplotly(stl_df_var2)

Extracción señales variable 3

# Cargar librerías necesarias
library(ggplot2)
library(plotly)

# Descomposición de la serie temporal
stl_decomp_VEH_CALI3_ts <- stl(VEH_CALI3_ts, s.window = "periodic")

# Convertir la descomposición a un data frame para graficar con ggplot2
stl_df_var3 <- data.frame(
  Time = rep(time(VEH_CALI3_ts), 12),  # Tiempo repetido para cada componente (son 12 componentes)
  Value = c(stl_decomp_VEH_CALI3_ts$time.series[, "seasonal"], 
            stl_decomp_VEH_CALI3_ts$time.series[, "trend"], 
            stl_decomp_VEH_CALI3_ts$time.series[, "remainder"], 
            VEH_CALI3_ts),
  Component = rep(c("Estacional", "Tendencia", "Residuo", "Serie Original"), each = length(VEH_CALI3_ts))
)

# Crear gráfico con ggplot2
 stl_df_var3 <- ggplot(stl_df_var3, aes(x = Time, y = Value, color = Component)) +
  geom_line() +
  facet_wrap(~Component, scales = "free_y", ncol = 1) + 
  theme_minimal() +
  labs(title = "Descomposición temporal de la variable VEH_CALI",
       x = "Tiempo",
       y = "Valor")

# Convertir a gráfico interactivo con plotly
ggplotly(stl_df_var3)

Después de la descomposición temporal de cada variable, se extrae la variable ajustada por estacionalidad para graficarla junto con la serie original:

Se crea la variable1 VEH_COLOMBIA ajustada por estacionalidad

# Extraer los componentes de la descomposición
VEH_COLOMBIA_sa <- VEH_COLOMBIA1_ts - stl_decomp_VEH_COLOMBIA1_ts$time.series[, "seasonal"]

Se crea la variable2 VEH_VALLE ajustada por estacionalidad

# Extraer los componentes de la descomposición
VEH_VALLE_sa <- VEH_VALLE2_ts - stl_decomp_VEH_VALLE2_ts$time.series[, "seasonal"]

Se crea la variable3 VEH_CALI ajustada por estacionalidad

# Extraer los componentes de la descomposición
VEH_CALI_sa <- VEH_CALI3_ts - stl_decomp_VEH_CALI3_ts$time.series[, "seasonal"]

Ahora si se puede graficar las series originales versus la ajustada por estacionalidad

Gráfico serie original VS ajustada Variable 1 VEH_COLOMBIA

# Crear vector de fechas correctamente alineado con la serie
fechas_VEH_COLOMBIA <- seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = length(VEH_COLOMBIA1_ts))

# Gráfico mejorado con fechas en el eje X
grafico_ajustada_VEH_COLOMBIA <- ggplot() +
  geom_line(aes(x = fechas_VEH_COLOMBIA, y = VEH_COLOMBIA1_ts), color = "grey", size = 0.5, linetype = "solid", name = "Serie Original") +
  geom_line(aes(x = fechas_VEH_COLOMBIA, y = VEH_COLOMBIA_sa), color = "black", size = 0.6, linetype = "solid", name = "Serie Ajustada") +
  ggtitle("Variable 1 VEH_COLOMBIA:Serie Original vs Serie Ajustada por Estacionalidad") +
  xlab("Tiempo") +
  ylab("Unidad de medida variable 1 VH_COLOMBIA") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas para mejor visualización

# Convertir a gráfico interactivo
ggplotly(grafico_ajustada_VEH_COLOMBIA)

Gráfico serie original VS ajustada Variable 2 VEH_VALLE

# Crear vector de fechas correctamente alineado con la serie
fechas_VEH_VALLE <- seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = length(VEH_VALLE2_ts))

# Gráfico mejorado con fechas en el eje X
grafico_ajustada_VEH_VALLE2 <- ggplot() +
  geom_line(aes(x = fechas_VEH_VALLE, y = VEH_VALLE2_ts), color = "grey", size = 0.5, linetype = "solid", name = "Serie Original") +
  geom_line(aes(x = fechas_VEH_VALLE, y = VEH_VALLE_sa), color = "black", size = 0.6, linetype = "solid", name = "Serie Ajustada") +
  ggtitle("Variable 2 VEH_VALLE:Serie Original vs Serie Ajustada por Estacionalidad") +
  xlab("Tiempo") +
  ylab("Unidad de medida variable 2 VEH_VALLE") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas para mejor visualización

# Convertir a gráfico interactivo
ggplotly(grafico_ajustada_VEH_VALLE2)

Gráfico serie original VS ajustada Variable 3 VEH_CALI

# Crear vector de fechas correctamente alineado con la serie
fechas_VEH_CALI <- seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = length(VEH_CALI3_ts))

# Gráfico mejorado con fechas en el eje X
grafico_ajustada_VEH_CALI3 <- ggplot() +
  geom_line(aes(x = fechas_VEH_CALI, y = VEH_CALI3_ts), color = "grey", size = 0.5, linetype = "solid", name = "Serie Original") +
  geom_line(aes(x = fechas_VEH_CALI, y = VEH_CALI_sa), color = "black", size = 0.6, linetype = "solid", name = "Serie Ajustada") +
  ggtitle("Variable 3 VEH_CALI:Serie Original vs Serie Ajustada por Estacionalidad") +
  xlab("Tiempo") +
  ylab("Unidad de medida variable 3") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas para mejor visualización

# Convertir a gráfico interactivo
ggplotly(grafico_ajustada_VEH_CALI3)

Ahora graficamos serie original vs tendencia

Primero se debe obtener la tendencia de cada variable y luego graficarla

Tendencia Variable 1 VEH_COLOMBIA

library(ggplot2)
library(plotly)

# Convertir la serie a un vector numérico
VEH_COLOMBIA_vec <- as.numeric(VEH_COLOMBIA1_ts)
tendencia_VEH_COLOMBIA <- as.numeric(stl_decomp_VEH_COLOMBIA1_ts$time.series[, "trend"])

# Asegurar que 'fechas' tenga la misma longitud
fechas <- seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = length(VEH_COLOMBIA1_ts))

# Gráfico interactivo de la serie original vs tendencia
grafico_tendencia_VEH_COLOMBIA <- ggplot() +
  geom_line(aes(x = fechas, y = VEH_COLOMBIA_vec, color = "Serie Original"), size = 0.7, linetype = "solid") +
  geom_line(aes(x = fechas, y = tendencia_VEH_COLOMBIA, color = "Tendencia"), size = 0.8, linetype = "solid") +
  scale_color_manual(values = c("Serie Original" = "grey", "Tendencia" = "black")) +
  ggtitle("Variable 1 VEH_COLOMBIA: Serie Original vs Tendencia") +
  xlab("Tiempo") +
  ylab("Unidad de medida Variable 1 VEH_COLOMBIA") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas del eje X

# Convertir a gráfico interactivo con plotly
ggplotly(grafico_tendencia_VEH_COLOMBIA)

Tendencia Variable 2 VEH_VALLE

library(ggplot2)
library(plotly)

# Convertir la serie a un vector numérico
VEH_VALLE_vec <- as.numeric(VEH_VALLE2_ts)
tendencia_VEH_VALLE <- as.numeric(stl_decomp_VEH_VALLE2_ts$time.series[, "trend"])

# Asegurar que 'fechas' tenga la misma longitud
fechas <- seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = length(VEH_VALLE2_ts))

# Gráfico interactivo de la serie original vs tendencia
grafico_tendencia_VEH_VALLE <- ggplot() +
  geom_line(aes(x = fechas, y = VEH_VALLE_vec, color = "Serie Original"), size = 0.7, linetype = "solid") +
  geom_line(aes(x = fechas, y = tendencia_VEH_VALLE, color = "Tendencia"), size = 0.8, linetype = "solid") +
  scale_color_manual(values = c("Serie Original" = "grey", "Tendencia" = "black")) +
  ggtitle("Variable 2 VEH_VALLE: Serie Original vs Tendencia") +
  xlab("Tiempo") +
  ylab("Unidad de medida Variable 2 VEH_VALLE") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas del eje X

# Convertir a gráfico interactivo con plotly
ggplotly(grafico_tendencia_VEH_VALLE)

Tendencia Variable 3 VEH_CALI

library(ggplot2)
library(plotly)

# Convertir la serie a un vector numérico
VEH_CALI_vec <- as.numeric(VEH_CALI3_ts)
tendencia_VEH_CALI <- as.numeric(stl_decomp_VEH_CALI3_ts$time.series[, "trend"])

# Asegurar que 'fechas' tenga la misma longitud
fechas <- seq.Date(from = as.Date("2014-01-01"), by = "month", length.out = length(VEH_CALI3_ts))

# Gráfico interactivo de la serie original vs tendencia
grafico_tendencia_VEH_CALI <- ggplot() +
  geom_line(aes(x = fechas, y = VEH_CALI_vec, color = "Serie Original"), size = 0.7, linetype = "solid") +
  geom_line(aes(x = fechas, y = tendencia_VEH_CALI, color = "Tendencia"), size = 0.8, linetype = "solid") +
  scale_color_manual(values = c("Serie Original" = "grey", "Tendencia" = "black")) +
  ggtitle("Variable 3 VEH_CALI: Serie Original vs Tendencia") +
  xlab("Tiempo") +
  ylab("Unidad de medida Variable 3 VEH_CALI") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas del eje X

# Convertir a gráfico interactivo con plotly
ggplotly(grafico_tendencia_VEH_CALI)

Ahora calculamos la tasa de crecimiento de la serie original vs tendencia:

Tasa de crecimiento de la serie de tendencia y original para la variable 1 VEH_COLOMBIA

#Cálculo de la tasa de crecimiento anual correctamente alineada
tasa_crecimiento_VEH_COLOMBIA <- (VEH_COLOMBIA1_ts[(13:length(VEH_COLOMBIA1_ts))] / VEH_COLOMBIA1_ts[1:(length(VEH_COLOMBIA1_ts) - 12)] - 1) * 100
tasa_tendencia_VEH_COLOMBIA <- (tendencia_VEH_COLOMBIA[(13:length(tendencia_VEH_COLOMBIA))] / tendencia_VEH_COLOMBIA[1:(length(tendencia_VEH_COLOMBIA) - 12)] - 1) * 100

# Crear vector de fechas corregido, es decir que inicie desde enero 2015
fechas_corregidas_VEH_COLOMBIA <- seq(from = as.Date("2015-01-01"), by = "month", length.out = length(tasa_crecimiento_VEH_COLOMBIA))

# Verificar longitudes
print(length(fechas_corregidas_VEH_COLOMBIA))

## [1] 120

print(length(tasa_crecimiento_VEH_COLOMBIA))

## [1] 120

print(length(tasa_tendencia_VEH_COLOMBIA))

## [1] 120

*Gráfico variable original y tendencia variable 1 VEH_COLOMBIA: tasa de crecimiento anual**

# Gráfico de la tasa de crecimiento anual variable 1 VEH_COLOMBIA
grafico_crecimiento_VEH_COLOMBIA <- ggplot() +
  geom_line(aes(x = fechas_corregidas_VEH_COLOMBIA, y = tasa_crecimiento_VEH_COLOMBIA), color = "grey", size = 0.7) +
  geom_line(aes(x = fechas_corregidas_VEH_COLOMBIA, y = tasa_tendencia_VEH_COLOMBIA), color = "black", size = 0.8, linetype = "dashed") +
  ggtitle("Variable1 VEH_COLOMBIA: Tasa de crecimiento anual %") +
  xlab("Tiempo") +
  ylab("% de Crecimiento Anual") +
  scale_y_continuous(labels = scales::label_comma()) +
  coord_cartesian(ylim = c(-150, 150)) +
  theme_minimal()
# Convertir a gráfico interactivo
ggplotly(grafico_crecimiento_VEH_COLOMBIA)

Ahora calculamos la tasa de crecimiento de la serie original vs tendencia: variable 2 VEH_VALLE

#Cálculo de la tasa de crecimiento anual correctamente alineada
tasa_crecimiento_VEH_VALLE <- (VEH_VALLE2_ts[(13:length(VEH_VALLE2_ts))] / VEH_VALLE2_ts[1:(length(VEH_VALLE2_ts) - 12)] - 1) * 100
tasa_tendencia_VEH_VALLE <- (tendencia_VEH_VALLE[(13:length(tendencia_VEH_VALLE))] / tendencia_VEH_VALLE[1:(length(tendencia_VEH_VALLE) - 12)] - 1) * 100

# Crear vector de fechas corregido
fechas_corregidas_VEH_VALLE <- seq(from = as.Date("2015-01-01"), by = "month", length.out = length(tasa_crecimiento_VEH_VALLE))

# Verificar longitudes
print(length(fechas_corregidas_VEH_VALLE))

## [1] 120

print(length(tasa_crecimiento_VEH_VALLE))

## [1] 120

print(length(tasa_tendencia_VEH_VALLE))

## [1] 120

# Gráfico de la tasa de crecimiento anual variable 2 VEH_VALLE
grafico_crecimiento_VEH_VALLE <- ggplot() +
  geom_line(aes(x = fechas_corregidas_VEH_VALLE, y = tasa_crecimiento_VEH_VALLE), color = "grey", size = 0.7) +
  geom_line(aes(x = fechas_corregidas_VEH_VALLE, y = tasa_tendencia_VEH_VALLE), color = "black", size = 0.8, linetype = "dashed") +
  ggtitle("Variable2 VEH_VALLE: Tasa de crecimiento anual %") +
  xlab("Tiempo") +
  ylab("% de Crecimiento Anual") +
  scale_y_continuous(labels = scales::label_comma()) +
  coord_cartesian(ylim = c(-150, 150)) +
  theme_minimal()
# Convertir a gráfico interactivo
ggplotly(grafico_crecimiento_VEH_VALLE)

Ahora calculamos la tasa de crecimiento de la serie original vs tendencia: variable 3 VEH_CALI

#Cálculo de la tasa de crecimiento anual correctamente alineada
tasa_crecimiento_VEH_CALI <- (VEH_CALI3_ts[(13:length(VEH_CALI3_ts))] / VEH_CALI3_ts[1:(length(VEH_CALI3_ts) - 12)] - 1) * 100
tasa_tendencia_VEH_CALI <- (tendencia_VEH_CALI[(13:length(tendencia_VEH_CALI))] / tendencia_VEH_CALI[1:(length(tendencia_VEH_CALI) - 12)] - 1) * 100

# Crear vector de fechas corregido
fechas_corregidas_VEH_CALI <- seq(from = as.Date("2015-01-01"), by = "month", length.out = length(tasa_crecimiento_VEH_CALI))

# Verificar longitudes
print(length(fechas_corregidas_VEH_CALI))

## [1] 120

print(length(tasa_crecimiento_VEH_CALI))

## [1] 120

print(length(tasa_tendencia_VEH_CALI))

## [1] 120

# Gráfico de la tasa de crecimiento anual variable 3 VEH_CALI
grafico_crecimiento_VEH_CALI <- ggplot() +
  geom_line(aes(x = fechas_corregidas_VEH_CALI, y = tasa_crecimiento_VEH_CALI), color = "grey", size = 0.7) +
  geom_line(aes(x = fechas_corregidas_VEH_CALI, y = tasa_tendencia_VEH_CALI), color = "black", size = 0.8, linetype = "dashed") +
  ggtitle("Variable3 VEH_CALI: Tasa de crecimiento anual %") +
  xlab("Tiempo") +
  ylab("% de Crecimiento Anual") +
  scale_y_continuous(labels = scales::label_comma()) +
  coord_cartesian(ylim = c(-150, 150)) +
  theme_minimal()
# Convertir a gráfico interactivo
ggplotly(grafico_crecimiento_VEH_CALI)

Analizar la tasa de crecimiento anual ayuda a detectar cambios en el entorno económico que afectan el sector. Se pueden prever crisis o períodos de auge y prepararse para ellos.

Modelo ARIMA

División en conjunto de entrenamiento y prueba para la variable 1 VEH_COLOMBIA que es la elegida para pronosticar

El código siguiente divide una serie temporal (variable1_ts) en dos subconjuntos:

Conjunto de entrenamiento (train): Datos desde enero de 2014 hasta septiembre de 2024. Conjunto de prueba (test): Datos desde octubre de 2024 hasta diciembre de 2024.

Esto se hace para evaluar el desempeño de modelos de predicción en datos no vistos.

# Esta división idealmente podria se 80%-70% de los datos para entrenamiento y 20%-30% para prueba o test

# En este ejemplo el conjunto de entrenamiento es: Enero 2014-Septiembre 2024 y  el conjunto de prueba o test: noviembre 2024-diciembre 2024 

train_size <- length(VEH_COLOMBIA1_ts) - 3 # Se deja fuera los últimos 3 valores para usarlos como set de prueba.
train_ts <- window(VEH_COLOMBIA1_ts, end = c(2024, 9))  # Entrenamiento hasta septiembre 2024
test_ts <- window(VEH_COLOMBIA1_ts, start = c(2024, 10))  # Prueba inicia desde oct2024

Modelo ARIMA automático normal (sin tener en cuenta el factor estacional)

Identificación automática del modelo ARIMA

library(forecast)

# Ajustar un modelo ARIMA automático sin estacionalidad, por eso se pone seasonal=FALSE
auto_arima_model_no_seasonal <- auto.arima(train_ts, seasonal = FALSE)

# Mostrar el modelo seleccionado
summary(auto_arima_model_no_seasonal)

## Series: train_ts 
## ARIMA(1,1,2) 
## 
## Coefficients:
##          ar1      ma1     ma2
##       0.7711  -1.3734  0.3966
## s.e.  0.1198   0.1463  0.1323
## 
## sigma^2 = 17741966:  log likelihood = -1248.95
## AIC=2505.91   AICc=2506.23   BIC=2517.32
## 
## Training set error measures:
##                     ME     RMSE      MAE       MPE     MAPE      MASE
## Training set -434.9487 4146.303 2840.497 -67.76767 77.01546 0.7694171
##                     ACF1
## Training set -0.03775687

Estimación del modelo identificado automatico y validación de Significancia de coeficientes

library(lmtest)

# Evaluar la significancia estadística de los coeficientes del modelo ARIMA
coeftest(auto_arima_model_no_seasonal)

## 
## z test of coefficients:
## 
##     Estimate Std. Error z value  Pr(>|z|)    
## ar1  0.77111    0.11976  6.4389 1.204e-10 ***
## ma1 -1.37338    0.14632 -9.3860 < 2.2e-16 ***
## ma2  0.39664    0.13232  2.9977  0.002721 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# Ajuste del modelo ARIMA(1,1,2) automático sin parte estacional y crearlo como variable darima_auto para luego poder graficarlo y crear la tabla
darima_auto <- Arima(train_ts, 
                order = c(1, 1, 2))  # Especificamos directamente (p=1, d=1, q=2)  

# Mostrar resumen del modelo ajustado
summary(darima_auto)

## Series: train_ts 
## ARIMA(1,1,2) 
## 
## Coefficients:
##          ar1      ma1     ma2
##       0.7711  -1.3734  0.3966
## s.e.  0.1198   0.1463  0.1323
## 
## sigma^2 = 17741966:  log likelihood = -1248.95
## AIC=2505.91   AICc=2506.23   BIC=2517.32
## 
## Training set error measures:
##                     ME     RMSE      MAE       MPE     MAPE      MASE
## Training set -434.9487 4146.303 2840.497 -67.76767 77.01546 0.7694171
##                     ACF1
## Training set -0.03775687

Validación de residuales o errores del modelo

# Diagnóstico del modelo (los residuos deben ser ruido blanco)
checkresiduals(darima_auto)  # Verificar si los residuos son aleatorios y no presentan patrones

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,1,2)
## Q* = 91.804, df = 21, p-value = 7.87e-11
## 
## Model df: 3.   Total lags used: 24

Pronóstico modelo ARIMA automático dentro de muestra o en el set de prueba

# Generar pronóstico para el conjunto de prueba
forecast_arima_auto <- forecast(darima_auto, h = length(test_ts))  # Predecir los valores futuros

# Crear dataframe para gráfico interactivo del pronóstico
forecast_data_auto <- data.frame(Tiempo = time(forecast_arima_auto$mean), 
                            Pronostico = as.numeric(forecast_arima_auto$mean),
                            Observado = as.numeric(test_ts))

# Graficar pronóstico junto con los valores observados reales
p1auto <- ggplot(forecast_data_auto, aes(x = Tiempo)) +
  geom_line(aes(y = Pronostico, color = "Pronóstico")) +
  geom_line(aes(y = Observado, color = "Observado")) +
  ggtitle("Pronóstico vs Observado") +
  xlab("Tiempo") + ylab("variable1 VEH_COLOMBIA")

ggplotly(p1auto)  # Convertir el gráfico en interactivo

Interpretación modelo automatico (1,1,2):El modelo automático (1,1,2) no es tan acertado, el modelo no ajusta bien a la serie observada en el horizonte del pronostico, la diferencia entre lo observado y el pronostico sugiere un sesgo sistematico, lo que podria ser resultado de una estacionalidad mal capturada, para tal caso el modelo arima no modela directamente estacionalidades.

Pronóstico automático dentro del set de prueba como tabla

# Cargar librerías necesarias
library(forecast)
library(dplyr)

# Generar pronóstico con el modelo ARIMA identificado
arima_forecast_auto <- forecast(auto_arima_model_no_seasonal, h = length(test_ts))

# Crear un dataframe con los valores observados y pronosticados
forecast_table_auto <- data.frame(
  Tiempo = time(arima_forecast_auto$mean),  # Extraer las fechas del pronóstico
  Observado = as.numeric(test_ts),  # Valores reales
  Pronosticado = as.numeric(arima_forecast_auto$mean)  # Valores pronosticados
)

# Mostrar la tabla
print(forecast_table_auto)

##     Tiempo Observado Pronosticado
## 1 2024.750     18521     17098.79
## 2 2024.833     21824     17230.50
## 3 2024.917     25331     17332.06

Ahora pronosticamos con el modelo automatico fuera del periodo de análisis, es decir enero 2025

Es decir, le sumamos al periodo de prueb auna observación más. Es decir, se estan pronosticando 4 observaciones o meses.

# Cargar librerías necesarias
library(forecast)

# Hacer un pronóstico para el siguiente mes (1 período adicional)
next_forecast_auto <- forecast(auto_arima_model_no_seasonal, h = length(test_ts) + 1)

# Extraer el pronóstico del próximo mes
next_month_forecast_auto <- data.frame(
  Tiempo = time(next_forecast_auto$mean),  # Extraer la fecha del pronóstico
  Pronostico = as.numeric(next_forecast_auto$mean)  # Valor pronosticado
)

# Mostrar el pronóstico completo
print(next_month_forecast_auto)

##     Tiempo Pronostico
## 1 2024.750   17098.79
## 2 2024.833   17230.50
## 3 2024.917   17332.06
## 4 2025.000   17410.38

# Extraer solo el valor del mes adicional (último de la tabla)
next_month <- tail(next_month_forecast_auto, 1)
print(paste("Pronóstico para enero 2025:", next_month$Tiempo, "=", next_month$Pronostico))

## [1] "Pronóstico para enero 2025: 2025 = 17410.3806092015"

Modelo SARIMA automático

Este modelo podria ser una solución o mejora al modelo arima tradicional ya que recoge el efecto estacional de las variables, es recomendable por tanto para datos que si tienen un componente estacional fuerte.

El modelo ajustado en este ejemplo es un SARIMA(0,1,1)(1,0,0)[12], lo que significa:

(0,1,1): Parte ARIMA no estacional: 0 términos autorregresivos (AR). 1 diferenciación (d), lo que indica que la serie fue diferenciada una vez para hacerla estacionaria. 1 término de media móvil (MA).

(1,0,0)[12]: Parte estacional con periodicidad 12 (mensual si los datos son mensuales): 1 término autorregresivo estacional (SAR). 0 diferenciaciones estacionales. 0 términos de media móvil estacionales (SMA).

El modelo SARIMA(0,1,1)(1,0,0)[12] sugiere que el comportamiento de la serie temporal depende de:

Una tendencia no estacionaria que fue estabilizada con una sola diferenciación.

Un componente de ruido (shock) reciente que se corrige con el término MA.

Un efecto estacional anual (por ejemplo, ventas similares cada diciembre), modelado por el término autorregresivo estacional.

Este tipo de modelo es comúnmente utilizado en variables que presentan patrones estacionales regulares y comportamiento persistente en el tiempo.

Identificación dautomática del modelo SARIMA

# Identificación automática modelo SARIMA
auto_arima_model <- auto.arima(train_ts)  # Busca automáticamente los mejores parámetros del modelo ARIMA
print(auto_arima_model)

## Series: train_ts 
## ARIMA(1,0,1)(0,1,1)[12] with drift 
## 
## Coefficients:
##          ar1      ma1     sma1     drift
##       0.7888  -0.2348  -0.7359  -70.6592
## s.e.  0.0862   0.1431   0.1037   29.5113
## 
## sigma^2 = 8678623:  log likelihood = -1103.54
## AIC=2217.07   AICc=2217.61   BIC=2230.88

A continuación, se crea el objeto darima para luegO poder graficar los valores reales y observados:

# Cargar el paquete necesario
library(forecast)

# Ajustar el modelo SARIMA(1,0,1)(0,1,1)[12] #Modelo identificado en el paso anterior
darima <- Arima(train_ts, 
                order = c(1, 0, 1),  # (p,d,q) -> (0,1,1)
                seasonal = list(order = c(0, 1, 1),  # (P,D,Q) -> (0,1,1)
                                period = 12))  # Periodicidad estacional de 12 meses

# Mostrar resumen del modelo ajustado
summary(darima)

## Series: train_ts 
## ARIMA(1,0,1)(0,1,1)[12] 
## 
## Coefficients:
##          ar1      ma1     sma1
##       0.8594  -0.2746  -0.6908
## s.e.  0.0685   0.1319   0.0972
## 
## sigma^2 = 8996361:  log likelihood = -1105.43
## AIC=2218.86   AICc=2219.21   BIC=2229.9
## 
## Training set error measures:
##                     ME     RMSE      MAE       MPE     MAPE      MASE
## Training set -460.2031 2819.623 1872.857 -56.37066 63.58204 0.5073084
##                      ACF1
## Training set -0.006174814

Validación de residuales del modelo automatico SARIMA

# Diagnóstico del modelo (los residuos deben ser ruido blanco)
checkresiduals(darima)  # Verificar si los residuos son aleatorios y no presentan patrones

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,0,1)(0,1,1)[12]
## Q* = 18.835, df = 21, p-value = 0.5958
## 
## Model df: 3.   Total lags used: 24

Pronóstico con el modelo SARIMA dentro del set de prueba-Gráfico líneas

# Generar pronóstico para el conjunto de prueba
forecast_arima <- forecast(darima, h = length(test_ts))  # Predecir los valores futuros

# Crear dataframe para gráfico interactivo del pronóstico
forecast_data <- data.frame(Tiempo = time(forecast_arima$mean), 
                            Pronostico = as.numeric(forecast_arima$mean),
                            Observado = as.numeric(test_ts))

# Graficar pronóstico junto con los valores observados reales
p4 <- ggplot(forecast_data, aes(x = Tiempo)) +
  geom_line(aes(y = Pronostico, color = "Pronóstico")) +
  geom_line(aes(y = Observado, color = "Observado")) +
  ggtitle("Pronóstico vs Observado") +
  xlab("Tiempo") + ylab("Unidad Variable 1")

ggplotly(p4)  # Convertir el gráfico en interactivo

Interpretación modelo SARIMA (1,0,1) (0,1,1): El gráfico compara los valores observados (línea roja) con los pronosticados (línea azul) para una variable temporal durante el último trimestre de 2024. La tendencia de ambos conjuntos es creciente, lo que indica que el modelo captó correctamente la dirección general del comportamiento de la variable.

Sin embargo, se observa una subestimación sistemática del modelo: los valores pronosticados son consistentemente más bajos que los valores reales.

Esto sugiere que, aunque el modelo es estructuralmente correcto (detecta la tendencia), podría haber un sesgo negativo o estar omitiendo otros factores influyentes que explica el mayor crecimiento observado.

Pronóstico del modelo automático SARIMA en el set de prueba-Tabla

# Cargar librerías necesarias
library(forecast)
library(dplyr)

# Generar pronóstico con el modelo ARIMA identificado
arima_forecast <- forecast(auto_arima_model, h = length(test_ts))

# Crear un dataframe con los valores observados y pronosticados
forecast_table <- data.frame(
  Tiempo = time(arima_forecast$mean),  # Extraer las fechas del pronóstico
  Observado = as.numeric(test_ts),  # Valores reales
  Pronosticado = as.numeric(arima_forecast$mean)  # Valores pronosticados
)

# Mostrar la tabla
print(forecast_table)

##     Tiempo Observado Pronosticado
## 1 2024.750     18521     16948.05
## 2 2024.833     21824     18731.03
## 3 2024.917     25331     21993.52

Pronóstico del modelo automático SARIMA fuera de muestra, es decir, en enero 2025

Es decir, le sumamos al periodo de prueba una observación más. Es decir, se estan pronosticando 4 meses.

# Cargar librerías necesarias
library(forecast)

# Hacer un pronóstico para el siguiente mes (1 período adicional)
next_forecast <- forecast(auto_arima_model, h = length(test_ts) + 1)

# Extraer el pronóstico del próximo mes
next_month_forecast <- data.frame(
  Tiempo = time(next_forecast$mean),  # Extraer la fecha del pronóstico
  Pronostico = as.numeric(next_forecast$mean)  # Valor pronosticado
)

# Mostrar el pronóstico completo
print(next_month_forecast)

##     Tiempo Pronostico
## 1 2024.750   16948.05
## 2 2024.833   18731.03
## 3 2024.917   21993.52
## 4 2025.000   11743.77

# Extraer solo el valor del trimestre adicional (último de la tabla)
next_month <- tail(next_month_forecast, 1)
print(paste("Pronóstico para enero 2025:", next_month$Tiempo, "=", next_month$Pronostico))

## [1] "Pronóstico para enero 2025: 2025 = 11743.76551736"

pronostico modelo SARIMA

Según el modelo, el valor estimado para enero de 2025 es 11,743.77, lo que representa una caída abrupta en comparación con diciembre de 2024 (21,993.52). Esta reducción de casi el 50% sugiere un comportamiento atípico, sin embargo en el analisis realizado en la base de datos de la variable, se denota que enero es un mes con una demanda descendiente posiblemente debido a factores estacionales, desgaste economico del mercado, fin del ciclo comercial de modelos, menor oferta promocional, habitos de consumo, entre otros.

Conclusión final del modelo estaditico escogido (SARIMA):

El modelo SARIMA(1,0,1)(0,1,1)[12] resulta ser el más apropiado para modelar y pronosticar las ventas de vehículos en colombia, ya que:

Captura explícitamente la estacionalidad mensual, fundamental en una variable como las ventas de vehículos, que típicamente presentan patrones cíclicos a lo largo del año (por promociones, cierres fiscales, variaciones de demanda, etc.).

Presenta residuos más aleatorios (tipo ruido blanco), con menor autocorrelación residual, lo que indica que el modelo recoge mejor la estructura temporal de la serie.

Aunque el ARIMA(1,1,2) tiene residuos con distribución más simétrica, no modela estacionalidad, lo que lo limita para series mensuales como esta.

El SARIMA también permite una mejor descomposición e interpretación de los efectos separados: tendencia, ruido y estacionalidad, lo cual es valioso para análisis económico y toma de decisiones públicas o empresariales.

unque el modelo ARIMA(1,1,2) muestra residuos más normalmente distribuidos, el modelo SARIMA es más apropiado pues la serie presenta estacionalidad anual, ya que:

El SARIMA reduce mejor la autocorrelación pues considera la estructura estacional, clave porque la serie tiene patrones anuales.

logo2

Recomendaciones para el concesionario de vehículos

Verificación y ajuste del modelo de pronóstico siempre y cuando se agreguen otros factores influyentes en el mercado automotriz. Revisar si hay cambios estructurales no capturados por el modelo (como cierre de año contable, políticas tributarias, escasez de inventario, politicas monetarias, etc.), considerar integrar variables influentes como promociones, tasas de interés, normativas fiscales o estacionalidad del mercado automotor.
Prepararse para una posible baja en ventas Si el modelo refleja una proyección realista, se debera plantear una reducción preventiva de inventario para evitar sobrecostos por acumulación de stock. Se debera realizar el ajuste de metas comerciales para enero de 2025, preparando al equipo de ventas con objetivos realistas, se deberan optimizar los gastos operativos, priorizando inversiones necesarias y aplazando aquellas que no generen valor inmediato.
Crear estrategias de ventas anticipadas; por lo cual se deberan ejecutar campañas de pre-venta y descuentos en diciembre para mitigar el bajón de enero, se debera incentivar la venta de unidades del 2024 con bonos de mantenimiento o financiación diferenciada.
Fortalecer canales digitales; debido a que enero suele ser un mes de bajo tráfico físico, por lo que potenciar las paginas web puede compensar la caída en ventas, al realizar promociones en redes sociales, cotizaciones virtuales, test drives a domicilio y agendamiento online pueden marcar diferencia.
Se recomienda realizar un análisis cualitativo del entorno, a manera de ejemplo, consultar factores del entorno económico relacionada a reformas tributarias, politicas monetarias, entre otros.

Conclusión general

La drástica caída proyectada para enero de 2025 debe ser tratada con cautela, debido a que puede ser un reflejo de una estacionalidad típica, sin embargo, para un concesionario de vehículos, anticiparse estratégicamente con campañas, ajustes de inventario y acciones comerciales inteligentes puede transformar un mes crítico en una oportunidad.

Caso 2. Extracción de señales y Pronóstico ARIMA_Vehiculos

Miguel Patiño Y Sebastian Soto

2025-06-07