INTRODUCCIÓN
El pronóstico de condiciones agrometeorológicas, en sentido general, se refiere a las predicciones esperadas que influyen sobre la producción agropecuaria (Caraza y Quintero, 1997). La utilidad de los mismos se orienta en proveer una guía para la adopción de decisiones, las que en el área de la producción agropecuaria resultan de interés dada la dependencia de esta actividad con las condiciones de la atmósfera.
En la ganadería, se han utilizado las salidas de los modelos globales empleados en la predicción del tiempo asociadas con indicadores como el índice de temperatura y humedad (THI) o rangos de temperatura del aire para la evaluación de las condiciones de los animales (Domínguez et al., 2010). Además de ello, muchas de las técnicas de predicción agrometeorológica se fundan en la relación estadística que existe entre las variables dependientes que han de ser estimadas y las variables agrometeorológicas independientes (precipitación, temperatura), o las variables deducidas (índice de humedad del suelo, influencia de la humedad atmosférica).
La estadística inferencial provee varias herramientas para el análisis que resultan de interés para las ciencias meteorológicas. Los modelos autorregresivos integrados de promedios móviles o ARIMA (acrónimo del inglés autoregressive integrated moving average, también referido como método de Box - Jenkins), donde las estimaciones futuras se sustentan únicamente en los datos previos, han sido aplicados en Cuba para la modelación y pronóstico en diferentes ramas del saber. Por ejemplo, en meteorología se encuentran los trabajos de Osés et al., (2017) y Cabanas (2005), mientras que en la rama agropecuaria se hallan, por ejemplo, los trabajos de Sánchez (2013) y Guevara et al., (2004).
El objetivo general de este trabajo es predecir las condiciones de los animales mediante modelos ARIMA de series temporales para el escenario nacional.
MATERIALES Y MÉTODOS
Se emplearon los datos de 64 estaciones meteorológicas pertenecientes a la red de estaciones de observación del Instituto de Meteorología. Igualmente se usaron los Boletines de Indicadores Fundamentales de la ganadería. De la primera, se utilizaron los datos de la empresa pecuaria La Vitrina.
A continuación se listan los indicadores seleccionados:
Producción de leche
Reproductivos:
Relativos a la masa total:
Para la evaluación del estado de los animales se empleó el Índice de Temperatura y Humedad (ITH) mediante la ecuación 1 (Nienabar y Hahn, 2007; Mujica, 2005 citado por Aguilar, 2006) para un periodo de 21 años (1993 2013).
Donde, t d es la temperatura del bulbo seco (oC) y HR representa la humedad relativa expresada en forma decimal. En la tabla 1 aparecen los intervalos adoptados para la valoración.
Tabla 1.
Escala para la evaluación del ITH.
| Valor del ITH | Evaluación | |
|---|---|---|
| Inferior | Superior | |
| 74 | Normal | |
| 75 | 78 | Alerta |
| 79 | 83 | Peligro |
| 84 | - | Emergencia |
Fuente: Nienabar y Hahn (2007) y Mujica (2005).
Se aplicaron técnicas exploratorias para la identificación de patrones y para ello, los datos fueron divididos y a partir de esto, se crearon series temporales como sigue:
Valores diarios (cada año con 365 datos por estación).
Valores semanales (cada año incluyó 52 datos por estación),
Por valores mensuales (cada año incluyó 12 datos por estación),
Valores anuales (un total de 21 puntos),
Por estaciones (valores diarios, semanales, mensuales y anuales para cada una de las 23 estaciones).
Para el análisis y representación se crearon tres grupos: ITH1, que incluyó las estaciones de las provincias más orientales (desde Camagüey hasta Guantánamo), ITH2 (Matanzas hasta Ciego de Ávila) e ITH3 (Pinar del Rio, Mayabeque, La Habana e Isla de la Juventud).
Se procedió al examen visual de los datos, mediante la representación gráfica con el propósito de identificar irregularidades en las series. Posteriormente se pasó a la descomposición de los mismos (tendencia, ciclicidad y estacionalidad). Se aplicaron los test de Dickey - Fuller y KPSS (Kwiatkowski- Phillips - Schmidt - Shin) para identificar la estacional de las series.
La prueba Dickey - Fuller se basó en asumir que la serie se puede aproximar por un proceso AR(1) con tres variantes: media cero, media diferente de cero y tendencia lineal. En este caso se partió de las hipótesis:
La hipótesis nula (2) asume que la serie no es estacionaria. El estadístico de la prueba se denota por p y su distribución bajo H 0 permite calcular los valores críticos, de tal forma que el criterio de rechazo es p < p0.05, con p es valor calculado del estadístico.
Se representó un diagrama estacional como se describe en Hyndman y Athanasopoulos (2013), es decir, semejante a un gráfico de tiempo excepto que los datos se trazaron contra las estaciones en años separados. La elección de un modelo ARIMA óptimo (p,d,q) se efectuó mediante una función AUTOARIMA. Se consideró a p como el parámetro auto regresivo, d al número de fases de diferenciación no estacionales y q al parámetro de medias móviles.
Los pronósticos de objetos STL (por sus siglas en inglés, Standard Triangle Language) se obtuvieron aplicando un método de predicción no estacional a los datos ajustados estacionalmente y re estadificando utilizando el último año del componente estacional. Se empleó el método Holt Winters de la familia de series temporales de suavizado exponencial.
Para la evaluación e iteración, se plotearon los ACF y PAFC para los residuales de los modelos resultantes. Todos los análisis se efectuaron mediante la aplicación de la herramienta R de R Core Team, (2018).
RESULTADOS Y DISCUSIÓN
Obtención de las series temporales para las diferentes escalas de trabajo
Se adicionaron los datos diarios, semanales y mensuales para conformar series de tiempo de la suma de valores de ITH. Para el ajuste de la serie de tiempo semanal, primeramente se truncaron los primeros y últimos registro de las series y se asignaron los nuevos valores de principio y fin.
Chequeo de la estacionalidad de los datos
Para el ajuste a un modelo tipo ARIMA se requiere que las series sean estacionarias, es decir, que su media, varianza y autocovarianza sean invariantes en el tiempo. En la Tabla 2 aparecen los resultados obtenidos para las pruebas seleccionadas.
Tabla 2.
Test de Dickey - Fuller
Como se puede apreciar en cualquiera de los casos, se rechaza la hipótesis nula y todas las series son estacionarias (p < 0.05). Lo que se corrobora con los resultados del test KPSS (Tabla 3).
Tabla 3.
Test de KPSS para el nivel de estacionalidad.
La autocorrelación en las series es un indicador de cómo cada valor de la serie temporal. En este caso se empleó la función acf() para la detección de esta característica.
El panel en la figura 1 muestra los autocorrelogramas, donde se muestran la correlación entre la serie y sus rezagos. Siempre se tiene que ρ(0)=1. Las líneas discontinuas en azul representan las bandas de confianza de ρ(k) de nivel 95% bajo la hipótesis de que la serie es un ruido blanco (incorrelada).
Figura 1.
Autocorrelogramas para las series.
Se observa que dentro de las series de tiempo hay una fuerte correlación entre cada valor y el valor anterior. Esto significa que los valores de las series temporales no son independientes entre sí. También hay una correlación estadísticamente significativa para el segundo y tercer rezagos. En este caso, el retardo predeterminado es donde N es el número de observaciones y m el número de series.
Una forma simple de lidiar con la autocorrelación es la diferenciación. Diferenciar significa restar valores previos de los valores actuales. La diferenciación de orden 1 significa que se sustrae el valor anterior del valor actual. En este caso, todavía hay autocorrelaciones estadísticamente significativas.
Los gráficos mostrados en la Figura 2 muestran las series temporales.
Figura 2.
Series temporales para los valores medios del ITH.
Una vez obtenidas las series temporales se confeccionaron los diagramas estacionales como se describe en Hyndman y Athanasopoulos (2018).
Los resultados de la aplicación del pronóstico automático mediante suavizado exponencial (método Holt Winters), tanto para la escala semanal como mensual, aparecen en el panel representado en la Figura 3. Para el primer lapso la predicción se realizó para 52 semanas (un año); mientras que en el segundo, se determinó usar 36 meses (tres años).