UNIVERSIDAD NACIONAL DE INGENIERÍA

FACULTAD DE INGENIERÍA ECONÓMICA,

ESTADÍSTICA Y CCSS
ESCUELA PROFESIONAL DE INGENIERÍA
ESTADÍSTICA

SEMINARIO CAPSTONE

TÍTULO
“Comparación de métodos Box-Jenkins y Holt-Winters para el
pronóstico de la producción de plata en el departamento de
Áncash-Perú de febrero a diciembre 2024”

POR:
Ccora Yupanqui, Steve
De la Rosa Hinostroza, Bridgyth
López Meléndez, Edwin

SUPERVISOR:
Leighton L. Estrada Rayme

LIMA-PERÚ
2024
 Abstract

In the present study, a comparative analysis was carried out between the Box-Jenkis
and Holt-Winters methods to forecast silver production in the Ancash region, Peru for
the period from February to December 2024. The predictions generated by each method
and its results were compared.
Firstly, it was found that the SARIMA model(0, 1, 1)(1, 0, 1)12 presented the lowest
forecast error among all the variants considered within the ARIMA/SARIMA family.
This model provided relevant information about the behavior pattern of the time se-
ries, which allowed improving forecast efficiency.
Second, the Holt-Winters model demonstrated good predictive capabilities when the
level, trend, and seasonality smoothing parameters were set to specific values. These
parameters (0.4345, 0 and 0, respectively) allowed us to obtain satisfactory predictions
with this approach.
In summary, the study contrasted the Box-Jenkins and Holt-Winters methods for fo-
recasting a specific time series. the Holt-Winters model provides more accurate and
reliable forecasts for this specific time series. Therefore, it is recommended to use the
Holt-Winters model to make predictions of silver production in the Ancash region, Pe-
ru for the period from February to December 2024.
vspace10pt

Keywords: Box-Jenkins, Holt-Winters, Forecasting

1
 Resumen

En el presente estudio, se realizó un análisis comparativo entre los métodos Box-

Jenkis y Holt-Winters para pronosticar la producción de plata en la region de Ancash,
Perú para el periodo de febrero a diciembre del 2024. Se evaluaron las predicciones
generadas por cada método y se contrastaron sus resultados.
En primer lugar, se encontró que el modelo SARIMA(0, 1, 1)(1, 0, 1)12 presentó el me-
nor error de pronóstico entre todas las variantes consideradas dentro de la familia
ARIMA/SARIMA. Este modelo proporcionó información relevante sobre el patrón de
comportamiento de la serie de tiempo, lo que permitió mejorar la eficiencia del pronósti-
co.
En segundo lugar, el modelo de Holt-Winters demostró buenas capacidades predictivas
cuando se establecieron los parámetros de suavización del nivel, tendencia y estaciona-
lidad en valores especı́ficos. Estos parámetros (0.4345, 0 y 0, respectivamente) permi-
tieron obtener predicciones satisfactorias con este enfoque.
En resumen, el estudio contrastó los métodos Box-Jenkins y Holt-Winters para el
pronóstico de una serie de tiempo especı́fica. el modelo Holt-Winters proporciona
pronósticos más precisos y confiables para esta serie temporal especı́fica. Por lo tanto,
se recomienda utilizar el modelo Holt-Winters para realizar predicciones de la produc-
ción de plata en la region de Ancash, Perú para el periodo de febrero a diciembre del
2024.

Palabras clave: Box-Jenkins, Holt-Winters, Pronóstico

2
Índice

1. INTRODUCCIÓN 5
1.1. Descripción de la situación problemática . . . . . . . . . . . . . . . . . 7
1.2. Formulación del problema . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.2.1. Problema general . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.2.2. Problema especı́fico . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.3. Objetivos de la investigación . . . . . . . . . . . . . . . . . . . . . . . . 7
1.3.1. Objetivo general . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.3.2. Objetivos especı́ficos . . . . . . . . . . . . . . . . . . . . . . . . 7
1.4. Justificación, alcances y limitaciones . . . . . . . . . . . . . . . . . . . 8
1.4.1. Justificación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.4.2. Alcances . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.4.3. Limitaciones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

2. FUNDAMENTOS TEÓRICOS 10
2.1. Antecedentes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.1.1. Antecedentes nacionales . . . . . . . . . . . . . . . . . . . . . . 10
2.1.2. Antecedentes internacionales . . . . . . . . . . . . . . . . . . . . 11
2.2. Bases teóricas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2.2.1. Conceptos Generales . . . . . . . . . . . . . . . . . . . . . . . . 12
2.2.2. Prueba de Dickey-Fuller . . . . . . . . . . . . . . . . . . . . . . 15
2.2.3. Prueba de Ljung-Box . . . . . . . . . . . . . . . . . . . . . . . . 16
2.2.4. Modelos lineales estacionarios . . . . . . . . . . . . . . . . . . . 17
2.2.5. Modelos lineales no Estacionarios . . . . . . . . . . . . . . . . . 17
2.3. Metodologı́a Box-Jenkins . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.3.1. Identificación del modelo . . . . . . . . . . . . . . . . . . . . . . 19
2.3.2. Estimación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.3.3. Validación o verificación . . . . . . . . . . . . . . . . . . . . . . 20
2.3.4. Pronóstico o predicción . . . . . . . . . . . . . . . . . . . . . . . 20
2.4. Método Holt-Winters . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.4.1. Holt-Winters aditivo . . . . . . . . . . . . . . . . . . . . . . . . 21
2.4.2. Holt-Winters multiplicativo . . . . . . . . . . . . . . . . . . . . 21
2.4.3. Estimación de α, β y γ e indicadores de error del modelo . . . . 22
2.5. Métricas de comparación . . . . . . . . . . . . . . . . . . . . . . . . . . 22

3. MATERIALES, MÉTODOS Y PROCEDIMIENTOS 23

3.1. Tipo, nivel y diseño de la investigación . . . . . . . . . . . . . . . . . . 23
3.1.1. Tipo de investigación . . . . . . . . . . . . . . . . . . . . . . . . 23
3.1.2. Nivel de investigación . . . . . . . . . . . . . . . . . . . . . . . . 23
3.1.3. Diseño de la investigación . . . . . . . . . . . . . . . . . . . . . 23
3.1.4. Población y muestra . . . . . . . . . . . . . . . . . . . . . . . . 23
3.1.5. Técnicas de recolección de datos . . . . . . . . . . . . . . . . . . 23
3.1.6. Descripción de la variable . . . . . . . . . . . . . . . . . . . . . 23
3.1.7. Herramientas y software . . . . . . . . . . . . . . . . . . . . . . 23
3.1.8. Metodologı́a Box-Jenkins . . . . . . . . . . . . . . . . . . . . . . 24
3.1.9. Método Holt-Winters . . . . . . . . . . . . . . . . . . . . . . . . 25

3
4. PRESENTACIÓN Y ANÁLISIS DE RESULTADOS 27
4.1. Análisis exploratorio . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
4.1.1. Suavizado de la serie . . . . . . . . . . . . . . . . . . . . . . . . 28
4.1.2. Gráfico de cajas . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
4.1.3. Gráficos de autocorrelación (ACF) y autocorrelación parcial (PACF) 30
4.1.4. Prueba de Estacionariedad . . . . . . . . . . . . . . . . . . . . . 32
4.1.5. Diferenciación . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
4.2. Modelamiento Box-Jenkins . . . . . . . . . . . . . . . . . . . . . . . . . 35
4.3. Modelamiento Holt-Winters . . . . . . . . . . . . . . . . . . . . . . . . 40
4.4. Comparación de modelos . . . . . . . . . . . . . . . . . . . . . . . . . . 44

5. CONCLUSIONES 45

6. BIBLIOGRAFIA 46

7. ANEXOS 47
7.1. NORMAS, ESTANDARES Y REGULACIONES . . . . . . . . . . . . 47
7.2. RESTRICCIONES MULTIPLES, LIMITACIONES DE LA INVESTI-
GACION . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
7.3. CÓDIGO DEL PROGRAMA . . . . . . . . . . . . . . . . . . . . . . . 49
7.4. HOJAS DE DATOS . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54

4
1. INTRODUCCIÓN
La minerı́a es una de las principales actividades económicas en Perú, destacándose
particularmente en la producción de plata.
Según Minerı́a & Energı́a (2023), Perú, reconocido globalmente como un gigante en la
producción de plata, enfrentó un panorama complejo en 2023. A pesar de un incremento
mensual en la producción, el paı́s experimentó una disminución anual, lo que planteó
interrogantes sobre las tendencias del mercado y la sustentabilidad de la industria.
Según datos proporcionados el Ministerio de Energı́a y Minas de Perú (MINEM) que
se pueden observar en la Fig. 1, la producción de plata en diciembre de 2023 fue de
296,800 kg.f, un 6.5 % más que los 278,608 kg de diciembre de 2022. Este aumento
puntual contrasta con la tendencia anual. En el periodo de enero a diciembre de 2023,
la producción totalizó 3,038,485 kg.f, marcando una caı́da en comparación con el mismo
periodo de 2022, que registró 3,059,412 kg.f.
Esta disminución se atribuye principalmente a la menor producción de dos gigantes
mineros: Compañı́a Minera Antamina y Compañı́a Minera Ares, con reducciones del
27.6 % y 8.8 % respectivamente.

Figura 1: Producción de plata total (kg. f) de enero 2022 a octubre 2023.

Fuente: Elaboración propia.

En términos regionales, Pasco lidera con un 20.1 % del total nacional, seguido por
Áncash con el 17.8 % y Junı́n con el 16.0 %. Estas cifras reflejan no solo la distribución
geográfica de la minerı́a de plata en Perú, sino también la importancia económica y
social de la minerı́a para estas regiones.
La producción de plata en Perú está influenciada por diversos factores, tanto in-
ternos como externos. Entre ellos, las polı́ticas gubernamentales, la inversión en tecno-
logı́a, los precios globales de la plata, y las condiciones geológicas y climáticas juegan
roles cruciales. Además, el contexto social y los desafı́os ambientales presentan retos
constantes para la industria.
El sector minero en Perú enfrenta desafı́os significativos, incluyendo la necesidad de
prácticas mineras sostenibles y responsables socialmente. La inversión en tecnologı́as
más limpias y eficientes, ası́ como en el desarrollo de relaciones comunitarias sólidas,
son fundamentales para el futuro de la industria.

5
 Por otro lado, las oportunidades de crecimiento y desarrollo son evidentes. La de-
manda mundial de plata, impulsada por su uso en diversas industrias como la electróni-
ca y la medicina, ofrece un panorama prometedor para el sector.
La producción de plata en Perú en 2023 refleja un sector en constante evolución,
enfrentando desafı́os pero también abriendo puertas a nuevas oportunidades. La clave
para el éxito sostenido de la industria radica en la adaptación a los cambios del mercado,
la inversión en tecnologı́a y prácticas sostenibles, y el compromiso con las comunidades
locales.
Especı́ficamente, la producción de plata en el departamento de Áncash, ha des-
empeñado un papel importante en la economı́a regional y nacional durante décadas.
Áncash se destaca como una de las principales regiones productoras de metales precio-
sos en el paı́s, contribuyendo significativamente al crecimiento económico y al desarrollo
industrial. En los últimos años, la producción de plata en esta región ha experimen-
tado fluctuaciones que pueden ser atribuidas a una variedad de factores, incluyendo
cambios en los precios internacionales de los metales, variaciones en la demanda global
y dinámicas locales de extracción minera.
Según datos proporcionados por el Ministerio de Energı́a y Minas de Perú (MI-
NEM), que se pueden observar en la Fig. 2, la producción de plata en el departamento
de Áncash durante diciembre de 2023 reportó un total de 51,660.7 kg.f, representando
un aumento en comparación al mismo mes del año 2022 con 45,715.8 kg.f.

Figura 2: Producción de plata Áncash, Perú (kg.f) de enero 2022 a diciembre 2023.

Fuente: Elaboración propia.

Esta tendencia alcista en la producción refleja el continuo crecimiento y la impor-

tancia de la industria minera en la región. Sin embargo, la predicción precisa de la
producción futura de plata en Áncash es fundamental para la planificación estratégica
y la toma de decisiones en el sector minero.
El presente trabajo de investigación se centra en evaluar y comparar la eficacia de
los métodos Box-Jenkins y Holt-Winters en la predicción de la producción de plata
en el departamento de Áncash para los meses restantes del 2024, con el objetivo de
proporcionar herramientas analı́ticas que faciliten la toma de decisiones en la industria
minera y contribuyan al crecimiento económico y social de la región.

6
1.1. Descripción de la situación problemática
La producción de plata en Áncash enfrenta desafı́os significativos relacionados con
la necesidad de anticipar y gestionar las fluctuaciones en la extracción de este metal
precioso. La variabilidad en los precios internacionales de los metales, los cambios en
la demanda global, la disponibilidad de recursos minerales y los desafı́os ambientales y
sociales asociados con la actividad minera son solo algunos de los factores que influyen
en la producción de plata en la región. Estas fluctuaciones dificultan la planificación
estratégica y la toma de decisiones tanto para las empresas mineras como para las
autoridades gubernamentales, lo que subraya la importancia de desarrollar modelos
de pronósticos precisos y confiables. Por lo tanto, surge la necesidad de investigar
y comparar diferentes enfoques de pronóstico, como los métodos Box-Jenkis y Holt-
Winters, para mejorar la capacidad predictiva y proporcionar información valiosa para
la gestión de la producción de plata en Áncash.

1.2. Formulación del problema

1.2.1. Problema general
¿Cuál de los modelos de los métodos Box-Jenkis o Holt-Winters ofrece el mejor
rendimiento en el pronóstico de la producción de plata en el departamento de Áncash,
Perú?

1.2.2. Problema especı́fico

¿Cuál es el modelo del método Box-Jenkis que mejor se ajusta a la serie de tiempo
de producción de plata en Áncash, Perú?
¿Cuál es el modelo del método Holt-Winters que mejor se ajusta a la serie de
tiempo de producción de plata en Áncash, Perú?

1.3. Objetivos de la investigación

1.3.1. Objetivo general
Determinar el mejor modelo entre los métodos Box-Jenkins y Holt-Winters para
la predicción de la producción de plata en el departamento de Áncash, Perú, para el
periodo de febrero a diciembre de 2024.

1.3.2. Objetivos especı́ficos

Implementar y ajustar un modelo del método Box-Jenkins para el pronóstico de
la producción de plata en el departamento de Áncash, Perú, para el periodo de
febrero a diciembre de 2024.
Implementar y ajustar un modelo del método Holt-Winters para el pronóstico de
la producción de plata en el departamento de Áncash, Perú, para el periodo de
febrero a diciembre de 2024.

7
1.4. Justificación, alcances y limitaciones
1.4.1. Justificación

La producción de plata en Áncash, Perú, ha sido una actividad económica importan-

te durante siglos. Sin embargo, la industria ha experimentado fluctuaciones significati-
vas en la producción de los últimos años. Esto ha llevado a la necesidad de desarrollar
modelos precisos para pronosticar la producción futura de plata, lo que puede ayu-
dar tanto las empresas mineras como al gobierno regional de Áncash a tomar mejores
decisiones de inversión y planificación respecto a la producción de este mineral.
Los métodos Box-Jenkins y Holt-Winters son dos métodos estadı́sticos comunes
utilizados para pronosticar series de tiempo. Los modelos del método Box-Jenkins
son adecuados para series de tiempo que son estacionarias, mientras que los modelos
Holt-Winters son más flexibles y pueden usarse para series de tiempo con tendencias
estacionales [Chatfield, 2015].
Ambos modelos han demostrado ser precisos para pronosticar series de tiempo de
naturaleza similar. Por ejemplo, un estudio aplicó con éxito modelos del método Box-
Jenkins, especı́ficamente modelos SARIMA para pronosticar la producción de cobre en
una mina chilena [Carmona et. al., 2018]. Otro estudio encontró que los modelos Holt-
Winters eran efectivos para pronosticar la demanda de plata a corto plazo [Mohammadi
et. al., 2020].

1.4.2. Alcances
Este trabajo se enfocará en la comparación de los modelos de los métodos Box-
Jenkins y Holt-Winters usando datos proporcionados por el BCRP de la serie Produc-
ción de productos mineros según departamentos - Plata - Áncash (enero de 2001 a
enero de 2024).

Los modelos de los métodos Box-Jenkins y Holt-Winters pueden generar pronósti-

cos precisos de la producción de plata, lo que puede ayudar a las empresas mineras
y al gobierno a tomar mejores decisiones de inversión y planificación respecto a
la producción de este mineral.
Las empresas mineras pueden usar pronósticos de producción de plata del mejor
modelo para mejorar la eficiencia de sus operaciones, lo que puede conducir a
menores costos y mayores ganancias.
Los pronósticos precisos de la producción de plata pueden alentar a los inversores
a invertir en la industria minera de Áncash, lo que puede crear empleos y estimular
el crecimiento económico.

8
1.4.3. Limitaciones
La precisión de los modelos de los métodos Box-Jenkins y Holt-Winters depende
de la calidad y disponibilidad de los datos históricos de producción de plata. Si los
datos son incompletos o inexactos, esto puede afectar la precisión de los modelos. A la
fecha en la que se recopilaron los datos de este trabajo, solo se dispone de los datos de
producción hasta enero de 2024.
La producción de plata puede verse afectada por una variedad de factores externos,
como los precios de los metales, las condiciones geológicas y las regulaciones guberna-
mentales. Los modelos de los métodos Box-Jenkins y Holt-Winters no pueden tener
en cuenta todos estos factores, ya que se trata de modelos univariados, por lo que los
pronósticos generados por estos modelos pueden estar sujetos a cierta incertidumbre.

9
2. FUNDAMENTOS TEÓRICOS
2.1. Antecedentes
2.1.1. Antecedentes nacionales
1. Mejı́a, E., Gonzales, S. (2019) en el estudio “Predicción del consumo de energı́a
eléctrica residencial de la Región de Cajamarca mediante modelos Holt-Winters”
realizaron un estudio con el objetivo de predecir el consumo de energı́a eléctrica
residencial de la Región Cajamarca de las tres empresas de esta región; Electro-
norte S.A., Hidrandina S.A., y Electro Oriente S.A., de enero 2005 a diciembre
2017.
La predicción se hizo usando modelos Holt-Winters ya que se consideran útiles
para la predicción a corto y mediano plazo de ventas de energı́a eléctrica. La
previsión del consumo de energı́a eléctrica tiene importancia en la planificación
energética regional y nacional; a partir de sus resultados los agentes del mercado
de energı́a eléctrica toman decisiones más adecuadas para su labor. El méto-
do Holt-Winters se estimó para diferentes constantes de suavización, que incluye
métodos para patrones estacionales aditivos y multiplicativos, considerando como
principal estimador de la capacidad del modelo el error porcentual medio (MA-
PE); el modelo de Holt-Winters aditivo fue seleccionado como el mejor modelo de
predicción, presentando menor error de ajuste entre el modelo y la serie de datos,
un MAPE de 1.7 y constantes de suavización estimadas de 0.4.
2. Chafloque, M., Chung, V., Cotrina, A. (2019) vieron la importancia de conocer
el comportamiento de las exportaciones mineras desde 1985 hasta el 2017, por lo
que la investigación “Modelación de las dinámicas, estimación y análisis de las
exportaciones mineras: caso peruano” tuvo como objetivo aplicar metodologı́as de
predicción y suavizamiento para pronosticar a 12 meses las exportaciones mineras
(en millones de dólares) con datos mensuales del Banco Central de Reserva del
Perú (BCRP) sobre la exportación minera desde 1985 a 2017. Se utilizó el análisis
estadı́stico descriptivo y la comparación de los métodos predictivos ARIMA y
SARIMA. Los resultados muestran que la mejor proyección a corto plazo se realiza
con el modelo SARIMA (0, 1, 1)(0, 1, 1)12 . Los resultados de esta investigación
pueden ser utilizados por instituciones tanto públicas como privadas para tomar
decisiones relacionadas con inversiones en el sector minero. Se concluyó que el
modelo SARIMA es óptimo para este tipo de predicciones.
3. Chung, A., e Inche, J. (2022) muestran en la “Aplicación de herramientas eco-
nométricas y machine learning para series de tiempo univariadas en una escuela
profesional de una universidad pública” un modelo de pronóstico para series de
tiempo que permite determinar la cantidad de posibles alumnos matriculados por
ciclo. Se tomó como muestra a los alumnos matriculados en el curso de Quı́mica
Orgánica de la Universidad Nacional Mayor de San Marcos (UNMSM) entre los
periodos académicos de 2009-II al 2021-I. Para ello, se probaron tres metodo-
logı́as: redes neuronales aplicadas a series de tiempos, ARIMA y Holt-Winters.
La primera de las tres metodologı́as mencionadas pertenece al campo del machine
learning y las dos restantes son ampliamente utilizadas en econometrı́a. Se con-
cluyó que el modelo que mejor se adapta a la serie de datos es Holt-Winters, ya
que posee un mejor ajuste. La principal contribución radica en presentar nuevos

10
 campos de aplicación de la econometrı́a y machine learning, esta vez orientados
a la gestión de las escuelas profesionales en las universidades públicas y cómo
se puede aprovechar el potencial de estas herramientas predictivas en diversos
campos.

2.1.2. Antecedentes internacionales

1. Marino, M., Arango, A., Lotero, L., y Jiménez, M. (2021) utilizaron metodologı́as
de análisis de series de tiempo para ajustar, validar y comparar ”Modelos de series
temporales para pronóstico de la demanda eléctrica del sector de explotación de
minas y canteras en Colombia”. Se tomaron 102 registros mensuales disponibles en
el portal XM, desde enero de 2010 hasta junio de 2018. Se ajustó un modelo aditivo
de componentes y errores autorregresivos de orden dos, un modelo SARIMA y un
modelo Holt-Winters. Los tres modelos cumplieron los supuestos sobre errores, y
el modelo Holt-Winters fue seleccionado como el más óptimo debido a su menor
nivel de error.
2. Garza, D., Ovalle, H., Tovar, O. (2022) realizaron un estudio para pronosticar la
producción de minerales y el prospecto productivo para desarrollo exploratorio,
ası́ como las regalı́as en los departamentos extractivos del paı́s a través del tiem-
po, con énfasis en la extracción de carbón. Se utilizó una base de datos histórica
proporcionada por la Agencia Nacional de Minerı́a (ANM), superpuesta con da-
tos de la Unidad de Planeación Minero-Energética (UPME), con una muestra
de tamaño (n=45215). Se aplicaron herramientas estadı́sticas descriptivo-analı́ti-
cas como SARIMA y SARIMAX para modelar y pronosticar adecuadamente la
exploración, producción y contraprestación del carbón en los departamentos de
Cesar y La Guajira.
3. Nina Yucra, S. M. (2021) desarrolló un modelo predictivo de producción de mi-
nerales de plata en Bolivia basado en datos históricos. Utilizó redes neuronales
artificiales para optimizar la producción de este mineral. El estudio aplicó in-
teligencia artificial para el pronóstico de la producción de minerales de plata,
ofreciendo resultados que pueden informar estrategias y polı́ticas para la planifi-
cación de la producción de este mineral, utilizando tanto herramientas de machine
learning como econometrı́a para el pronóstico.

11
2.2. Bases teóricas
2.2.1. Conceptos Generales
1. Serie de tiempo. Una serie de tiempo es un conjunto de observaciones generadas
secuencialmente a lo largo del tiempo.
Si el conjunto es continuo, se dice que la serie temporal es continua. Si el conjunto
es discreto, se dice que la serie temporal es discreta. Por lo tanto, las observaciones
de una serie de tiempo discreta realizadas en los momentos τ1 , τ2 , . . . , τt , . . . , τN
pueden denotarse por Z(τ1 ), Z(τ2 ), . . . , Z(τt ), . . . , Z(τN ). En las series de tiem-
po discretas, las observaciones se realizan en un intervalo fijo h. Cuando tene-
mos N valores sucesivos de dicha serie disponibles para el análisis, escribimos
Z1 , Z2 , . . . , Zt , Zt+1 , . . . , ZN para denotar observaciones realizadas en intervalos
de tiempo equidistantes τ0 + h, τ0 + 2h, . . . , τ0 + th, τ0 + (t + 1)h, . . . , τ0 + N h.
Es necesario definir exactamente los tiempos de observación, lo cual puede espe-
cificarse mediante estos dos valores. Si adoptamos τ0 como el origen y h como
la unidad de tiempo, podemos considerar Zt como la observación en el tiempo t
(Box, Jenkins, Reinsel, & Ljung, 2015).
2. Estacionariedad. Las caracterı́sticas de un proceso estocástico pueden ser Gaus-
siano (variables aleatorias que tienen una distribución normal) o Estacionario
(está asociada con la estabilidad en el tiempo). Se conocen dos tipos de estacio-
nariedad:
Estacionariedad en sentido estricto: Una serie de tiempo se considera estric-
tamente estacionaria si la distribución conjunta de Zt1 , Zt2 , . . . , Zts es la misma
que la de Zt1 +k , Zt2 +k , . . . , Zts +k , para todo t1 , t2 , . . . , ts y k.
Estacionariedad en sentido débil: Se considera débilmente estacionaria si
cumple con las siguientes condiciones:
La media de la serie de tiempo es constante a lo largo del tiempo, es decir,
se cumple que µt = µ para todo t.
La varianza de la serie de tiempo es constante a lo largo del tiempo, es decir,
se cumple que σt2 = σ 2 para todo t.
La covarianza entre dos observaciones de la serie de tiempo depende solo de
la diferencia en sus tiempos, es decir, Cov(Zt , Zt+k ) = γ(k), donde γ(k) es
una función que solo depende de k y no de t.
3. Descomposición en Componentes no observables
La serie de tiempo se puede descomponer en cuatro componentes:
a) Tendencia (T ): Movimientos o variaciones de larga duración que se man-
tienen durante todo el periodo de observación.
b) Ciclo (C): Oscilaciones periódicas a largo plazo alrededor de la tendencia
producidas por perı́odos irregulares y alternativos de prosperidad y depre-
sión.
c) Estacionalidad (S): Movimientos o variaciones periódicas regulares que se
producen dentro de un periodo anual por motivos estrictamente estacionales.
d ) Aleatoriedad (A): Movimientos erráticos generados por causas ajenas a la
variación sistemática, es aleatoria.

12
4. Esquemas de descomposición en componentes no observables
La descomposición de una serie tiene el siguiente modelo estadı́stico: Zt = f (Tt , St , ϵt ),
donde la función f puede ser lineal o no lineal. Las estructuras más usuales son:
a) Aditivo: Zt = Tt + Ct + St + At
b) Multiplicativo: Zt = Tt · Ct · St · At
c) Mixto: Zt = Tt · Ct · St + At
5. Procesos Estocásticos. Un fenómeno estadı́stico que evoluciona en el tiempo
según leyes probabilı́sticas se denomina proceso estocástico. La serie de tiempo
que se va a analizar puede considerarse entonces como una realización particular,
producida por el mecanismo de probabilidad subyacente, del sistema en estudio.
En otras palabras, las series temporales se definen como un caso particular de los
procesos estocásticos.
Procesos estocásticos estacionarios. Un proceso estocástico es estrictamente
estacionario si sus propiedades no se ven afectadas por un cambio de origen tem-
poral, es decir, si la distribución de probabilidad conjunta asociada con m observa-
ciones Zt1 , Zt2 , . . . , Ztm , realizadas en cualquier conjunto de tiempos t1 , t2 , . . . , tm ,
es la misma que la asociada con m observaciones Zt1 +k , Zt2 +k , . . . , Ztm +k , hechas
en los momentos t1 + k, t2 + k, . . . , tm + k. Por lo tanto, para que un proceso
discreto sea estrictamente estacionario, la distribución conjunta de cualquier con-
junto de observaciones no debe verse afectada al desplazar todos los tiempos de
observación hacia adelante o hacia atrás en cualquier cantidad entera k.
Media y varianza de un proceso estacionario.
Función media: Cuando m = 1, el supuesto de estacionariedad implica que
la distribución de probabilidad p(zt ) es la misma para cualquier t y puede
escribirse como p(z). Por lo tanto, el proceso estocástico tiene una media
constante. Z ∞
µ = E[Zt ] = zp(z) dz (1)
−∞

que define el nivel sobre el cual fluctúa.

Función Varianza: Una varianza constante que mide su dispersión alrededor
de este nivel.
Z ∞
2 2
σZ = E[(Zt − µ) ] = (Z − µ)2 p(z) dz (2)
−∞

13
6. Funciones de autocovarianza y autocorrelación
Autocovarianza. El supuesto de estacionariedad también implica que la dis-
tribución de probabilidad conjunta p(zt1 , zt2 ) es la misma para todos los tiempos
t1 , t2 , que están separados por un intervalo constante. En particular, se deduce que
la covarianza entre los valores Zt y Zt+k , separados por k intervalos de tiempo o
por el rezago k, debe ser la misma para todo t bajo el supuesto de estacionariedad.
Esta covarianza se llama autocovarianza en el retraso k y se define como:

γk = cov[Zt , Zt+k ] = E[(Zt − µ)(Zt+k − µ)] (3)

Autocorrelación. Según Villavicencio (2010), en ocasiones en una serie de tiem-

po los valores que toma una variable en el tiempo no son independientes entre
sı́, sino que un valor determinado depende de los valores anteriores. Existen dos
formas de medir esta dependencia de las variables:
Función de autocorrelación (ACF): La autocorrelación mide la correlación
entre dos variables separadas por k periodos.
cov(Zj , Zj−k )
ρj = corr(Zj , Zj−k ) = p p (4)
Var(Zj ) Var(Zj−k )

La función de autocorrelación simple tiene las siguientes propiedades:

ρ0 = 1
−1 ≤ ρj ≤ 1
ρj = ρ−j
Función de autocorrelación parcial (PACF): La autocorrelación parcial mi-
de la correlación entre dos variables separadas por k periodos cuando no se con-
sidera la dependencia creada por los retardos intermedios existentes entre ambas.

πj = corr (Zj , Zj−k /Zj−1 , Zj−2 , . . . , Zj−k+1 ) (5)

cov(Zj − Ẑj , Zj−k − Ẑj−k )

πj = q q (6)
V (Zj − Ẑj ) V (Zj−k − Ẑj−k )

7. Ruido Blanco
El ruido blanco es un caso simple de los procesos estocásticos que son indepen-
dientes e idénticamente distribuidas, se denota por εt . Formalmente, una serie de
tiempo se considera ruido blanco si cumple con las condiciones:

E(εt ) = 0 (7)

V ar(εt ) = σ 2 (8)
Cov(εt , εt+j ) = 0, ∀j ̸= 0 (9)

14
 8. Medias Móviles
Para Parra y Vicente una media móvil trata, sencillamente de una media aritméti-
ca que se caracteriza porque toma un valor para cada momento del tiempo y
porque en su cálculo no entran todas las observaciones de la muestra disponible.
Entre los distintos tipos de medias móviles que se pueden construir nos vamos
a referir a dos tipos: medias móviles centradas y medias móviles asimétricas. El
primer tipo se utiliza para la representación de la tendencia, mientras que el se-
gundo lo aplicaremos para la predicción en modelos con media constante.
Entre los distintos tipos de medias móviles que se pueden construir nos vamos
a referir a dos tipos: medias móviles centradas y medias móviles asimétricas. El
primer tipo se utiliza para la representación de la tendencia, mientras que el se-
gundo lo aplicaremos para la predicción en modelos con media constante.

Media móvil centrada. Se caracterizan porque el número de observaciones

que entran en su cálculo es impar, asignándose cada media móvil a la observación
central. Ası́, una media móvil centrada en t de longitud 2n + 1 viene dada por la
siguiente expresión:
n
1 X
M (Zt ) = Zt+j (10)
2n + 1 j=−n

El subı́ndice asignado a la media móvil, t, es el mismo que el de la observación

central, Yt. Obsérvese también que, por construcción, no se pueden calcular las
medias móviles correspondientes a las n primeras y a las n últimas observaciones.

La utilización de medias móviles implica la elección arbitraria de su longitud

u orden, es decir, del número de observaciones que intervienen en el cálculo de
cada media móvil. Cuanto mayor sea la longitud, mejor se eliminarán las irre-
gularidades de la serie, ya que al intervenir más observaciones en su cálculo se
compensarán las fluctuaciones de este tipo, pero por el contrario, el coste informa-
tivo será mayor. Por el contrario, cuando la longitud es pequeña, la media móvil
refleja con mayor rapidez los cambios que puedan producirse en la evolución de
la serie. Es conveniente, pues, sopesar estos factores al decidir la longitud de la
media móvil.

2.2.2. Prueba de Dickey-Fuller

David Dickey y Wayne Fuller (1979) desarrollaron un procedimiento estadı́stico uti-
lizado para determinar la presencia de raı́ces unitarias en una serie temporal. La prue-
ba evalúa la hipótesis nula de que una serie temporal tiene una raı́z unitaria contra la
hipótesis alternativa de que es estacionaria después de eliminar dicha raı́z unitaria. Una
serie temporal estacionaria es aquella cuyas propiedades estadı́sticas, como la media y
la varianza, son constantes a lo largo del tiempo.
La forma más común de la prueba de Dickey-Fuller es la versión aumentada (ADF),
propuesta por Said E. Said y David A. Dickey en 1984, que incluye términos adicionales
en el modelo de regresión para mejorar la robustez y la precisión de la prueba.
La prueba se realiza a través de un proceso de regresión donde la serie temporal Zt se

15
regresa contra su rezago Zt−1 , con la inclusión de rezagos adicionales y términos deter-
minı́sticos si es necesario. El estadı́stico de prueba resultante se compara con valores
crı́ticos tabulados para determinar si se puede rechazar la hipótesis nula de raı́z unitaria.

Si consideramos un proceso estocástico de la forma:

Zt = ϕZt−1 + εt , −1 ≤ ϕ ≤ 1 (11)
donde yt es la serie de tiempo y εt es un término de error de ruido blanco.

H0 : La serie temporal tiene una raı́z unitaria, no es estacionaria.

H1 : La serie temporal es estacionaria, no tiene una raı́z unitaria.

2.2.3. Prueba de Ljung-Box

La prueba de Ljung-Box es una técnica estadı́stica utilizada para examinar la pre-
sencia de autocorrelación en una serie temporal. Desarrollada por Greta M. Ljung y
George E. P. Box en 1978, esta prueba evalúa si los primeros m autocorrelaciones de
una serie temporal son conjuntamente diferentes de cero. La prueba es útil para iden-
tificar patrones en los datos que no son aleatorios, lo que puede indicar la presencia de
estructura en la serie temporal, como tendencias o ciclos.

H0 : La serie temporal no presenta autocorrelación para los retardos 1, 2, . . . , m.

H1 : Al menos una de las autocorrelaciones de los retardos 1, 2, . . . , m es diferente
de cero.

El estadı́stico de prueba de Ljung-Box, Q, se calcula de la siguiente manera:

m
X ρ̂2k
Q = n(n + 2) (12)
k=1
n−k
donde n es el tamaño de la muestra y ρ̂k es la autocorrelación muestral en el re-
tardo k. El estadı́stico Q sigue aproximadamente una distribución χ2 con m grados de
libertad.

16
2.2.4. Modelos lineales estacionarios
1. Modelo Autorregresivo AR(p)
Un modelo denotado como AR(p) es una serie temporal que se basa en la idea de
que el valor actual Zt se puede explicar por los valores anteriores hasta p periodos
atrás. El caso más simple es un AR(1) o ARIMA(1, 0, 0).

Zt = ϕ0 + ϕ1 Zt−1 + ϕ2 Zt−2 + . . . + ϕp Zt−p + εt (13)

Donde:
Zt es el valor de la serie temporal en el tiempo t.
ϕ0 es el término constante de la serie.
ϕ1 , ϕ2 , . . . , ϕp son los parámetros autoregresivos del modelo.
εt es el término de error o ruido blanco en el tiempo t, con media cero y
varianza constante σ 2 .
2. Modelo de Medias Móviles MA(q)
El modelo estacionario de medias móviles de orden q, denotado como MA(q),
representa el valor de la serie temporal Zt como la suma ponderada de los errores
de periodos precedentes hasta q.

Zt = θ0 − θ1 εt−1 − θ2 εt−2 − . . . − θq εt−q + εt (14)

Donde:
Zt es el valor de la serie temporal en el tiempo t.
θ0 es el término constante de la serie.
θ1 , θ2 , . . . , θq son los parámetros de la media móvil del modelo.
εt es el término de error o ruido blanco en el tiempo t, con media cero y
varianza constante σ 2 .
3. Modelo mixto ARMA(p,q)
Es una serie temporal denotada por su abreviatura ARMA(p, q), donde el valor de
la variable está representado en función de sus observaciones pasadas y los valores
de los errores actuales y rezagados. En este proceso, se utilizan “p” términos
autorregresivos y “q” términos de media móvil.

Zt = ϕ0 + ϕ1 Zt−1 + ϕ2 Zt−2 + . . . + ϕp Zt−p + θ1 εt−1 + θ2 εt−2 + . . . + θq εt−q + εt (15)

2.2.5. Modelos lineales no Estacionarios

1. Modelo Autorregresivo Integrado de Media Móvil ARIMA (p, d. q)
Son considerados una extensión del modelo ARMA, que se ocupan de series con
problemas de estacionariedad. En otros términos, se refieren a series donde la me-
dia y la varianza no se mantienen constantes a lo largo del tiempo. Estos modelos
pertenecen a una clase particular de procesos, donde la cantidad de términos au-
torregresivos está denotada por “p”, mientras que “d” indica la cantidad de veces
que se debe diferenciar la serie para lograr la estacionariedad, y “q” la cantidad

17
 de términos de medias móviles. Es decir, se basan primordialmente en convertir
una serie no estacionaria en una serie estacionaria mediante el uso de un modelo
ARIMA(p, d, q).
Su expresión algebraica es:

∆d Φ(B)Zt = Θ(B)εt (16)

(1 − ϕ1 B − ϕ2 B 2 − · · · − ϕp B p )Zt = (1 − θ1 B − θ2 B 2 − · · · − θq B q )εt (17)
Donde:
Φ(B) = 1 − ϕ1 B − ϕ2 B 2 − · · · − ϕp B p es el polinomio autoregresivo de orden
p.
Θ(B) = 1 − θ1 B − θ2 B 2 − · · · − θq B q es el polinomio de media móvil de orden
q.
∆d es el operador de diferencia y d es el número de diferencias necesarias
para retirar la tendencia de la serie.
B es el operador de rezago.
2. Modelo Autorregresivo Integrado de Media Móvil con Parte Estacional
SARIMA (p, d, q) (P, D, Q)
Estos procesos son caracterizados por sus fluctuaciones cı́clicas, conocidas también
como variaciones estacionales, estas series pueden ser representadas mediante un
modelo SARIMA (p, d, q) (P, D, Q), donde el primer parte minúscula indica la
parte regular y la mayúscula indica la parte estacional.

Φp (B)ΦP (B s )∆d ∆D s
s Zt = Θq (B)ΘQ (B )εt (18)
Donde:
Φp (B) = 1 − ϕ1 B − ϕ2 B 2 − · · · − ϕp B p es el polinomio autoregresivo no
estacional de orden p.
Θq (B) = 1 − θ1 B − θ2 B 2 − · · · − θq B q es el polinomio de media móvil no
estacional de orden q.
ΦP (B s ) = 1 − Φ1 B s − Φ2 B 2s − · · · − ΦP B P s es el polinomio autoregresivo
estacional de orden P.
ΘQ (B s ) = 1 − Θ1 B s − Θ2 B 2s − · · · − ΘQ B Qs es el polinomio de media móvil
estacional de orden Q.
∆ds es el operador de diferencia generalizado, cuando dos observaciones están
distantes entre sı́ de s intervalos de tiempo que presentan alguna semejan-
za y D es el numero de diferencias de rezagos s necesarias para retirar la
estacionalidad de la serie.
εt es el término de error o ruido blanco en el tiempo t.

18
2.3. Metodologı́a Box-Jenkins
La metodologı́a de los modelos ARIMA fue formalizada por Box y Jenkins en 1976,
por lo que también se les denomina modelos Box- Jenkins. Este enfoque parte del
hecho de que la serie temporal que se trata de predecir es geneerada por un proceso
estocástico cuya naturaleza puede ser caracterizada mediante un modelo. Para efec-
tuar la estimación de un modelo ARIMA se requiere de una serie temporal mensual
o trimestral que cuente con un elevado número de observaciones. Básicamente, la me-
todologı́a Box-Jenkins consiste en encontrar un modelo matemático que represente el
comportamiento de una serie temporal de datos, y permita hacer previsiones única-
mente introduciendo el perı́odo de tiempo correspondiente (Chatfield, 1989).
Su metodologı́a consta de cuatro etapas: identificación del modelo, estimación, verifi-
cación o validación y pronóstico.

2.3.1. Identificación del modelo

En esta etapa se identifican los parámetros más apropiados para los posbiles mode-
los. Es crucial que la serie sea estacionaria para proceder. Las herramientas principales
utilizadas en la identificación son:
Gráfico de la serie temporal: Proporciona una visualización estructurada del
comportamiento de la serie a lo largo del tiempo. Este gráfico indica si la serie es
estacionaria. Una serie se considera estacionaria si su oscilación alrededor de algún
valor no se aleja permanentemente.
Correlogramas: Se emplean los correlogramas de la autocorrelación (ACF) y au-
tocorrelación parcial (PACF) para identificar los posibles órdenes p y q del modelo
ARIMA. Los picos en estos gráficos ayudan a determinar los valores adecuados para
los componentes autoregresivos y de media móvil.
Durante la identificación del modelo:

Si la serie tiene tendencia, se aplican diferencias regulares hasta que desaparezca.

Normalmente, el orden de diferencia es 1, raramente mayor a 3.
Si la serie presenta estacionalidad, se aplican diferencias estacionales hasta que
desaparezca el patrón estacional. Usualmente, no es necesario aplicar más de una
diferencia estacional.
Si la serie no tiene varianza constante (es heterocedástica), se debe transformar la
serie. A menudo, basta con aplicar el logaritmo, aunque existen transformaciones
más avanzadas como las de Box-Cox.

Una vez que la transformación de la serie indica que es estacionaria, se procede a

deducir la estructura del modelo observando los correlogramas ACF y PACF de la serie
transformada.

2.3.2. Estimación
Después de identificar los valores más apropiados para (p, d y q) de los modelos
seleccionados, el siguiente paso es estimar los parámetros autorregresivos (AR) y de
media móvil (MA) incluidos en el modelo.
El proceso de estimación generalmente implica utilizar algoritmos de optimización para
encontrar los valores de los parámetros que maximizan la función de verosimilitud.

19
Este proceso puede ser computacionalmente intensivo y puede requerir la exploración
de múltiples combinaciones de parámetros.

2.3.3. Validación o verificación

Después de seleccionar el mejor modelo según los criterios previamente mencionados
y estimar los parámetros, se procede a verificar si dicho modelo seleccionado se ajusta
de manera razonable a los datos. Es posible que existan otros modelos que también se
ajusten.
Un modelo será válido cuando los residuos estimados se comporten como ruido blanco,
lo cual implica que deben seguir una distribución normal con media nula, covarianza
nula y varianza constante.
En los gráficos de las funciones de autocorrelación (AFC) y autocorrelación parcial
(FACP), los valores de los residuos deben ser no significativos, es decir, pequeños y
dentro de las bandas o lı́neas de confianza.

2.3.4. Pronóstico o predicción

Después de identificar y validar el modelo más adecuado, se puede utilizar para
realizar predicciones futuras.

2.4. Método Holt-Winters

El método Holt-Winters, también conocido como suavización exponencial triple, es
una técnica ampliamente utilizada para el pronóstico de series temporales. Este método
es particularmente útil cuando los datos muestran tanto tendencia como estacionalidad.
Fue desarrollado por Charles Holt en 1957 y refinado por Peter Winters en 1960.
Existen dos variantes principales del método Holt-Winters: aditiva y multiplicativa,
dependiendo de la naturaleza de la estacionalidad (constante o proporcional). Este
método tiene los siguientes componentes:
1. Nivel (at ): Representa el valor promedio de la serie en el tiempo t. Se estima
utilizando la suavización exponencial del último nivel, ajustado por la tendencia
y la estacionalidad.
2. Tendencia (bt ): Captura la tasa de cambio del nivel en el tiempo t. Se estima
utilizando la suavización exponencial de la última tendencia, ajustada por la
tendencia anterior y la estacionalidad.
3. Estacionalidad (st ): Representa los efectos estacionales que se repiten en interva-
los regulares.
Dependiendo de si la estacionalidad es aditiva o multiplicativa, se estima utilizan-
do la suavización exponencial de la última estacionalidad, ajustada por el nivel y
la tendencia.
Estos componentes se actualizan y suavizan iterativamente utilizando ecuaciones es-
pecı́ficas que dependen de los parámetros de suavización (α, β y γ), los cuales controlan
la velocidad de adaptación de cada componente a los datos observados. La elección en-
tre el modelo aditivo y el multiplicativo depende de cómo varı́a la estacionalidad con
respecto al nivel de la serie temporal.

20
2.4.1. Holt-Winters aditivo
Para una serie de tiempo con tendencia, y una estación constante, se establecen las
siguientes ecuaciones de actualización de sus componentes:

Pronóstico:
Ẑt+h|t = at + hbt + st−m+h(k+1) (19)

Nivel:
at = α(Zt − st−m ) + (1 − α)(at−1 + bt−1 ) (20)

Tendencia:
bt = β(at − at−1 ) + (1 − β)bt−1 (21)

Estacionalidad:
st = γ(Zt − at−1 − bt−1 ) + (1 − γ)st−m (22)

En esta expresión, la variable ”k” corresponde a la parte entera de la ecuación (18),

asegurando ası́ que las estimaciones de los ı́ndices estacionales empleadas en el proceso
de pronóstico estén fundamentadas en los datos del último año de la muestra.

En la ecuación (19), se observa un promedio ponderado entre la observación ajustada

estacionalmente (Zt − st−m ) y el pronóstico no estacional (at−1 + bt−1 ) para el tiempo t.
La ecuación (20) es idéntica al método lineal de Holt. En la ecuación (21), se muestra
un promedio ponderado entre el ı́ndice estacional actual (Zt − at−1 − bt−1 ) y el ı́ndice
estacional correspondiente al mismo periodo de la temporada del año anterior (es decir,
m periodos de tiempo atrás).

2.4.2. Holt-Winters multiplicativo

Adecuado cuando la serie tiene una tendencia, al menos localmente, y un patrón
estacional creciente. Las ecuaciones de actualización son:
Pronóstico:
Ẑt+h|t = (at + hbt )st−m+h(k+1) (23)

Nivel:  
Zt
at = α + (1 − α)(at−1 + bt−1 ) (24)
st−m
Tendencia:
bt = β(at − at−1 ) + (1 − β)bt−1 (25)

Estacionalidad:  
Zt
st = γ + (1 − γ)st−m (26)
at−1 + bt−1

21
2.4.3. Estimación de α, β y γ e indicadores de error del modelo
En el estudio realizado por Hudson y Cohen (2016), se menciona que los coeficien-
tes de ponderación α, β y γ se estiman inicialmente estableciendo valores iniciales y
luego se ajustan mediante un proceso de optimización iterativo para lograr un ajuste
adecuado. Comúnmente, se utiliza la minimización del Error Cuadrático Medio (MSE)
o el Error Estándar de Estimación (σ̂) como objetivo de optimización. Estos indicado-
res se utilizan para comparar la efectividad del modelo utilizado, como se describe en
el estudio de Yang et al. (2017). La precisión de la predicción es un criterio de gran
importancia para evaluar la validez del modelo predictivo.

2.5. Métricas de comparación

1. Error Absoluto Medio (MAE): Se define como la media de los errores abso-
lutos de predicción, dado por la ecuación (26):
T
1X
M AE = eM AE = |Zt − Ẑt | (27)
T t=1

2. Raı́z del Error Cuadrático Medio (RMSE): Se define como la raı́z cuadrada
del error cuadrático medio, calculado como por la ecuación (28):
T
1X
M SE = (Zt − Ẑt )2 (28)
T t=1
v
u
u1 X T
RM SE = eRM SE = t (Zt − Ẑt )2 (29)
T t=1

3. Error Porcentual Absoluto Medio (MAPE): Se define como el error absoluto

medio expresado como porcentaje, dado por la ecuación (29):
T
1 X Zt − Ẑt
M AP E = eM AP E = (30)
T t=1 Zt

NOTA: El valor ideal de estas métricas serı́a 0, lo que indicarı́a que no hay error de
predicción.

22
3. MATERIALES, MÉTODOS Y PROCEDIMIENTOS
3.1. Tipo, nivel y diseño de la investigación
3.1.1. Tipo de investigación
Este trabajo de investigación es de tipo aplicado. Según Murillo (2008), la investiga-
ción aplicada también se conoce como investigación práctica o empı́rica, caracterizada
por la aplicación de conocimientos adquiridos y la generación de nuevos conocimientos
a través de la implementación y estructuración práctica basada en la investigación.

3.1.2. Nivel de investigación

El nivel de investigación es explicativo, dado que se emplean técnicas cuantitativas
para el análisis de datos históricos con el fin de predecir el comportamiento futuro
de una variable. Schwarz (2017) destaca que este enfoque permite una comprensión
completa del problema investigado.

3.1.3. Diseño de la investigación

El diseño de la presente investigación es cuantitativo, ya que se basa en modelos
estadı́sticos y matemáticos para el pronóstico de series de tiempo. Es un diseño no
experimental, según Álvarez Risco (2021), porque no se manipulan las variables por
parte del investigador.

3.1.4. Población y muestra

La población de estudio comprende la producción de plata (en Kg) del departamento
de Áncash, con datos mensuales desde enero de 2001 hasta enero de 2024.

3.1.5. Técnicas de recolección de datos

La recolección de datos se realizó mediante la recopilación de datos secundarios
disponibles en la página web del BCRP. La base de datos se puede consultar en el si-
guiente enlace: BCRP - Producción de Plata en Áncash ([Link]
.[Link]/estadisticas/series/mensuales/resultados/RD12997DM/html). El data-
set consta de 289 observaciones de una única variable.

3.1.6. Descripción de la variable

La variable de estudio es la producción total de plata del departamento de Áncash,

medida en kilogramos finos (kg.f). Es una variable de tipo cuantitativo.

3.1.7. Herramientas y software

Se usó el software R versión 4.0.3, para el tratamiento de datos e implementación
del método Box-Jenkis, por otro lado, se usó el programa estadı́stico Excel 2019 para
la implementación del método Holt-Winters.

23
3.1.8. Metodologı́a Box-Jenkins
1. Identificación del modelo:
En esta fase se analizarán los gráficos de la serie y su función de autocorre-
lación para verificar indicios de tendencia y estacionalidad. En caso de que la
serie no sea estacionaria, se aplicarán diferencias para eliminar la tendencia y la
estacionalidad. Para la identificación del modelo se analizarán las funciones de
autocorrelación (ACF) y autocorrelación parcial (PACF) de la serie diferenciada
(libre de tendencia y estacionalidad).
2. Estimación de parámetros:
Después de identificar los valores más apropiados para (p, d y q) de los mode-
los seleccionados, se procede a la estimación de los parámetros autorregresivos
(AR) y de media móvil (MA) incluidos en el modelo. El proceso de estimación
generalmente implica utilizar algoritmos de optimización para encontrar los valo-
res de los parámetros que maximizan la función de verosimilitud. Se analizará la
adecuación del modelo mediante la prueba de significancia de los coeficientes del
modelo y el análisis de los residuos.
3. Validez:
Se verifica si el modelo seleccionado se ajusta de manera razonable a los datos.
Esto implica verificar si los residuos estimados se comportan como ruido blanco,
es decir, que cumplan con la distribución normal, media y covarianza nulas, y
varianza constante. En caso necesario, se ajusta o refina el modelo, lo que puede
implicar probar diferentes órdenes de AR, MA e I en el modelo ARIMA o ajustar
los términos estacionales en el modelo SARIMA.
4. Pronóstico:
Se comparan los pronósticos del modelo estimado con los datos reales en la fase de
prueba utilizando métricas de errores de pronóstico como el error absoluto medio
(MAE), la raı́z del error cuadrático medio (RMSE) y el error porcentual absoluto
medio (MAPE).

24
 Figura 3: Diagrama de flujo de la metodologı́a Box-Jenkins

Fuente: Elaboración propia

3.1.9. Método Holt-Winters

1. Tipo de suavización:
Se evaluará si se optará por el método aditivo o multiplicativo de Holt-Winters
dependiendo de las caracterı́sticas de las fluctuaciones estacionales en los datos.
El enfoque aditivo se prefiere cuando las fluctuaciones estacionales mantienen
una constancia en términos absolutos, mientras que el enfoque multiplicativo es
apropiado cuando las fluctuaciones estacionales cambian en proporción al nivel
de la serie.
2. Estimación de los coeficientes:
Se asignarán valores iniciales adecuados a los coeficientes de suavizado α, β y
γ. Estos parámetros regulan la importancia de los elementos relacionados con el
nivel, la tendencia y la estacionalidad en el modelo. La elección de estos valores
iniciales puede realizarse de forma heurı́stica o a través de técnicas de optimiza-
ción.

25
3. Ajuste:
Se optimizarán los coeficientes de suavizado α, β y γ a través de un proceso ite-
rativo orientado a reducir el error de ajuste con respecto a los datos históricos.
Esto puede llevarse a cabo empleando técnicas de optimización tales como la mi-
nimización del Error Cuadrático Medio (MSE) o el Error Estándar de Estimación
(SEE).
4. Predicción y validez:
Se compararán los pronósticos del modelo estimado con los datos reales en la fase
de prueba utilizando métricas de errores de pronóstico como el error absoluto
promedio (MAE), la raı́z del error cuadrático medio (RMSE) y el error porcentual
absoluto medio (MAPE).

Figura 4: Diagrama de flujo de la metodologı́a Holt-Winters

Fuente: Elaboración propia.

26
4. PRESENTACIÓN Y ANÁLISIS DE RESULTADOS
4.1. Análisis exploratorio
Realizamos un análisis exploratorio de los datos comenzando con la gráfica de la
serie de tiempo de la producción de plata (en Kg) desde enero de 2001 hasta enero de
2024 en Áncash:

Figura 5: Serie de tiempo de la producción de plata en Áncash, Perú

Fuente: Elaboración propia.

En la Figura 5, observamos un incremento significativo a mediados de 2001, seguido

de un perı́odo de estabilidad entre los 30,000 y 50,000 kilogramos aproximadamente
hasta un marcado decrecimiento a finales de 2007, con puntos de decremento muy
pronunciados. El descenso más notable se observa en el primer semestre de 2020, seguido
de un incremento de 50,000 kilogramos.

27
 Tabla 1: Resumen de estadı́sticas descriptivas de la serie de tiempo

Descripción Valor
Número de observaciones 277
Media aritmética 46,713.35
Desviación estándar 13,062.47
Valor mı́nimo 9,724.15
Primer cuartil (Q1) 37,815.44
Mediana (Q2) 46,116.52
Tercer cuartil (Q3) 54,881.65
Valor máximo 82,306.29

Fuente: Elaboración propia.

En la Tabla 1 vemos el resumen de la serie de tiempo , el menor valor es de 9,724.15
kilogramos en enero del 2001, la media de la serie es de 46,713.35 kilogramos.

4.1.1. Suavizado de la serie

Para realizar el suavizado de la serie, usaremos medias moviles centrales, tomando
las observaciones en bloques de 5, es decir, el primer valor que tendremos en la serie
suavizada será el promedio de enero de 2001 a mayo de 2001, y este se vera en marzo
de 2001.
En la Figura 6, podemos ver lo mencionado, las producción suavizada, es la de color
rojo, vemos que todos los picos se reducen, la más notoria es de mayo de 2020.
Figura 6: Comparación de la serie de tiempo y la suavización

Fuente: Elaboración propia.

28
 Usaremos estos valores obtenidos desde marzo de 2020 a julio de 2020 en la serie
original para reducir el pico que tenemos entre esos puntos, graficamente se veria como
en la Figura 7, esta será la serie que usaremos para realizar los modelos.

Figura 7: Serie de tiempo suavizada de la producción de plata en Áncash, Perú

Fuente: Elaboración propia.

En la Tabla 2, tenemos los nuevos valores de la serie, los cuartiles se mantienen, pero
la media, mediana y desviación estándar varian ligeramente.
Tabla 2: Resumen de estadı́sticas descriptivas de la serie de tiempo suavizada

Descripción Valor
Número de observaciones 277
Media aritmética 46,879.42
Desviación estándar 12,920.58
Valor mı́nimo 9,724.15
Primer cuartil (Q1) 37,815.44
Mediana (Q2) 46,022.97
Tercer cuartil (Q3) 54,826.35
Valor máximo 82,306.29

Fuente: Elaboración propia.

29
4.1.2. Gráfico de cajas

En la Figura 8 se observa la distribución de la producción de plata en Áncash a

lo largo de los meses. Se puede notar una clara estacionalidad, con picos más altos
de producción en junio y agosto, mientras que los otros meses muestran cifras de
producción más regulares.

Figura 8: Gráfico de cajas para la distribución mensual de la producción de plata en

Áncash, Perú

Fuente: Elaboración propia

4.1.3. Gráficos de autocorrelación (ACF) y autocorrelación parcial (PACF)

Para realizar el análisis de estacionalidad, se examinaron los correlogramas de la
Función de Autocorrelación (ACF) en la Figura 9 y de la Función de Autocorrelación
Parcial (PACF) en la Figura 10.

Dado que la serie consta de 277 observaciones

√ mensuales, inicialmente se consideró
un rango de rezagos de aproximadamente 277 ≈ 16. Sin embargo, para el modelado
ARIMA, se sugiere ampliar este rango hasta 24 rezagos debido a la periodicidad men-
sual de los datos.

En la Figura 9, se observa que hay una ligera presencia de tendencia en la serie,

como lo indica la autocorrelación decreciente. La primera autocorrelación es cercana a
uno y disminuye notablemente alrededor de 0.3 para el segundo rezago. A partir de allı́,
la autocorrelación muestra decaimiento lento, por lo que es un indicador para realizar
una diferenciación d=1 a la serie de tiempo.

30
 Figura 9: Correlograma de la producción de plata en Áncash, Perú

Fuente: Elaboración propia.

En la Figura 10, se identifican tres rezagos significativos que sobresalen de las ban-
das de confianza, lo que sugiere la posible inclusión de términos autorregresivos en el
modelo. Por lo tanto, se considera apropiado incorporar un componente autorregresivo
(AR) de hasta orden 2 en el modelo.

Figura 10: Correlograma parcial de la producción de plata en Áncash, Perú

Fuente: Elaboración propia.

31
 Realizamos la descomposición de la serie de tiempo para poder visualizarlo mejor.

Figura 11: Descomposición de los componentes de la serie de tiempo

Fuente: Elaboración propia.

En la descomposición de la serie (Figura 11), vemos que se presenta una ligera

tendencia al igual que en el ACF, se ve tambien que hasta inicios del año 2019 tenemos
una tendencia creciente, que en ese año empieza decrecer; esto poniendolo en el contexto
mundial, nos indica que el COVID-19 afecto a la producción de plata en Áncash. La
serie presenta tambien estacionalidad.
Para los métodos de predicción de modelos que desarrollaremos no se debe de eliminar
la estacionalidad, por lo que se realizara la prueba de estacionariedad, sin eliminar la
estacionalidad.

4.1.4. Prueba de Estacionariedad

Para realizar los modelos predictivos es necesario que la serie sea estacionaria, pa-
ra verificarlo realizaremos el test de Dickey-Fuller (ADF), para poder saber cuántas
diferenciaciones necesita la serie hasta que se cumpla la estacionariedad. Planteamos
nuestras hipótesis:

Ho : La serie tiene una raı́z unitaria, es no estacionaria.

H1 : La serie no tiene una raı́z unitaria, es estacionaria.

32
 En la Tabla 3, se observara 6 rezagos y se aplicó la Prueba ADF (Dickey-Fuller Au-
mentada), cuyo resultado mostró un P-valor menor que 0.05. Esto indica que podemos
rechazar la hipótesis nula y concluir que la serie es estacionaria.

Tabla 3: Resultados Test Dickey-Fuller

Descripción Valor
Dickey - Fuller -3.6641
Orden de rezagos 6
P-valor 0.02739

Fuente: Elaboración propia

Por lo que se ve en la Figura 10, la Función de Autocorrelación Parcial nos da 3

rezagos más significativos identificados, asi que también realizaremos el test de Dickey-
Fuller:

Tabla 4: Resultados Test Dickey-Fuller con tres rezagos

Descripción Valor
Dickey-Fuller -4.7406
Orden de rezagos 3
P-valor 0.01

Fuente: Elaboración propia.

Los resultados muestran que el estadı́stico de prueba de Dickey-Fuller es negativo,

lo cual confirma la presencia de estacionariedad en la serie. Además, el P-valor obtenido
es significativamente menor que 0.05, se rechaza la hipótesis nula (H0 ). Por lo tanto,
podemos concluir que la serie es estacionaria.

4.1.5. Diferenciación
Se observa un decaimiento lento en el gráfico ACF (Figura 9), esto podrı́a ser una
indicación de que la serie es estacionaria en varianza pero no en media. La diferencia-
ción sigue siendo necesaria para eliminar la tendencia de la media y asegurar que la
serie sea completamente estacionaria.
Al usar la función ndiff() en RStudio corroboramos que se puede realizar una diferen-
ciación.

En la Figura 12, presentamos como se grafica la serie de tiempo después de aplicar

la primera diferencia.

33
 Figura 12: Gráfico lineal de la serie de tiempo diferenciada

Fuente: Elaboración propia.

En la Figura 13 se muestran las funciones de autocorrelaciones y autocorrelaciones

parciales donde se corrobora que la serie es estacionaria, en esta figura se puede apreciar
que la ACF cae rápidamente a cero después del primer rezago, indicando un componente
MA(1). La PACF muestra un solo pico significativo en el primer rezago, indicando que
no hay correlaciones adicionales significativas después del primer rezago.

Figura 13: Funciones de autocorrelación y autocorrelación parcial de la serie diferenciada

Fuente: Elaboración propia.

34
4.2. Modelamiento Box-Jenkins
En base a la metodologı́a de Box-Jenkins, seleccionamos cinco modelos, resultado
de numerosas pruebas, obteniéndose los parámetros significativos para tales modelos:

Tabla 5: Estimación de parámetros de los modelos ARIMA y SARIMA

Modelo Parámetro Estimación Error Estándar P-valor

ARIMA(0,1,1) θ1 -0.6167 0.05401 2,2 × 10−16
SARIMA(0,1,1)(0,0,1)12 θ1 -0.6167 0.05354 2 × 10−16
Θ1 0.0459 0.06336 0.4689
SARIMA(0,1,1)(1,0,0)12 θ1 -0.6149 0.05357 2 × 10−16
Φ1 0.0421 0.06085 0.4889
SARIMA(0,1,1)(1,0,1)12 θ1 -0.6159 0.05371 2 × 10−16
Φ1 -0.2223 0.61735 0.7188
Θ1 0.2713 0.61219 0.6576
ARIMA(1,1,1) ϕ1 0.0113 0.0854 0.9319
θ1 -0.6259 0.1166 7,91 × 10−8

Fuente: Elaboración propia.

Después de la estimación de los parámetros aplicaremos la prueba de Ljung-Box

test de autocorrelación para verificar si los residuos del modelo propuestos tienen una
distribución de ruido blanco.
La prueba de Ljung-Box se puede definir de la siguiente manera:

Ho : Los datos se distribuyen de forma independiente.

H1 : Los datos no se distribuyen de forma independiente.

Tabla 6: Resultados Test DE Ljung-Box

Modelos P-valor
ARIM A(0, 1, 1) 0,118
SARIM A(0, 1, 1)(0, 0, 1)12 0.1325
SARIM A(0, 1, 1)(1, 0, 0)12 0.1325
SARIM A(0, 1, 1)(1, 0, 1)12 0.1289
ARIM A(1, 1, 1) 0.1191

Fuente: Elaboración propia.

En la Tabla 6, observamos el p-valor y concluimos que se rechaza la Ho , en los 5

modelos, por lo que sus residuos serán ruido blanco.
Esta afirmación se puede corroborar con el análisis de los residuales de los modelos.

35
 Figura 14: Residuales del modelo ARIMA(0, 1, 1).

Fuente: Elaboración propia.

Figura 15: Residuales del modelo SARIMA(0, 1, 1)(0, 0, 1)12 .

Fuente: Elaboración propia.

36
Figura 16: Residuales del modelo SARIMA(0, 1, 1)(1, 0, 0)12 .

Fuente: Elaboración propia.

Figura 17: Residuales del modelo SARIMA(0, 1, 1)(1, 0, 1)12 .

Fuente: Elaboración propia.

37
 Figura 18: Residuales del modelo ARIMA(1, 1, 1).

Fuente: Elaboración propia.

De acuerdo con las Figuras 14, 15, 16, 17 y 18 se observa que los errores se dis-
tribuyen aleatoriamente. Por tanto, es posible asegurar que tienen varianza constante.
Las funciones de autocorrelación de los modelos no son significativas. El Q-Q de los
residuos muestra una distribución que se aproxima a una normal. Según la prueba de
Box-Ljung, los p-valor son mayores a 0.05, por lo tanto, no se rechaza la Ho con un
nivel de significancia de 0.01, concluyendo que los residuos son independientes, lo que
indica que los residuos se comportan como ruido blanco

Después de realizar los pronósticos, aplicamos medidas de error para evaluar la

precisión de cada modelo SARIMA. Las medidas fueron el MAE, RMSE y MAPE. En
la Tabla 7 se muestran los resultados obtenidos para cada medida en cada modelo.

Tabla 7: Medida de los errores de pronósticos de los modelos ARIMA/SARIMA ajustados

Modelos MAE RSME MAPE

ARIMA(0, 1, 1) 6134.035 7975.458 13.3188
SARIMA(0, 1, 1)(0, 0, 1)12 6116.994 7967.562 13.28295
SARIMA(0, 1, 1)(1, 0, 0)12 6119.375 7968.24 13.28884
SARIMA(0, 1, 1)(1, 0, 1)12 6109.007 7965.297 13.26147
ARIMA(1, 1, 1) 6136.562 7975.289 13.32286

Fuente: Elaboración propia.

38
 El modelo SARIMA(0, 1, 1)(1, 0, 1)12 fue el modelo que obtuvo los menores valores
para las medidas de error de pronóstico MAE, RMSE y MAPE.

Luego de haber validado los residuos de los modelos e identificado cual es el mejor
modelo para el pronóstivo de la producción de plata en el departamento de Áncash-Perú
en el Perú para los meses de febrero a diciembre del 2024, se mostrarán los resultados
del pronóstico en la Tabla 8 y su gráfico en la Figura 19.

Tabla 8: Pronósticos para los meses de febrero a diciembre del 2024

Mes del 2024 Pronóstico

Febrero 44134.87
Marzo 43426.53
Abril 43993.20
Mayo 44211.89
Junio 44367.99
Julio 44025.38
Agosto 43962.44
Setiembre 44040.01
Octubre 44064.37
Noviembre 44483.38
Diciembre 44493.79

Fuente: Elaboración propia

Figura 19: Pronóstico del modelo SARIMA(0, 1, 1)(1, 0, 1).

Fuente: Elaboración propia.

39
4.3. Modelamiento Holt-Winters
En base a la Figura 7 donde se encuentra la serie de tiempo, que muestra la pro-
ducción de plata fina en Ancash entre enero de 2001 y enero de 2024, se planteo aplicar
un modelo Holt-Winters Multiplicativo. Esta elección se fundamenta en las siguientes
observaciones:

Según la Figura 11 de la descomposición de la serie, la tendencia es decreciente a

largo plazo pero no es tan pronunciada como en la serie original
Según la Figura 11 en la descomposición de la serie se presenta un patrón esta-
cional, donde la producción de plata aumenta en determinados periodos del año
y disminuye en otros. Sin embargo, este patrón estacional no es constante a lo
largo del tiempo, sino que aumenta a medida que la producción general también
aumenta.
El patrón estacional creciente sugiere que la estacionalidad no es un efecto aditivo
sobre la tendencia, sino que multiplica el valor de la tendencia. En otras palabras,
la variación estacional es proporcional al nivel general de producción. Esto sig-
nifica que las fluctuaciones estacionales aumentan o disminuyen en proporción al
nivel general de la serie temporal.

El modelamiento comienza con parámetros de suavización aleatorios los cuales se

irán ajustando con iteraciones mientras se minimiza el MAPE:

Tabla 9: Parámetros iniciales

Parámetros Valor
α 0.3
β 0.3
γ 0.3
a1 9724.15
b1 0

Fuente: Elaboración propia.

La preparación de los datos es fundamental para aplicar correctamente las ecuacio-

nes de actualización de los parámetros de Holt-Winters.

40
 Figura 20: Pronósticos iniciales

Fuente: Elaboración propia

Figura 21: Gráfico de lı́neas producción de plata y su pronóstico

Fuente: Elaboración propia

41
 En la Figura 21 se observa que no existe un buen ajuste a los datos por lo que sı́
se tendrá que modificar los valores de los parámetros para obtener los óptimos para el
modelo.
En la Figura 20 visualizamos los valores iniciales marcados en amarillo, sirvieron
como punto de partida para actualizar las fórmulas recursivas de nivel, tendencia y es-
tacionalidad. Tras completar las actualizaciones, se calcularon los errores de pronóstico
y las métricas MAE, RSME y MAPE , considerando también los valores iniciales. Es-
tos valores nos permiten optimizar los parámetros de suavizado, con el objetivo de
minimizar el MAPE bajo las siguientes restricciones:

α ≥ 0 ,β ≥ 0 ,γ ≥ 0
α ≤ 1 ,β ≤ 1 ,γ ≤ 1
Utilizando la función SOLVER de Excel y la técnica de optimización GRG Nonli-
near, se consiguió minimizar el MAPE. La configuración del SOLVER para este proceso
se detalla en la Figura 22.

Figura 22: Optimzación de parámetros con SOLVER

Fuente: Elaboración propia.

42
 Tras minimizar el MAPE, se actualizaron los parámetros de Holt-Winters para el
nivel, la tendencia y la estacionalidad. Los nuevas parámetros actualizados se muestran
a continuación además los pronosticos calculados para los meses restantes de 2024:

Tabla 10: Parámetros optimizados

Parámetros Valor
α 0.4345
β 0
γ 0

Fuente: Elaboración propia.

Figura 23: FEB-DIC 2024 con parámetros optimizados

Fuente: Elaboración propia.

43
Figura 24: Grafico de lineas pronósitcos FEB-DIC 2024 con parámetros optimizados

Fuente: Elaboración propia.

La optimización de los parámetros no solo ha dado como resultado una actualización

de los mismos, sino que también ha traı́do consigo una reducción considerable en los
errores de predicción obteniendo los siguentes MAE, RSME y MAPE.

Tabla 11: Medida de los errores de pronósticos Holt-Winters

Métricas Valor
MAE 5905.7036
RSME 7890.2678
MAPE 12.8895

Fuente: Elaboración propia.

4.4. Comparación de modelos

Tabla 12: Medida de los errores de pronósticos de los modelos

Modelo MAE RSME MAPE

SARIMA(0, 1, 1)(1, 0, 1)12 6109.007 7965.297 13.26147
Holt-Winters 5905.7036 7890.2678 12.8895

Fuente: Elaboración propia.

Se puede observar en la Tabla 12 que el mejor modelo es el Holt-Winters ya que

tiene un menor valor en el MAE, RMSE y MAPE.

44
5. CONCLUSIONES
El mejor modelo del método Box-Jenkins es el SARIMA(0, 1, 1)(1, 0, 1)12 ya que
tiene menores valores de error de pronóstico entre toda su familia, el análisis con
dicha metodologı́a ofrece información relevante sobre el patrón de comportamiento
de la serie de tiempo, dicha información ayuda en la elaboración de los modelos
más eficientes.
El modelo Holt-Winters en su versión multiplicativa, ha mostrado ser efectivo
para capturar las tendencias y estacionalidades presentes en la producción de
plata fina en Ancash con los parámetros de suavizacion del nivel, tendencia y
estación con valores 0.4345, 0 y 0 respectivamente. Con estos valores el resultado
calculado para la métrica de comparación MAPE = 12.8895.
Al realizar la comparación entre los 2 métodos propuestos, se identificó el modelo
del método Holt-Winters multiplicativo demuestra un mejor desempeño en com-
paración con el modelo SARIMA(0, 1, 1)(1, 0, 1)12 en todos los aspectos evaluados:
MAE, RMSE y MAPE. Estos resultados sugieren que el modelo Holt-Winters pro-
porciona pronósticos más precisos y confiables para esta serie temporal especı́fica.
Por lo tanto, se recomienda utilizar el modelo Holt-Winters para futuras predic-
ciones de la producción de plata en la región de Áncash, Perú para el periodo de
febrero a diciembre del 2024.

45
6. BIBLIOGRAFIA
Box, G. E., Jenkins, G. M., Reinsel, G. C., & Ljung, G. M. (2015). Time series
analysis: forecasting and control. John Wiley & Sons.
Chafloque, M., Chung, V., & Cotrina, A. (2018). Modelación de las dinámicas, esti-
mación y análisis de las exportaciones mineras: caso peruano. Alternativas, 19(3),
37-45. ISSN: 1390-1915.
Chatfield, C. (2015). The analysis of time series: an introduction (6th ed.). Chapman
and Hall/CRC.
Chung, A., & Inche, J. (2022). Aplicación de técnicas econométricas y machine lear-
ning para series de tiempo univariadas en una escuela profesional de una univer-
sidad pública. Véritas, 23(1), 47-51.
Garza, D. C., Ovalle, H., & Tovar, O. C. (2022). Modelos estadı́sticos para el pronóstico
de las regalı́as producto de la extracción de minerales en Colombia.
Marino, M., Arango, A., Lotero, L., & Jiménez, M. (2021). Modelos de series tempo-
rales para pronóstico de la demanda eléctrica del sector de explotación de minas
y canteras en Colombia. Revista EIA, Escuela de Ingenierı́a de Antioquia, 18(35),
1-17.
Mejı́a, E., & Gonzales, S. (2019). Predicción del consumo de energı́a eléctrica resi-
dencial de la Región de Cajamarca mediante modelos Holt-Winters. Revista de
Ingenierı́a Eléctrica (RIE).
Minerı́a & Energı́a. (2023). Producción de Plata en Perú 2023: Tendencias, Retos y
Lı́deres del Sector. Recuperado de [Link]
-de-plata-en-peru-2023-tendencias-retos-y-lideres-del-sector/
Mohammadi, H., Nezhad, M. S., & Rafiei, M. (2020). Short-term silver demand fo-
recasting: A Holt-Winters approach with seasonal decomposition. International
Journal of Forecasting, 36(4), 1271-1282.
Nina Yucra, S. M. (2021). Modelo predictivo de producción de minerales de plata
mediante redes neuronales artificiales. INF-FCPN-PGI Revista PGI, (8), 176-179.
Villavicencio, J. (2010). Introducción a series de tiempo. Puerto Rico.

46
7. ANEXOS
7.1. NORMAS, ESTANDARES Y REGULACIONES
Estándares Internacionales de Orientación sobre las Técnicas Estadı́sticas para
la Norma NTP-ISO 9001:2001.

¿Qué es? El análisis de series de tiempo es una familia de métodos para estudiar
una colección de observaciones realizadas secuencialmente en el tiempo. Se utiliza para:
- La búsqueda de patrones de ’retraso’ mediante el análisis estadı́stico de cómo se corre-
laciona cada observación con la observación inmediatamente anterior, y repitiendo esto
para cada periodo de separación sucesivo. - La búsqueda de patrones cı́clicos o estacio-
nales, para entender cómo factores causales en el pasado podrı́an influir repetidamente
en el futuro. - La utilización de herramientas estadı́sticas para predecir observaciones
futuras o para entender cuáles son los factores que han contribuido más a las variacio-
nes en una serie de tiempo.
Aunque las técnicas empleadas en el análisis de series de tiempo pueden incluir gráfi-
cos de tendencias simples, en esta Norma Técnica Peruana tales gráficos elementales
se enumeran entre los métodos gráficos simples citados en la estadı́stica descriptiva.

¿Para qué se utiliza? - El análisis de series de tiempo se utiliza para describir

patrones en los datos de series de tiempo, para identificar valores atı́picos (es decir,
valores extremos cuya validez deberı́a ser investigada), bien para ayudar a entender los
patrones o para hacer ajustes, y para detectar los puntos de cambio en una tendencia.
Otro uso que se le da es para explicar los patrones en una serie de tiempo con los de
otra serie de tiempo, con todos los objetivos inherentes en el análisis de regresión.
- El análisis de series de tiempo se utiliza para predecir valores futuros en serie de
tiempo, habitualmente con algunos lı́mites superiores e inferiores, conocidos como el
intervalo de predicción. Se utiliza ampliamente en el área de control y frecuentemente
se aplica en procesos automatizados. En ese caso, se ajusta un modelo de probabilidad
al historial de la serie de tiempo, se predicen los valores futuros y luego se ajustan
parámetros del proceso especı́ficos para mantener el proceso de acuerdo con lo previs-
to, con tan poca variación como sea posible.
SARIMA (Seasonal Autoregressive Integrated Moving Average): Este mo-
delo es adecuado para series temporales con componentes estacionales y no estacio-
nales. Es robusto para manejar datos con tendencias y estacionalidades complejas.
Holt-Winters: Es ideal para series temporales que muestran patrones estacionales y
tendencias, ajustando los datos mediante tres componentes: nivel, tendencia y estacio-
nalidad.
La comparación entre los modelos SARIMA y Holt-Winters para pronosticar la pro-
ducción de plata en el departamento de Ancash debe considerar los siguientes criterios:
Precisión: La precisión de los pronósticos se mide mediante métricas como el
error medio absoluto (MAE) y el error cuadrático medio (RMSE).
Robustez: La robustez de los modelos se refiere a su capacidad para manejar
datos incompletos o ruidosos.
Interpretabilidad: La interpretabilidad de los modelos se refiere a la facilidad
con la que se pueden entender los resultados y las relaciones entre las variables.

47
 Beneficios y Limitaciones:
Beneficios: Identificación de patrones cı́clicos y estacionales, pronósticos precisos
para planificación y control.
Limitaciones: Selección adecuada del modelo y precisión influenciada por la
inclusión u omisión de datos.

7.2. RESTRICCIONES MULTIPLES, LIMITACIONES DE LA INVES-

TIGACION
La precisión de los modelos predictivos depende en gran medida de la calidad
y cantidad de los datos disponibles. En esta investigación, la disponibilidad de
datos históricos de alta frecuencia y precisión fue limitada,ya que los datos se
encuentran hasta Enero de 2024 lo que podrı́a haber afectado los resultados.
El periodo de predicción se limitó a febrero a diciembre de 2024. Predicciones a
largo plazo podrı́an requerir ajustes adicionales a los modelos o la incorporación
de factores externos no considerados en este estudio.
En esta investigación, se ha demostrado que los modelos SARIMA y Holt-Winters
son herramientas valiosas para predecir la producción de plata en Áncash. Sin
embargo, es importante reconocer que la producción de plata está sujeta a una
variedad de factores externos que no fueron explı́citamente incorporados en estos
modelos.
Cada modelo predictivo tiene sus propias limitaciones metodológicas que deben
ser consideradas al interpretar los resultados. SARIMA este modelo puede no ser
efectivo en series temporales con cambios estructurales abruptos, como cambios
repentinos en las tendencias o la estacionalidad. Holt-Winters este modelo puede
ser menos preciso en la presencia de tendencias no lineales complejas o cuando la
estacionalidad no es constante a lo largo del tiempo.
Es importante destacar que los hallazgos de esta investigación tienen un alcance
especı́fico y no pueden generalizarse directamente a otras regiones o a la produc-
ción de otros metales.

48
7.3. CÓDIGO DEL PROGRAMA

#Cargamos las librerias necesarias

library (tidyverse)
library ( forecast )
library ( tseries )
library (ggplot2)
library (readxl)
library (tidyr)
library (zoo)
library (urca)
library ( viridis )

#Cargamos la base de datos

datos plata <− read excel(”D:/Datos [Link]”,
col types = c(”date”, ”numeric”))
datos plata <− datos plata %> %
mutate(Fecha = [Link](Fecha))

###Analisis Exploratorio
#Descripcion general
summary(datos plata)
head(datos plata)
tail (datos plata)
str (datos plata)
any(is .na(datos plata))

#Convertimos a una serie de tiempo

ts plata <− ts(datos plata$Produccion, start(2001,1), frequency = 12)
#Graficamos la serie de tiempo
ggplot(datos plata, aes(x = Fecha, y = Produccion)) +
geom line() + labs( title = ”Produccion de Plata (Enero 2001 − Enero 2024)”, x
= ”Fecha”, y = ”Produccion”) +
scale x date(date breaks = ”1 years”, date labels = ” %Y”) + theme minimal()

###Suavizado Central
#Aplicamos el suavizado central tomando la media de 5 observaciones
datos plata <− datos plata %> %
mutate(Produccion suavizada = rollmean(Produccion, k = 5, fill = NA, align =
”center”))

#Graficamos la serie de tiempo original y la suavizada

ggplot(datos plata, aes(x = Fecha)) +
geom line(aes(y = Produccion), color = ”blue”, size = 1) +
geom line(aes(y = Produccion suavizada), color = ”red”, size = 1, linetype =
”dashed”) +
labs( title = ”Produccion de Plata con Suavizado Central (Enero 2001 − Enero
2024)”, x = ”Fecha”, y = ”Produccion”) +

49
 scale x date(date breaks = ”1 years”, date labels = ” %Y”) +
theme minimal() +
scale color viridis ( discrete = TRUE) +
theme([Link] = ”bottom”) +
scale y continuous(labels = scales :: comma)

#Reemplazamos los datos de Produccion por los datos suavizados de marzo 2020 a
julio 2020
datos plata <− datos plata %> %
mutate(Produccion = ifelse(Fecha >= [Link](”2020−03−01”) & Fecha <=
[Link](”2020−07−01”), Produccion suavizada, Produccion))

#Creamos una nueva tabla con las columnas Fecha y Produccion

plata suavizado <− datos plata %> %
select (Fecha, Produccion)
# Visualizamos la nueva tabla
print(plata suavizado)

#Convertimos la nueva tabla a una serie de tiempo

plata ts <− ts(plata suavizado$Produccion, start = c(2001, 1), frequency = 12)
#Graficamos la nueva serie de tiempo
ggplot(plata suavizado, aes(x = Fecha, y = Produccion)) +
geom line() +
labs( title = ”Serie de Tiempo de Produccion de Plata Suavizada (Enero 2001 −
Enero 2024)”,
x = ”Fecha”, y = ”Produccion”) +
scale x date(date breaks = ”1 years”, date labels = ” %Y”) +
theme minimal()

### Analisis Exploratorio de la serie suavizada

#Descripcion general
summary(plata suavizado)
head(plata suavizado)
tail (plata suavizado)
str (plata suavizado)
any(is .na(plata suavizado))

# Boxplots
colors <− viridis(12, option = ”C”) #definimos el color
boxplot(plata ts ˜ cycle(plata ts) , xlab = ”Mes”, ylab = ”Produccion de Plata”,
main = ”Boxplot de Produccion de Plata por Mes”, col = colors)
##Funciones de autocorrelacion
#ACF
acf vals <− acf(plata ts, lag .max = 24, plot = FALSE)
plot(acf vals , main = ”Funcion de Autocorrelacion (ACF) hasta 24 rezagos”,
col = ”steelblue”, lwd = 2, ylim = c(−0.2, 1), xlab = ”Rezago”,
ylab = ”Autocorrelacion”)
abline(h = 0, col = ”red”, lwd = 2)

50
abline(h = c(−1.96/sqrt(length(plata ts)) , 1.96/sqrt(length(plata ts))) , col =
”blue”,
lty = 2, lwd = 2)
#PACF
pacf vals <− pacf(plata ts, lag .max = 24, plot = FALSE)
plot(pacf vals , main = ”Funcion de Autocorrelacion Parcial (PACF) hasta 24
rezagos”,
col = ”steelblue”, lwd = 2, ylim = c(−0.2, 1), xlab = ”Rezago”,
ylab = ”Autocorrelacion”)
abline(h = 0, col = ”red”, lwd = 2)
abline(h = c(−1.96/sqrt(length(plata ts)) , 1.96/sqrt(length(plata ts))) , col =
”blue”,
lty = 2, lwd = 2)

#Descomposicion de la serie
plata ts components <− decompose(plata ts)
plot(plata ts components,
xlab = ”Fecha”,
ylab = ”Valor”)

#Test de Dickey − Fuller

adf test <− [Link](plata ts)
print(adf test )
# Test de Dickey − Fuller con 3 rezagos
adf test 3 <− [Link](plata ts , k=3)
print(adf test 3)

#########
## Diferenciacin de la serie suavizada
ndiffs (plata ts) #diferenciacion simple d=1
nsdiffs (plata ts) #diferenciacion estacional D=0

#Graficamos la serie de tiempo suavizada y diferenciada

dfplata <− diff(plata ts , differences =1)
dfplata df <− [Link](
Date = time(dfplata),
Value = [Link](dfplata)
)

ggplot(dfplata df , aes(x = Date, y = Value)) +

geom line(color = ”black”) +
labs( title = ”Serie de tiempo suavizada y diferenciada”,
x = ”Fecha”,
y = ” Produccin ”) +
theme minimal()

51
##ACF Y PACF de la serie suavizada y diferenciada
#ACF
par(mfrow=c(2,2),col=’blue’,
col . axis=’darkblue’,bty=’l’ ,lwd=2,font=4,[Link]=1.2)
acf df<− acf(dfplata, lag .max = 24, plot = FALSE)
plot(acf vals , main = ”Funcin de Autocorrelacin (ACF)”,
col = ”steelblue”, lwd = 2, ylim = c(−0.5, 1), xlab = ”Rezago”,
ylab = ”ACF”)
abline(h = 0, col = ”red”, lwd = 2)
abline(h = c(−1.96/sqrt(length(dfplata)), 1.96/sqrt(length(dfplata))) , col =
”blue”,
lty = 2, lwd = 2)

#PACF
pacf df <− pacf(dfplata, [Link] = 24, plot = FALSE)
plot(pacf vals , main = ”Funcin de Autocorrelacin Parcial (PACF)”,
col = ”steelblue”, lwd = 2, ylim = c(−0.5, 0.2), xlab = ”Rezago”,
ylab = ”PACF”)
abline(h = 0, col = ”red”, lwd = 2)
abline(h = c(−1.96/sqrt(length(dfplata)), 1.96/sqrt(length(dfplata))) , col =
”blue”,
lty = 2, lwd = 2)

###METODO BOX − JENKINS

#Ajustamos los modelos ARIMA/SARIMA
modelo <− [Link](plata ts, seasonal = TRUE, stepwise = TRUE, trace =
TRUE)
print(modelo)
summary(modelo)
coeftest (modelo,df=Inf)
## A partir de lo obtenido generamos nuestros mejor 5 modelos:
modelo 1 <− Arima(plata ts, order = c(0, 1, 1))
print(modelo 1)
summary(modelo 1)
coeftest (modelo 1,df=Inf)
f1 <− forecast(modelo 1,h=11)
accuracy(f1)
modelo 2 <− Arima(plata ts, order = c(0, 1, 1), seasonal = c(0, 0, 1))
print(modelo 2)
summary(modelo 2)
coeftest (modelo 2,df=Inf)
f2 <− forecast(modelo 2,h=11)
accuracy(f2)
modelo 3 <− Arima(plata ts, order = c(0, 1, 1), seasonal = c(1, 0, 0))
print(modelo 3)
summary(modelo 3)
coeftest (modelo 3,df=Inf)
f3 <− forecast(modelo 3,h=11)

52
accuracy(f3)
modelo 4 <− Arima(plata ts, order = c(0, 1, 1), seasonal = c(1, 0, 1))
print(modelo 4)
summary(modelo 4)
coeftest (modelo 4,df=Inf)
f4 <− forecast(modelo 4,h=11)
accuracy(f4)
modelo 5 <− Arima(plata ts, order = c(1, 1, 1))
print(modelo 5)
summary(modelo 5)
coeftest (modelo 5,df=Inf)
f5 <− forecast(modelo 5,h=11)
accuracy(f5)
## Realizamos la prueba de Ljung Box, para los 5 modelos:
ljung box <− [Link](modelo 1$residuals, lag = 12, type = ”Ljung−Box”)
print(”Prueba Ljung−Box:”)
print(ljung box)
ljung box <− [Link](modelo 2$residuals, lag = 12, type = ”Ljung−Box”)
print(”Prueba Ljung−Box:”)
print(ljung box)
ljung box <− [Link](modelo 3$residuals, lag = 12, type = ”Ljung−Box”)
print(”Prueba Ljung−Box:”)
print(ljung box)
ljung box <− [Link](modelo 4$residuals, lag = 12, type = ”Ljung−Box”)
print(”Prueba Ljung−Box:”)
print(ljung box)
ljung box <− [Link](modelo 5$residuals, lag = 12, type = ”Ljung−Box”)
print(”Prueba Ljung−Box:”)
print(ljung box)
# Funcion para generar los graficos 2x2 con titulos diferenciados por modelo
plot residuales <− function(modelo, modelo nombre) {
par(mfrow=c(2,2))
# Grafico de Residuos
plot( residuals (modelo), main=paste(”Residuales del ”, modelo nombre),
ylab=”Residuos”)
# ACF de Residuos
acf( residuals (modelo), main=paste(”ACF de Residuos − ”, modelo nombre))
# PACF de Residuos
pacf( residuals (modelo), main=paste(”PACF de Residuos − ”, modelo nombre))
# Grafico Q−Q de Residuos
qqnorm(residuals(modelo), main=paste(”Q−Q de Residuos − ”, modelo nombre))
qqline( residuals (modelo), col=”red”)
# Prueba de Normalidad de Residuos (Shapiro−Wilk)
shapiro. test ( residuals (modelo))
par(mfrow=c(1,1)) # Restaurar la configuracion original de graficos
}
# Llama a la funcion para cada modelo, especificando el nombre del modelo
plot residuales (modelo 1, ”Modelo 1”)

53
plot residuales (modelo 2, ”Modelo 2”)
plot residuales (modelo 3, ”Modelo 3”)
plot residuales (modelo 4, ”Modelo 4”)
plot residuales (modelo 5, ”Modelo 5”)

#Grafico del pronostico del modelo SARIMA

f4 df <− [Link](
Date = time(f4$mean),
Forecast = [Link](f4$mean),
Lower80 = [Link](f4$lower[, 1]),
Upper80 = [Link](f4$upper[, 1]),
Lower95 = [Link](f4$lower[, 2]),
Upper95 = [Link](f4$upper[, 2])
)
original df <− [Link](
Date = time(plata ts),
Value = [Link](plata ts)
)
ggplot() +
geom line(data = original df , aes(x = Date, y = Value), color = ’black’) +
geom line(data = f4 df, aes(x = Date, y = Forecast), color = ’red’) +
geom ribbon(data = f4 df, aes(x = Date, ymin = Lower80, ymax = Upper80),
alpha = 0.2, fill = ’orange’) +
geom ribbon(data = f4 df, aes(x = Date, ymin = Lower95, ymax = Upper95),
alpha = 0.1, fill = ’red’) +
labs( title = ” Pronstico del modelo SARIMA(0, 1, 1)(1, 0, 1)”, x = ”Fecha”, y
= ”Produccin”) +
theme minimal()

7.4. HOJAS DE DATOS

54
 Fecha Producción en Kg. Nov-2004 30,390.3
Dic-2004 30,474.3
Ene-2001 9,724.2 Ene-2005 51,879.5
Feb-2001 14,146.2 Feb-2005 48,629.7
Mar-2001 12,280.7 Mar-2005 31,855.5
Abr-2001 13,555.8 Abr-2005 46,313.3
May-2001 14,667.5 May-2005 30,497.1
Jun-2001 13,483.5 Jun-2005 45,082.8
Jul-2001 26,463.4 Jul-2005 27,587.2
Ago-2001 34,252.6 Ago-2005 44,514.6
Set-2001 26,429.6 Set-2005 31,623.6
Oct-2001 36,949.5 Oct-2005 32,884.5
Nov-2001 28,567.1 Nov-2005 43,550.5
Dic-2001 38,437.6 Dic-2005 35,078.6
Ene-2002 38,681.6 Ene-2006 34,939.3
Feb-2002 25,657.8 Feb-2006 33,977.3
Mar-2002 35,048.9 Mar-2006 42,945.1
Abr-2002 27,231.4 Abr-2006 37,624.6
May-2002 37,730.7 May-2006 39,823.5
Jun-2002 35,507.6 Jun-2006 38,018.0
Jul-2002 34,068.9 Jul-2006 39,165.4
Ago-2002 38,035.7 Ago-2006 36,182.5
Set-2002 33,366.0 Set-2006 29,386.1
Oct-2002 37,219.9 Oct-2006 35,845.4
Nov-2002 32,717.6 Nov-2006 36,211.0
Dic-2002 36,496.8 Dic-2006 28,539.1
Ene-2003 46,104.7 Ene-2007 30,978.7
Feb-2003 39,384.5 Feb-2007 33,002.9
Mar-2003 39,859.2 Mar-2007 44,413.7
Abr-2003 37,646.1 Abr-2007 43,856.4
May-2003 45,618.0 May-2007 41,820.1
Jun-2003 34,573.3 Jun-2007 35,024.1
Jul-2003 30,300.7 Jul-2007 47,165.2
Ago-2003 34,540.1 Ago-2007 39,016.4
Set-2003 37,815.4 Set-2007 46,829.6
Oct-2003 35,177.6 Oct-2007 36,439.8
Nov-2003 34,194.7 Nov-2007 21,640.1
Dic-2003 33,186.9 Dic-2007 31,524.6
Ene-2004 32,516.8 Ene-2008 31,602.5
Feb-2004 33,892.3 Feb-2008 33,267.4
Mar-2004 36,730.6 Mar-2008 49,550.7
Abr-2004 44,163.7 Abr-2008 53,721.5
May-2004 35,134.1 May-2008 48,395.5
Jun-2004 28,698.4 Jun-2008 49,207.0
Jul-2004 28,668.7 Jul-2008 43,743.8
Ago-2004 48,116.7 Ago-2008 50,203.4
Set-2004 31,610.0 Set-2008 43,622.6
Oct-2004 27,890.8 Oct-2008 41,569.2

55
 Fecha Producción en Kg. Set-2012 52,582.7
Oct-2012 41,386.0
Nov-2008 42,461.5 Nov-2012 52,373.5
Dic-2008 43,733.9 Dic-2012 39,446.2
Ene-2009 45,853.6 Ene-2013 43,800.4
Feb-2009 43,575.6 Feb-2013 40,944.9
Mar-2009 49,695.0 Mar-2013 54,206.4
Abr-2009 52,404.1 Abr-2013 54,881.6
May-2009 55,229.4 May-2013 59,710.4
Jun-2009 57,865.4 Jun-2013 61,069.2
Jul-2009 47,614.9 Jul-2013 58,890.1
Ago-2009 54,826.3 Ago-2013 66,224.9
Set-2009 46,176.4 Set-2013 37,827.6
Oct-2009 54,292.3 Oct-2013 48,402.0
Nov-2009 67,267.9 Nov-2013 74,369.2
Dic-2009 53,293.2 Dic-2013 51,639.1
Ene-2010 57,397.9 Ene-2014 39,386.7
Feb-2010 48,265.9 Feb-2014 47,704.6
Mar-2010 46,766.2 Mar-2014 42,740.1
Abr-2010 58,020.5 Abr-2014 39,117.4
May-2010 51,423.8 May-2014 43,014.1
Jun-2010 60,264.0 Jun-2014 43,032.9
Jul-2010 49,262.7 Jul-2014 46,230.3
Ago-2010 55,036.2 Ago-2014 51,191.8
Set-2010 41,994.1 Set-2014 45,383.9
Oct-2010 56,771.9 Oct-2014 43,811.8
Nov-2010 41,708.2 Nov-2014 46,116.5
Dic-2010 47,284.2 Dic-2014 37,509.7
Ene-2011 40,385.4 Ene-2015 38,175.5
Feb-2011 43,104.6 Feb-2015 37,164.8
Mar-2011 38,076.9 Mar-2015 39,306.3
Abr-2011 36,520.6 Abr-2015 38,049.5
May-2011 44,366.5 May-2015 39,796.6
Jun-2011 32,988.9 Jun-2015 63,329.9
Jul-2011 40,549.2 Jul-2015 75,007.1
Ago-2011 40,040.4 Ago-2015 66,549.1
Set-2011 41,089.2 Set-2015 79,908.7
Oct-2011 56,506.5 Oct-2015 65,160.5
Nov-2011 40,735.3 Nov-2015 60,053.0
Dic-2011 48,814.6 Dic-2015 82,306.3
Ene-2012 46,023.0 Ene-2016 53,864.6
Feb-2012 41,915.4 Feb-2016 81,816.7
Mar-2012 51,968.5 Mar-2016 78,777.1
Abr-2012 45,254.5 Abr-2016 51,035.7
May-2012 43,514.0 May-2016 68,618.3
Jun-2012 43,052.6 Jun-2016 75,343.2
Jul-2012 35,672.2 Jul-2016 57,166.0
Ago-2012 37,823.8 Ago-2016 62,447.5

56
 Fecha Producción en Kg. May-2020 19,128.9
Jun-2020 42,011.5
Set-2016 57,121.5 Jul-2020 46,597.9
Oct-2016 79,654.9 Ago-2020 47,772.8
Nov-2016 56,898.4 Set-2020 53,281.2
Dic-2016 57,761.9 Oct-2020 69,328.0
Ene-2017 63,125.9 Nov-2020 50,839.8
Feb-2017 63,690.5 Dic-2020 57,385.5
Mar-2017 52,731.5 Ene-2021 53,984.9
Abr-2017 66,963.2 Feb-2021 47,443.3
May-2017 78,319.1 Mar-2021 42,803.2
Jun-2017 76,324.5 Abr-2021 43,821.1
Jul-2017 75,155.1 May-2021 51,725.9
Ago-2017 57,672.0 Jun-2021 58,778.4
Set-2017 81,719.2 Jul-2021 59,793.4
Oct-2017 57,973.4 Ago-2021 61,868.4
Nov-2017 60,906.7 Set-2021 60,572.8
Dic-2017 67,249.8 Oct-2021 61,727.1
Ene-2018 53,142.3 Nov-2021 51,080.5
Feb-2018 55,758.1 Dic-2021 60,675.4
Mar-2018 54,355.4 Ene-2022 54,550.9
Abr-2018 55,027.6 Feb-2022 51,490.3
May-2018 60,857.6 Mar-2022 58,563.0
Jun-2018 61,381.0 Abr-2022 53,216.8
Jul-2018 56,208.9 May-2022 59,460.3
Ago-2018 61,183.2 Jun-2022 53,000.6
Set-2018 60,669.9 Jul-2022 51,222.9
Oct-2018 53,730.2 Ago-2022 60,850.3
Nov-2018 52,685.9 Set-2022 54,387.8
Dic-2018 56,316.4 Oct-2022 50,612.4
Ene-2019 46,435.5 Nov-2022 48,533.9
Feb-2019 49,597.0 Dic-2022 45,715.8
Mar-2019 58,209.3 Ene-2023 41,123.3
Abr-2019 48,149.5 Feb-2023 46,831.8
May-2019 63,030.3 Mar-2023 34,440.2
Jun-2019 59,017.3 Abr-2023 44,427.2
Jul-2019 47,991.6 May-2023 49,905.2
Ago-2019 53,938.5 Jun-2023 50,808.6
Set-2019 56,628.0 Jul-2023 43,425.0
Oct-2019 51,169.3 Ago-2023 44,704.0
Nov-2019 54,899.5 Set-2023 44,711.9
Dic-2019 61,484.5 Oct-2023 44,451.6
Ene-2020 64,735.8 Nov-2023 52,822.1
Feb-2020 51,425.5 Dic-2023 51,660.7
Mar-2020 55,644.9 Ene-2024 35,524.6
Abr-2020 25,680.4

57