Aplicación del análisis de series temporales a variables atmosféricas

Universidad de Valladolid Facultad de Ciencias Trabajo de fin de grado Aplicación del análisis de series temporales a variables atmosféricas Autora: Lucía Alonso Losa Tutoras: María Pilar Rodríguez del Tío Ana María Burgos Pérez Índice general Resumen 5 Abstract 5 1. Introducción 6 2. Contexto 7 2.1. Tendencia ..................................... 7 2.1.1. Test de Mann-Kendall . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2. Metodología de Box-Jenkins. Modelos ARIMA. . . . . . . . . . . . . . . . . 9 3. Serie Temperatura 11 3.1. Temperatura Promedio Diaria . . . . . . . . . . . . . . . . . . . . . . . . . . 11 3.2. Temperatura Promedio Mensual . . . . . . . . . . . . . . . . . . . . . . . . . 17 3.3. Temperatura Máxima Mensual . . . . . . . . . . . . . . . . . . . . . . . . . . 22 3.4. Temperatura Mínima Mensual . . . . . . . . . . . . . . . . . . . . . . . . . . 25 4. Serie Humedad Relativa 29 4.1. Humedad Relativa Promedio Diaria . . . . . . . . . . . . . . . . . . . . . . . 29 4.2. Humedad Relativa Promedio Mensual . . . . . . . . . . . . . . . . . . . . . . 34 5. Serie Radiación Solar Incidente 39 5.1. Serie Radiación Solar Incidente Acumulada Diaria . . . . . . . . . . . . . . . 39 5.2. Serie Radiación Solar Incidente Promedio Mensual . . . . . . . . . . . . . . . 45 6. Conclusiones 49 7. Anexo 50 Bibliografía 57 1 Índice de figuras 3.1. Diagrama de cajas múltiple de la Temperatura promedio anual . . . . . . . . 13 3.2. Descomposición de la serie diaria Temperatura . . . . . . . . . . . . . . . . . 14 3.3. Gráfico de la serie diaria Temperatura . . . . . . . . . . . . . . . . . . . . . 14 3.4. ACF de la serie diaria Temperatura . . . . . . . . . . . . . . . . . . . . . . . 15 3.5. Diagrama de cajas múltiple de la serie Temperatura media mensual . . . . . 18 3.6. Descomposición de la serie Temperatura media mensual . . . . . . . . . . . . 19 3.7. Gráfico de la serie Temperatura media mensual . . . . . . . . . . . . . . . . 19 3.8. ACF de la serie Temperatura media mensual . . . . . . . . . . . . . . . . . . 20 3.9. Predicciones de la serie Temperatura media mensual . . . . . . . . . . . . . . 21 3.10. Descomposición de la serie mensual Temperatura máxima . . . . . . . . . . . 23 3.11. Serie Temperatura máxima mensual . . . . . . . . . . . . . . . . . . . . . . . 23 3.12. Predicciones de la serie Temperatura máxima mensual . . . . . . . . . . . . 24 3.13. Descomposición de la serie mensual Temperatura mínima . . . . . . . . . . . 26 3.14. Serie Temperatura mínima mensual . . . . . . . . . . . . . . . . . . . . . . . 26 3.15. Predicciones de la serie Temperatura máxima mensual . . . . . . . . . . . . 27 4.1. Diagrama de cajas múltiple del %Humedad promedio anual . . . . . . . . . 31 4.2. Descomposición clásica de la serie diaria %Humedad............. 32 4.3. Gráfico de la serie diaria %Humedad . . . . . . . . . . . . . . . . . . . . . . 32 4.4. Diagrama de cajas múltiple de la serie %Humedad media mensual . . . . . 35 4.5. Descomposición clásica de la serie %Humedad media mensual . . . . . . . . 36 4.6. Gráfico de la serie %Humedad media mensual . . . . . . . . . . . . . . . . . 36 5.1. Diagrama de cajas múltiple de la Radición solar promedio anual . . . . . . . 41 5.2. Descomposición de la serie diaria Radiación . . . . . . . . . . . . . . . . . . 42 5.3. Gráfico de la serie Radiación diaria . . . . . . . . . . . . . . . . . . . . . . . 42 5.4. ACF de la serie diaria Radiación . . . . . . . . . . . . . . . . . . . . . . . . 43 5.5. Normalidad de los residuales de la serie diaria Radiación . . . . . . . . . . . 44 5.6. Diagrama de cajas múltiple de la serie Radición media mensual . . . . . . . 46 5.7. Descomposición de la serie Radiación media mensual . . . . . . . . . . . . . 47 5.8. Gráfico de la serie Radiación media mensual . . . . . . . . . . . . . . . . . . 47 5.9. ACF de la serie Radiación media mensual . . . . . . . . . . . . . . . . . . . 48 7.1. Test de Ljung-Box de la serie diaria Temperatura . . . . . . . . . . . . . . . 50 7.2. Qqplot de la serie diaria Temperatura . . . . . . . . . . . . . . . . . . . . . . 50 7.3. Test de Ljung-Box de la serie Temperatura media mensual . . . . . . . . . . 51 2 7.4. Diagrama de cajas múltiple de la serie mensual Temperatura máxima . . . . 51 7.5. ACF de la serie Temperatura máxima mensual . . . . . . . . . . . . . . . . . 52 7.6. Test de Ljung-Box de la serie Temperatura máxima mensual . . . . . . . . . 52 7.7. Diagrama de cajas múltiple de la serie mensual Temperatura mínima . . . . 53 7.8. ACF de la serie Temperatura mínima mensual . . . . . . . . . . . . . . . . . 53 7.9. Test de Ljung-Box de la serie Temperatura mínima mensual . . . . . . . . . 54 7.10. ACF de la serie diaria %Humedad....................... 54 7.11. Test de Ljung-Box de la serie % Humedad diaria . . . . . . . . . . . . . . . . 54 7.12. Qqplot de la serie % Humedad diaria . . . . . . . . . . . . . . . . . . . . . . 55 7.13. ACF de la serie %Humedad media mensual . . . . . . . . . . . . . . . . . . 55 7.14. Test de Ljung-Box de la serie%HUmedad media mensual modelo 1 . . . . . 55 7.15. Test de Ljung-Box de la serie%Humedad media mensual modelo 2 . . . . . . 56 7.16. Test de Ljung-Box de la serie Radiación diaria . . . . . . . . . . . . . . . . . 56 7.17. Test de Ljung-Box de la serie Radiación media mensual . . . . . . . . . . . . 56 3 Índice de tablas 2.1. Resultados del Test de Mann-Kendall para los promedios anuales . . . . . . 8 3.1. Estadísticos resumen de la serie anual Temperatura . . . . . . . . . . . . . . 12 3.2. Estadísticos resumen de la serie Temperatura media mensual . . . . . . . . . 17 3.3. Estadísticos resumen de la serie mensual Temperatura máxima . . . . . . . . 22 3.4. Estadísticos resumen de la serie mensual Temperatura mínima . . . . . . . . 25 4.1. Estadísticos resumen de la serie anual %Humedad . . . . . . . . . . . . . . . 30 4.2. Estadísticos resumen de la serie %Humedad media mensual . . . . . . . . . 34 4.3. Estadísticos para la comparación de modelos . . . . . . . . . . . . . . . . . . 37 4.4. Comparación de modelos ARIMA . . . . . . . . . . . . . . . . . . . . . . . . 38 5.1. Estadísticos resumen de la serie diaria Radiación . . . . . . . . . . . . . . . . 40 5.2. Estadísticos resumen de la serie Radiación media mensual . . . . . . . . . . 45 4 Resumen A lo largo de este trabajo se van a modelizar, a través de la metodología ARIMA, series climáticas. En concreto, temperatura, humedad relativa y radiación solar incidente. Los datos de los que se parte son diarios, el análisis de las series se hará tanto con estos datos como con resúmenes mensuales. Para estudiar la posible tendencia, o no, de estas series se utilizará el test de Mann- Kendall, test no paramétrico recomendado por la Organización Mundial de Metereología en el estudio de series ambientales. Los programas utilizados para el tratamiento estadístico de los datos son R y SAS. Palabras claves: ARIMA, tendencia, climatología, Mann-Kendall, R, SAS. Abstract Throughout this document, several climate series will be modelized by using ARIMA methodology. More specifically, parameters such as temperature, relative huidity and incident solar radiation lead these models. The series analysis is based, mainly, on daily data. However, this research will also consider both daily and monthly synthesis. In order to study the hypothetically tendency of these series, the Mann-Kendall test will be conducted, a non-parametric test recommended by World Meteorological Organization for climate series researches. Software used for data and statistical processing will be both R and SAS. Keywords: ARIMA, tendency, climatology, Mann-Kendall, R, SAS. 5 Capítulo 1 Introducción Los datos utilizados en este trabajo consisten en tres series temporales diaras de temperatura (ºC), humedad relativa (%) y radiación solar incidente (decenas de kJ/m2 día) medidas a una altura estándar de 2 metros sobre el suelo. Las tres variables climáticas han sido suministradas por la estación metereológica de la Central Nuclear de Almaraz gracias a un proyecto de transferencia e innovación con la Universidad de Valladolid. La estación meteorológica, ubicada en la provincia de Cáceres, tiene como coordenadas geográficas 39º 48’ 29”N, 5º41’49”W. Los periodos temporales de las series aparecen en la siguiente tabla: Variable Inicio Fin Años Transformación datos horarios Temperatura 01/01/2001 31/12/2019 19 Promedio Humedad Relativa 01/03/2001 31/12/2019 17 Promedio Radiación solar incidente 01/01/1990 31/12/2017 27 Suma Los datos se recogieron cada hora con algunos valores ausentes que fueron imputados. A la autora de este trabajo se le han proporcionado los datos diarios ya resumidos e imputados. En cuanto a la motivación, diversos estudios buscan describir, analizar y predecir cómo será el clima utilizando distintos métodos estadísticos que permitirán conocer variabilidad, tendencias climáticas, y detectar posibles alteraciones climáticas entre otros. Específicamente, el análisis de tendencias de series temporales, es un método bastante aplicado en el estudio de variables climáticas. Entre los contrastes existentes para detectar tendencia está el propuesto por Mann y Kendall, recomendado por la Organización Mundial de Meteorología (OMM) [1]. El propósito de este Trabajo de Fin de Grado es llevar a cabo un análisis de series temporales, utilizando modelos ARIMA, de las variables antes mencionadas, prestando especial atención al análisis de tendencia, como se ha dicho antes, debido al carácter climático de los datos. 6 Capítulo 2 Contexto La Organización Meteorológica Mundial considera importante la aplicación de métodos estadísticos en meteorología, entre ellos el análisis de tendencias y frecuencias, el análisis univariado y multivariado de series temporales y la metodología de Box-Jenkins, entre otros. Según su Guía de instrumentos y métodos de observación [2] las variables meteorológicas son mediciones en puntos discretos en idénticos intervalos temporales en un periodo finito de tiempo. Es decir, situaciones climáticas continuas son discretizadas para formar una serie temporal equiespaciada. 2.1. Tendencia Si estas mediciones se grafican en función del tiempo suponen una herramienta cualitativa muy útil para detectar variaciones, dicho de otra forma, dectectar tendencias. En términos estadísticos podemos definir la tendencia como el comportamiento o movimiento de una serie a largo plazo, si los datos presentan crecimiento o decrecimiento general a lo largo del tiempo. Volviendo a la climatología, y según la Guía de prácticas climatológicas [3], la tendencia es una característica interesante puesto que confiere un resumen del comportamiento histórico de las observaciones. Puntualizar, que en periodos de tiempo más reducidos, lo que pueden parecer tendencias sostenidas podrían pertenecer a oscilaciones más lentas relacionadas con variaciones multi-decadales. En el análisis de tendencia en series ambientales se utiliza frecuentemente el test no paramétrico de Mann-Kendall. En los sucesivos análisis se va a utilizar el comando de R SeasonalMannKendall() para series mensuales del paquete Kendall [4] basado en el test de Mann-Kendall. 7 2.1.1. Test de Mann-Kendall El estadístico utilizado [5] en el contraste es S= n−1 X k=1 n X j=k+1 signo(xj−xk)donde signo(xj−xk) =      1si (xj−xk)>0 0si (xj−xk)=0 −1si (xj−xk)<0 Valores positivos del estadístico S apoyan la hipótesis de tendencia creciente y valores muy negativos apoyan la hipótesis de tendencia decreciente. Su varianza puede aproximarse del siguiente modo cuando se dispone del suficiente número de observaciones, n, y sin tener en cuenta coincidencias. V ar(S) = 1 18[n(n−1)(2n+ 5)]. Las hipótesis del contraste son H0: no hay tendencia vs. H1: hay tendencia creciente/decreciente, y la regla de decisión se establece a partir del valor de ZMk. ZMK =        S−1 √(V ar(S)) si S > 0 0si S= 0 S+1 √(V ar(S)) si S < 0. Si el test de Mann-Kendall proporciona un valor significativo, es decir un p-valor, inferior a 0.05/ 0.1 se puede decir que se rechazada la hipótesis nula por lo que las probabilididades de que sea cierta la presencia de tendencia son altas. A continuación, se probará dicho test mediante el comando MannKendall(), del paquete Kendall también, sobre los promedios anuales calculados a partir de los datos diarios primarios de las variables propuestas. Variable Años Estadístico P-valor Temperatura 19 0.0526 0.77957 Humedad Relativa 17 -0.412 0.023476 Radiación solar incidente 27 0.45 0.00084126 Tabla 2.1: Resultados del Test de Mann-Kendall para los promedios anuales A la vista de la anterior tabla se puede afirmar que la temperatura es la única serie que no tiene tendencia. En cuanto a la humedad relativa y a la radiación solar sí tienen tendencia, y fijando la atención en los signos estadísticos se puede asegurar que la tendencia de la humedad es decreciente y la tendencia de la radiación solar es creciente. Por tanto, se puede concluir en base a estos resultados, que no se observan tendencias significativas en la temperatura dentro del periodo de tiempo estudiado, en cambio, cada vez la humedad es menor y la radiación solar mayor. 8 de la tendencia del apartado anterior sí se vea cierta tendencia creciente. Esta hipótesis de existencia de tendencia se comprobará en el siguiente apartado una vez se haya modelizado la serie mediante la metodología ARIMA. La necesidad de incluir, o no, la constante en el supuesto modelo propuesto determinará, o no, la existencia de cierta tendencia en la serie. Estacionaridad y Modelo ARIMA propuesto Con motivo de poder estudiar la estacionaridad de la serie se va a representar la función de autocorrelación (ACF). Figura 3.4: ACF de la serie diaria Temperatura En vista del ACF se puede observar que por la parte regular decrece muy rápidamente por lo que es estacionaria por la parte regular. En cambio, si nos fijamos en los retardos múltiplos de trescientos sesenta y cinco podemos ver que decrecen muy lentamente por tanto por la parte estacional la serie no es estacionaria. Para solucionarlo se va a diferenciar la serie de forma estacional, ∇365xt=xt−xt−365, y posteriormente realizando el test de Dickey–Fuller, que resulta un p-valor significativo, se confirmará que la serie diferenciada es estacionaria por la parte regular y por la parte estacional. Después de estimar distintos modelos sugeridos observando la ACF y la PACF de la serie diferenciada, y comparando los modelos que mejor se validaban, se propone como modelo ARIMA(2,0,2)(0,1,1)365 sin constante y ecuación (1 −φ1B−φ2B2−B365 +φ1B366 + φ2B367)Xt= (1 −θ1B−θ2B2−ΘB365 +θ1ΘB366 +θ2ΘB367)at. Estimando los parámetros con SAS obtenemos la ecuación Xt= 1,65237Xt−1−0,65802Xt−2+Xt−365 −1,65237Xt−366 + 0,65802Xt−367 +at−069991at−1−0,26089at−2−0,96247at−365 +0,67364at−366 +0,2511at−367. Volviendo al apartado anterior, debido a la ausencia de constante en el modelo propuesto, así como en otros posibles modelos alternativos, se concluye la ausencia de tendencia en la temperatura promedio diaria. 15 Como muestra del grado de validación del modelo se representan los test de Ljung-Box (figura 7.1 del anexo). Para algunos retardos se puede ver que los p-valores están entre 0,01 y 0,05, por lo que la hipótesis de incorrelación de los residuales no se satisface todo lo bien que sería deseable. Para estudiar la distribución de los residuos, se va a tener en cuenta su qqplot (figura 7.2) para comprobar su normalidad, y debido a la ausencia de cambios bruscos en la pendiente de los puntos se va a suponer que los residuales son normales. Para determianar si dicha normal tiene media cero se realizará el contraste H0 : µ= 0, obteniendo un p-valor de 0.7622 por lo que se concluye que los residuales siguen una normal de media cero. En cuanto a la capacidad de predicción, se va a utilizar como medida relativa el error porcentual absoluto medio calculado como %EAM = p X i=1 |ˆyi−yi| p∗yi donde pcorresponde al número de periodos predichos, en este caso solo se va a predecir un periodo por lo que p= 365. Realizando una función muy sencilla en R y exportando los valores predichos desde SAS se obtiene un error del 21.90337 %. 16 3.2. Temperatura Promedio Mensual A continuación, se va a analizar la serie media mensual. Para dicha transformación se ha utilizado el comando de R aggregate(), la función divide los datos en subconjuntos, calcula estadísticas de resumen para cada subconjunto y devuelve el resultado en un grupo por formulario. En este caso se han aprovechado estas características para conseguir la serie mensual promediando en función de las variables mes y año la serie diaria. Por último se ha eliminado la variable día ya que ahora no tiene ningún interés. Análisis Descriptivo La variable temperatura promedio mensual oscila entre 4.988 ºC y 30.528 ºC con una media de 16.361 ºC, si se tienen en cuenta las mensualidades se obtienen los siguientes estadísticos. Mes Mínimo Media Mediana Máximo Rango Enero 4.988 7.925 7.835 10.55 5.56 Febrero 7.326 9.361 9.338 11.23 3.90 Marzo 10.97 12.45 12.55 14.83 3.86 Abril 13.07 15.3 15.11 17.83 4.76 Mayo 16.49 19.51 20.22 21.76 5.27 Junio 21.56 25.18 25.33 28.12 6.56 Julio 25.9 28.04 28.1 30.53 4.63 Agosto 25.52 27.8 28.2 29.73 4.21 Septimbre 21.57 23.4 23.49 26.13 4.56 Octubre 15.99 17.8 17.89 20.3 4.31 Noviembre 8.86 11.53 11.97 13.62 4.76 Diciembre 5.87 8.08 8.1 10.82 4.95 Tabla 3.2: Estadísticos resumen de la serie Temperatura media mensual Fijando la atención en la tabla anterior destaca la disminución tan importante que presentan los rangos en comparación con los de la serie diaria, esto informa de la estabilidad de las temperaturas medias mensuales independientemente de la existencia de años más o menos calurosos. Respecto a los valores de los máximos y los mínimos se puede ver reflejados en ellos el avance de las estaciones, es decir, el ciclo estacional de la temperatura. A continución se hará una representación gráfica donde se puedan ver reflejados los resultados anteriores. 17 Figura 3.5: Diagrama de cajas múltiple de la serie Temperatura media mensual En el diagrama se puede ver claramente lo antes comentado, el ciclo estacional de la temperatura. Destacar como desde los primeros meses del año hasta los meses de verano la temperatura aumenta más lentamente pero, en cambio, desde el fin del verano hasta el comienzo del invierno las temperaturas bajan de manera más abrupta. Respecto a los rangos intercuartiles se puede ver que son más pequeños que los del anterior diagrama, indicando una menor dispersion en las mediciones. Resaltar mayo y junio, donde su rango intercualtílico y sus bigotes respectivamente son especialmente diferentes en comparación con el resto. Terminando los comentario sobre el diagrama de cajas, ver la escasez de valores atípicos. Únicamente presentes en septiembre y dos en diciembre, correspondientes a los años 2018, 2002 y 2019 respectivamente. Por último, en este apartado, se realizará la representación de la descomposición de la serie temperatura promedio mensual. 18 Figura 3.6: Descomposición de la serie Temperatura media mensual Análisis de Tendencia Para comenzar el análisis de la serie se va a representar. Figura 3.7: Gráfico de la serie Temperatura media mensual A la vista del gráfico de la serie no hay una tendencia que se vea a simple vista. Fijándose en la tendencia descrita en la descomposición del apartado anterior se puede ver una forma similar a la descrita en la serie diaria por lo que se sospecha una ausencia de tendencia. A 19 modo de comprobación se va a realizar la variación del test de Mann-Kendall antes descrita, trás su ejecución se obtiene un p-valor de 0.1083 por lo que no se puede rechazar la ausencia de tendencia. Una vez analizada la tendencia de la serie promedio mensual, se ha apreciado que la tendencia de la parte superior y de la parte inferior del gráfico difieren algo, en concreto, se sospecha que las mediciones máximas no tienen una tendencia evidente, en cambio las mínimas sí que se puede apreciar una leve tendencia creciente por lo que las temperaturas mínimas han podido sufrir un aumento. Estacionaridad y Modelo ARIMA propuesto El siguiente punto a observar es la estacionaridad, para ello, primeramente, se va a calcular la función de autocorrelación (ACF). Figura 3.8: ACF de la serie Temperatura media mensual El ACF de esta serie a simple vista puede parecer diferente al de la temperatura media diaria pero prestando atención se pueden extraer las mismas conclusiones. Por la parte regular los retardos decrecen rápidamente por lo que sí es estacionaria, en cambio, por la parte estacional no es estacionaria, ya que en los retardos múltiplos de doce se puede ver claramente decrecimiento lineal, y no exponencial. Para hacer la serie estacionaria por ambas partes se va a diferenciar estacionalmente como en el caso anterior, ∇12xt=xt−xt−12. Una vez diferenciada la serie se aplica el test de Dickey–Fuller y se obtiene un p-valor significativo por lo que se puede concluir que la serie es estacionaria tanto por la parte regular como por la parte estacional. Estas observaciones también se pueden extraer del ACF de la serie diferenciada donde se ve un decrecimiento rápido de los retardos por ambas partes. En este caso se propone el modelo ARIMA(1,0,0)(1,1,1)12 sin constante cuya ecuación sería (1 −φB −(1 + Φ)B12 + (φ+φΦ)B13 + ΦB24 −φΦB25)Xt= (1 −ΘB12)aty estimando sus coeficientes Xt−0,39Xt−1= 0,8736039Xt−12 −0,3888372Xt−13 + 0,1263961Xt−24 − 20 0,04915Xt−25 +at+ 0,9354965at−12. En este modelo, como en el anterior, la constante no es significativa por lo que se tiene otra evidencia de que la serie no tiene tendencia. Con respecto a la validación del modelo, puede verse el gráfico resúmen de los contrastes de Ljung-Box (figura 7.3 del anexo), donde todos los p-valores son notoriamente no significativos por lo que no se rechaza la hipótesis nula y consideramos los residuos como incorrelados. Seguidamente, se comprobará que los residuales sigan una distribución normal de media cero. Por un lado se hará el test de normalidad de Shapiro-Wilk obteniendo un p-valor igual a 0.5177. Por otro lado se contrastará H0 : µ= 0 obteniendose un p-valor de 0.8435. Con estos resultados no podemos rechazar que los residuales sigan una normal de media cero y se concluye que el modelo es válido, y tanto las predicciones como las bandas de predicción podrían considerarse fiables. Por último, con el fin de medir la capacidad de predicción, y reservando dos años completos de datos (24 observaciones), se va a calcular el error porcentual absoluto medio. Esta vez todos los calculos se han realizado con R y con la ayuda del comando forecast() del paquete Forecast para realizar las predicciones. Finalmente se obtiene un error del 6.718701 %. Representando los datos reservados junto a las predicciones en dos colores distintos se obtiene el siguiente gráfico. Figura 3.9: Predicciones de la serie Temperatura media mensual 21 3.3. Temperatura Máxima Mensual Como se ha dicho en la anterior sección, se van a analizar los extremos mensuales. Primeramente se va a estudiar la serie de máximos mensuales, para ello se han creado los datos utilizando el comando aggregate() anteriormente explicado. La diferencia radica en la especificación de la función max() como argumento, además de indicar su aplicación respecto de las variables mes y año de la serie diaria. Análisis Descriptivo Los límites de esta variable son 9.37 ºC y 34.82 ºC, además su media toma un valor de 22.30 ºC. A continuación, se comenzará con la tabla resumen de los datos creados en función de los meses. Mes Mínimo Media Mediana Máximo Rango Enero 9.37 13.39 14.40 15.53 6.16 Febrero 11.72 13.74 13.65 15.70 3.98 Marzo 14.58 16.77 16.55 20.01 5.43 Abril 16.58 20.4 20.68 22.86 6.28 Mayo 21.01 25.91 26.42 30.12 9.11 Junio 25.27 31.08 31.64 34.12 8.85 Julio 30.60 32.31 32.14 34.12 3.52 Agosto 29.83 32.00 31.70 38.82 4.99 Septiembre 25.37 28.32 28.23 32.43 7.06 Octubre 19.62 22.99 23.11 25.98 6.36 Noviembre 11.84 16.67 16.79 20.00 8.16 Diciembre 11.27 14.01 13.82 16.37 5.10 Tabla 3.3: Estadísticos resumen de la serie mensual Temperatura máxima En la anterior tabla se vuelven a ver las mismas características que en la mensual, únicamente destacable las diferencias entre algunos valores de los rangos. Realizando un diagrama de cajas múltiples como en otros casos, dicha variabilidad en el rango supone un aumento en la dispersión de los valores que se verá reflejado en la amplitud de los intervalos intercuartílicos o en la amplitud de los bigotes (figura 7.4 del anexo). Como viene siendo costumbre, se terminará este apartado realizando la descomposición de la serie. 22 Figura 3.10: Descomposición de la serie mensual Temperatura máxima Análisis de Tendencia Para realizar el análisis de tendencia se va a proceder a la representación de la serie en función del tiempo. Figura 3.11: Serie Temperatura máxima mensual 23 En el gráfico de la serie no se aprecia una tendencia general clara. Volviendo a la descomposición del apartado anterior, la descripción de la tendencia vuelve a ser muy parecida a las de otras secciones de este capítulo por lo que probablemente la conclusión sea la misma: no hay tendencia. Realizando la variante del test de Mann-Kendall descrita se obtiene un p-valor de 0.1515 por lo que no se rechaza la hipótesis nula de ausencia de tendencia como se había previsto. Estacionaridad y Modelo ARIMA propuesto Representando la función de autocorrelación (figura 7.5 del anexo) se ha podido ver claramente la misma composición que en la ACF de la serie de promedios mensuales, por lo que las conclusiones sobre la estacionaridad de la serie son las mismas. Es estacionaria regularmente pero no estacionalmente. Habiendo diferenciado por la parte estacioanl la serie se propone el modelo ARIMA(0,0,2) (1,1,1)12 sin constante y con ecuación (1 −(1 + Φ)B12 + ΦB24)Xt= (1 −θ1B−θ2B2− ΘB12 +φ1ΘB13 +θ2ΘB14)at. Estimando los parámetros se obtiene Xt=−0,1597689t−12+ 0,1597689Xt −24 + at−0,1662011at−1−0,1972976at−2+ 0,7910377at−12 −0,7910377at−13 − 0,5937401at−14. Para validar el modelo se va a comenzar por ver los resultados del test de Ljung-Box (figura 7.6 del anexo), que al obtenerse p-valores no significativos se concluye que los residuales del modelo son incorrelados. Pasando a comprobar la distribución de los residuos se realiza el test de Shapiro-Wilks, al obtenerse un p-valor de 0.00258 se puede rechazar que los residuales sigan una normal. Por último, se comprueba si la media de la normal es igual a cero mediante un contraste que resulta 0.9654 por tanto los residuales no siguen una normal pero su distribución sí que tiene media cero. El error porcentual absoluto medio cometido por este modelo es del 5.299361 %, y representando las predicciones junto a sus intervalos se obtiene Figura 3.12: Predicciones de la serie Temperatura máxima mensual 24 Figura 4.1: Diagrama de cajas múltiple del %Humedad promedio anual A la vista del diagrama, llama la atención lo comentado anteriormente, las similitudes en los límites superiores y lo cambiantes que son los límites inferiores. Respecto a los máximos, resaltar la excepción del año 2012, que volviendo a la tabla, se puede ver un valor máximo más bajo. Viendo el diagrama de manera más general, se puede apreciar como los rangos intercuartilicos están dispuestos, ligeramente, en una trayectoria decreciente. Esta observación podría indicar tendencia, más adelante, en esta misma sección, se estudiará con más detenimiento. Seguidamente, se representará la descomposición de la serie. 31 Figura 4.2: Descomposición clásica de la serie diaria %Humedad Análisis de Tendencia A continuación, se representará gráficamente la serie. Figura 4.3: Gráfico de la serie diaria %Humedad En el gráfico de la serie, se aprecia una leve tendencia decreciente. Fijándose en la tendencia que describe la descomposición se ve más claramente ese decrecimiento. Para encontrar alguna evidencia de esta tendencia vista gráficamente habrá que esperar al siguiente apartado. Como en la anterior serie diaria, la necesidad de incluir, o no, constante 32 en el modelo nos proporcionará la confirmación, o la rectificación. Estacionaridad y Modelo ARIMA propuesto La función de autocorrelación de la serie humedad relativa media diaria (figura 7.10 del anexo) es muy similar al de la serie media diaria antes vista, la temperatura, por lo que se pueden extraer las mismas conclusiones: la serie es estacionaria por la parte regular pero no lo es por la parte estacional. Para corregir esa falta de estacionaridad estacional se va a diferenciar la serie por dicha parte. Haciendo uso del test de Dickey–Fuller se ha comprobando que la diferenciación sí ha conseguido el objetivo, debido al p-valor resultante muy significativo, se concluye que la serie es estacionaria regular y estacionalmente. Acto seguido, se propone el modelo ARIMA(2,0,3)(2,1,0)365 sin constante y ecuación (1 −φ1B−φ2B2−(1 + Φ1)B365 + (φ1+φ1Φ1)B366 + (φ2+φ2Φ1)B367 + (Φ1−Φ2)B730 + (φ1Φ2−φ1Φ1)B731)+(φ2Φ2−φ2Φ1)B732 +Φ2B1095 −φ1Φ2B1096 −φ2Φ2B1097)Xt= (1−θ1B− θ2B2−θ3B3)at. Estimando sus parámetros Xt= 1,68578Xt−1−0,68974Xt−2+0,35441Xt−365− 0,59746Xt−366+0,24445Xt−367+0,3076Xt−730−0,51854Xt−731+0,21216X732+0,33799Xt−1095− 0,56977Xt−10960,23313Xt−1097 +at−0,93966at−1−0,13318at−2+ 0,09953at−3. Volviendo al apartado anterior, debido a que la constante no es significativa en el modelo se concluye que la hipótesis de ausencia de tendencia vista se puede suponer cierta. Validando el modelo se ha visto que los test de Ljung-Box (figura 7.11 del anexo) presentan p-valores entre 0,05 y 0,01 por lo que algunos de los residuos pueden tener cierta correlación, por tanto la hipótesis de incorrelación del modelo no satisface como debería. En cuanto a la normalidad de los residuos, el qqplot (figure 7.12 del anexo) no presenta cambios bruscos en la pendiente que describen los puntos, por lo que se van a considerar normales. Junto con el contraste de la media, cuyo p-valor resulta 0,4183, se puede concluir que los residuales son normales de media cero. Reservando un año de observaciones y prediciendo, el error porcentual medio absoluto asciende al 21.3649 %. 33 4.2. Humedad Relativa Promedio Mensual En esta sección se analizará la serie de promedios mensuales de la humedad relativa siguiendo el mismo esquema, y calculada de igual forma, que la serie promedio mensual anterior. Análisis descriptivo En este caso, la variable oscila entre 29.87 %y el 93.15 %con un valor medio del 62.4 %. De igual forma, se van a visualizar estadísticos descriptivos, esta vez mensuales. Mes Mínimo Media Mediana Máximo Rango Enero 78.69 84.32 83.43 93.15 14.46 Febrero 53.03 73.64 74.28 83.79 30.76 Marzo 50.89 66.84 69.13 76.51 25.62 Abril 48.76 65.49 65.31 73.45 24.69 Mayo 43.07 56.38 58.76 76.04 32.97 Junio 35.75 43.59 42.19 56.62 20.87 Julio 29.87 35.06 35.52 39.41 9.54 Agosto 31.08 36.25 35.85 49.72 18.64 Septimbre 36.2 47.20 45.91 63.85 27.65 Octubre 49.06 66.96 69.81 74.71 25.65 Noviembre 73.72 79.74 79.62 88.77 15.05 Diciembre 77.97 85.62 85.81 93.1 15.13 Tabla 4.2: Estadísticos resumen de la serie %Humedad media mensual En vista de la tabla, los rangos presentan diferencias notables entre ellos, registrándose los valores más bajos en los meses de verano y los más altos en los meses de primavera y otoño. Fijando la atención en los máximos y en los mínimos, se pueden ver los ciclos estacionales que se comentaron en los promedios mensuales anteriores, en este caso invertidos. A continuación, se hará una representación gráfica donde se vean reflejados los estadísticos anteriores. 34 Figura 4.4: Diagrama de cajas múltiple de la serie %Humedad media mensual Como ya se ha dicho, se ven claramente los ciclos estacionales de la humedad. Dichos ciclos son contrarios a los vistos en los anteriores promedios mensuales, lo que es razonable ya que los meses de verano suelen ser propias las sequías y las lluvias escasas. También cabe destacar los bigotes tan grandes de los meses de marzo y septiembre, y el rango intercuartílico de mayo especialmente amplio. Por último, destacar la gran cantidad de valores atípicos presentes sobretodo en los meses de primavera y otoño, que volviendo a la tabla, concuerdan con los meses que poseen los rangos más grandes. Como en otras ocasiones, y a modo de comienzo del análisis de la serie temporal, se va a representar la descomposición de la serie. 35 Figura 4.5: Descomposición clásica de la serie %Humedad media mensual Análisis de Tendencia La siguiente figura corresponde con la representación gráfica de la serie. Figura 4.6: Gráfico de la serie %Humedad media mensual A la vista del gráfico se aprecia una leve tendencia decreciente al igual que en la serie de promedios diarios. Fijándose en el subgráfico tendencia de la descomposición sí que se puede ver más evidente esa tendencia decreciente. A modo de comprobación, se va a volver a realizar la variación estacional del test para tendencia de Mann-Kendall. Trás realizar el test 36 se obtiene un p-valor aproximadamente cero y un estadístico negativo, por lo que se puede concluir que existe tendencia, y además decreciente. Estacionaridad y Modelo ARIMA propuesto Representando el ACF (figura 7.13 del anexo) se puede observar que presenta las mismas características que la ACF de la temperatura promedio mensual. Por tanto, igulamente se concluye que la serie es estacionaria regularmente pero no estacionalmente. Realizando el test de Dickey-Fuller, y habiendo diferenciado por la parte estacional, se afirma la estacionaridad de la serie. A continuación, debido a que se obtienen resultados parecidos, se proponen dos modelos. El primero de ellos ARIMA(1,0,0)(2,1,0)12 sin constante y con ecuación (1 −φB − (1 + Φ1)B12 + (φ+φΦ1)B13 + (Φ1−Φ2)B24 + (φΦ2−φΦ1)B25 + Φ2B36 −φΦ2B37)Xt= at. Estimando sus parámetros se obtiene la ecuación Xt= 0,35409Xt−1+ 0,36205Xt−12 − 0,1282Xt−13 +0,24074X24 −0,3665Xt−25 +0,39721X36 −0,14065Xt−37 +at. El segundo modelo es ARIMA(1,0,1)(2,1,0)12 también sin constante y con ecuación (1−φB −(1+Φ1)B12 + (φ+φΦ1)B13 + (Φ1−Φ2)B24 + (φΦ2−φΦ1)B25 + Φ2B36 −φΦ2B37)Xt= (1 −θ)at. Estimando sus parámetros Xt= 0,71013Xt−1+ 0,36326Xt−12 −0,25796Xt−13 + 0,25168Xt−24 − 0,17872Xt−25 + 0,38506Xt−36 −0,27344Xt−37 +at−0,42376at−1. Ambos modelos han sido validados, por lo que sus residuales son incorrelados según los test de Ljung-Box (figuras 7.14 y 7.14 del anexo), y sus bandas de predicción fiables, ya que los residuales se pueden suponer normales de media cero. Algunos de los criterios que se han tenido en cuenta para comparar los modelos son el AIC (Akaike information criterion) y el SBC (Schwarz’s Bayesian Criterion). Se considerará como modelo 1 ARIMA(1,0,0)(2,1,0)12 y como modelo 2 ARIMA(1,0,1)(2,1,0)12. AIC SBC Modelo 1 1268.244 1278.017 Modelo 2 1266.535 1279.565 Tabla 4.3: Estadísticos para la comparación de modelos Estos estadísticos dan una visión muy general de cuál de los dos modelos es mejor, además vemos que la diferencia entre ambos es bastante pequeña, por lo que con la intención de seguir comparando los modelos, se estudiará su capacidad de predicción. Para esta tarea, se han recogido en una tabla las predicciones del último año, los valores de la serie en ese año, los errores absolutos y el error cuadrático acumulado (SSE). 37 Pred1 Pred2 HR Eabs1 Eabs2 Eqabs1 Eqabs2 83.5838 83.7050 81.5119 2.0718 2.1930 4.293 4.809 71.2830 71.7375 69.2475 2.0355 2.4900 8.436 11.010 68.2208 68.6967 55.1371 13.083 13.5596 179.620 194.873 63.6744 63.9427 64.2253 0.5509 0.2827 179.923 194.953 59.6125 59.8238 43.0674 16.545 16.7563 453.662 475.729 43.4699 43.6820 36.4940 6.9759 7.1880 502.325 527.396 36.0979 36.2709 36.2974 0.1995 0.0265 502.365 527.396 31.9494 32.0742 36.9819 5.0326 4.9077 527.692 551.482 41.0957 41.1995 47.4877 6.3920 6.2881 568.550 591.023 59.8716 59.8306 61.6274 1.7558 1.7968 571.633 594.251 78.8658 78.8716 73.8197 5.0461 5.0519 597.096 619.773 85.6801 85.6515 82.0345 3.6456 3.6170 610.386 632.856 Tabla 4.4: Comparación de modelos ARIMA En la tabla anterior se pueden ver las predicciones junto al SSE, como se ha dicho antes. Fijándose en los valores que toma el error, se puede ver un claro aumento en la tercera predicción de ambos modelos, obviando esa observación y llegando a la última predicción, los errores han ido aumentando poco a poco y han terminado con una diferencia de veinte puntos. Aunque la diferencia no sea alarmante, el modelo que menor error comete es el 2, ARIMA(1,0,1)(2,1,0)12. Además es el que posee valores más bajos de AIC y SBC, como se ha visto antes. Todo los resultados anteriores apuntan a que el modelo 2 es el más acertado pero, si se presta atención a las correlaciones de los componentes del modelo, entre la parte de media móvil regular y la parte autorregresiva regular existe un correlación de 0.913, por lo que en ese aspecto podría ser más acertado el modelo 1, ARIMA(1,0,0)(2,1,0)12. 38 Capítulo 5 Serie Radiación Solar Incidente La última serie a analizar, siguiendo el esquema propuesto, es la radiación solar incidente sobre la superficie terrestre, medida en decenas de kJ/m2. El dato de partida es la radiación solar que llega a una superficie de 1 m2 acumulada durante un día, en el punto donde está colocada la estación. A partir de estos datos, se van a analizar los promedios mensuales de la misma además de la serie mencionada de datos de radiación solar acumulada diaria. 5.1. Serie Radiación Solar Incidente Acumulada Diaria En esta sección se va a analizar de igual forma que en otras ocasiones la serie media diaria radiación solar incidente. En este caso, la serie es algo más larga que las dos anteriores, abarcando seis años más (de 1990 a 2017). De manera semejante a la humedad, al haber detectado errores en los primeros años de datos, se ha optado por eliminar los datos comprendidos entre 1987 y 1990. A partir de ese año, los datos erróneos puntuales se sustituyeron por el promedio de los laterales. Análisis Descriptivo Esta variable oscila, en el periodo que ha sido observada, entre 0 kJ/m2 y 14.23 kJ/m2, y su media es 6.882 kJ/m2. Comenzando la descripción de esta serie, se presenta la siguiente tabla con estadísticos anuales. 39 Año Mínimo Media Mediana Máximo Rango 1990 0.00 6.600 6.710 12.32 12.32 1991 0.00 6.974 7.090 13.05 13.05 1992 0.39 6.725 6.805 12.30 11.91 1993 0.38 6.509 6.340 12.67 12.29 1994 0.65 6.776 6.760 12.76 12.11 1995 0.52 7.032 7.280 12.78 12.26 1996 0.24 6.959 7.000 14.23 13.99 1997 0.16 6.546 6.740 12.83 12.67 1998 0.20 6.693 6.490 12.64 12.44 1999 0.43 6.764 6.380 12.49 12.06 2000 0.24 6.675 6.580 12.61 12.37 2001 0.64 6.780 6.320 12.69 12.05 2002 0.53 6.732 6.390 12.68 12.15 2003 0.55 6.738 6.750 12.50 11.95 2004 0.86 6.730 6.500 12.76 11.90 2005 0.66 6.994 7.060 12.70 12.04 2006 0.78 6.628 6.400 12.74 11.96 2007 0.90 6.602 6.560 12.79 11.89 2008 0.69 6.801 6.805 12.67 11.98 2009 0.91 7.017 7.240 12.69 11.78 2010 0.87 6.943 6.860 13.23 12.36 2011 0.84 7.243 7.250 12.68 11.84 2012 0.91 7.311 6.995 13.08 12.17 2013 0.95 7.236 7.210 13.15 12.20 2014 0.97 7.145 7.080 12.91 11.94 2015 0.81 7.236 7.250 12.84 12.03 2016 0.61 7.063 7.115 12.86 12.25 2017 0.88 7.237 7.180 12.64 11.76 Tabla 5.1: Estadísticos resumen de la serie diaria Radiación Observando la tabla resumen llama la atención que, de entre las series diarias, la radiación solar tiene las medias y las medianas anuales más similares entre sí. Comentando los rangos, destacar la similitud en los valores, consecuencia de la similitud entre los valores máximos y los valores mínimos. A lo largo de los casi 30 años de los que se disponen datos las mediciones de radiación no han variado significativamente. 40 Figura 5.7: Descomposición de la serie Radiación media mensual Análisis de Tendencia En el siguiente gráfico se representa la serie mensual radiación. Figura 5.8: Gráfico de la serie Radiación media mensual A la vista del gráfico de la serie se puede ver una sutil tendencia creciente en la serie, fijando la atención en la tendecia propuesta en la descomposición se puede ver más claramente esa tendencia creciente. Realizando la variante del test de tendencia de Mann-Kendall se 47 obtiene un p-valor de 7.2915e-08 lo que indica que el contraste sí detecta esa tendencia creciente en la serie. Estacionaridad y Modelo ARIMA propuesto Con el propósito de estudiar la estacionaridad de la serie se va a representar la función de autocorrelación. Figura 5.9: ACF de la serie Radiación media mensual En el ACF de la serie se puede distinguir, debido al rápido decrecimiento de los primeros retardos, que la serie es estacionaria por la parte regular. En cambio, por la parte estacional no, debido al decrecimiento lineal que presentan los retardos múltiplos de doce. Las anteriores características de estacionaridad han sido vistas en repetidas ocasiones y es bien sabido que la solución es la diferenciación estacional consiguiendo así que la serie sea estacionaria por ambas partes, observación comprobada mediante el test de Dickey–Fuller. Seguidamente se propone el modelo ARIMA(0,0,1)(0,1,1)12 con constante y ecuación (1 −B12)(Xt−µ) = (1 −θB)(1 −ΘB12)at. Estimando los parámetros Xt=Xt−12 +at+ 0,10426at−1−0,92828at−12 −0,0968at−13 + 0,0192. Puesto que la constante del modelo es significativa y coincide con la media de la serie diferenciada estacionalmente, por ser un modelo de media móvil sin parte autorregresiva, podemos decir que la radiación aumenta 0.0192 unidades al año. Resultado coherente con el obtenido tras el análisis diario de la serie. Pasando a la validación, según los test de Ljung-Box (figura 7.17 del anexo) los residuales son incorrelados. Pero, como en el caso anterior, no son normales aunque su media sí se pueda suponer cero. 48 Capítulo 6 Conclusiones En este trabajo se ha tratado tanto el estudio de tendencia como la modelización y predicción mediante la metodología ARIMA de variables atmosféricas. Debido al carácter climático de los datos, es de vital importancia hacer hincapié en el análisis de tendencia, con el fin de observar la evolución del cambio climático. En estas dos últimas décadas de datos, se ha concluido, después de tener en cuenta diferentes escenarios, que la temperatura en España no ha variado significativamente. Según la AEMET [11] sí que hubo un incremento de las temperaturas en los últimos 60 años, en concreto de 1.3 ºC, más evidente a partir de 2010. Volviendo a los datos de los que se dispone, probablemente, debido al periodo de observación más reducido, ese incremento no ha sido detectado, aunque en los gráficos de las series sí se puede ver ese crecimiento más fuerte a partir de 2010. En cambio, en las otras dos mediciones, se ha observado un cambio. La humedad relativa en los últimos años ha disminuido ligeramente, y la radiación solar ha aumentado en contraposición. Decir que el incremento en esta última se relaciona con una menor presencia de nubes, lo que aumenta la radiación directa del sol. Por último, una posible linea futura consistiría en mejorar los modelos para las series diarias utilizando otro tipo de metodología menos restrictiva que los modelos ARIMA y establecer modelos que permitan relacionar de un modo estocástico la temperatura y la radiación solar. 49 Capítulo 7 Anexo Serie Temperatura Figura 7.1: Test de Ljung-Box de la serie diaria Temperatura Figura 7.2: Qqplot de la serie diaria Temperatura 50 Figura 7.3: Test de Ljung-Box de la serie Temperatura media mensual Figura 7.4: Diagrama de cajas múltiple de la serie mensual Temperatura máxima 51 Figura 7.5: ACF de la serie Temperatura máxima mensual Figura 7.6: Test de Ljung-Box de la serie Temperatura máxima mensual 52 Figura 7.7: Diagrama de cajas múltiple de la serie mensual Temperatura mínima Figura 7.8: ACF de la serie Temperatura mínima mensual 53 Figura 7.9: Test de Ljung-Box de la serie Temperatura mínima mensual Serie Humedad Relativa Figura 7.10: ACF de la serie diaria %Humedad Figura 7.11: Test de Ljung-Box de la serie % Humedad diaria 54 Figura 7.12: Qqplot de la serie % Humedad diaria Figura 7.13: ACF de la serie %Humedad media mensual Figura 7.14: Test de Ljung-Box de la serie %HUmedad media mensual modelo 1 55 Figura 7.15: Test de Ljung-Box de la serie %Humedad media mensual modelo 2 Serie Radiación Solar Incidente Figura 7.16: Test de Ljung-Box de la serie Radiación diaria Figura 7.17: Test de Ljung-Box de la serie Radiación media mensual 56