Generalidades

Al parecer, los eventos extremos del tiempo atmosférico son cada vez más severos y frecuentes, con sus respectivos impactos devastadores; además, quizás son originados por el cambio climático. Estas dos aseveraciones requieren que los análisis estadísticos orientados a su predicción y pronóstico estén basados en datos meteorológicos y climáticos de calidad, es decir, que sean precisos y homogéneos. En la realidad, obtener registros o series amplias de datos climáticos perfectamente homogéneos es casi imposible debido al cambio inevitable del área que rodea el sitio de la estación climatológica, el cual afecta los datos observados. La pérdida de homogeneidad ocurre cuando hay cambios en el registro generados por causas no climáticas (^{Guijarro, 2014}; ^{Yozgatligil & Yazici, 2016}).

En resumen, todos los registros o series de datos que se miden en las estaciones climatológicas, con frecuencia sufren alteraciones que no están relacionadas con las variaciones del clima. Tales modificaciones se originan en los cambios de ubicación de la estación, al reemplazar sus instrumentos o cambiar la técnica de medición, y al alterarse las condiciones físicas que rodean al sitio. La homogenización es un procedimiento para detectar y corregir las alteraciones artificiales citadas. Respecto a las técnicas de homogenización, ^{Guijarro (2014)} las describe someramente y se enfoca en los paquetes de cómputo disponibles. Una serie de datos climatológicos se considera homogénea cuando las condiciones de medición en la estación no han variado con el tiempo (^{Beaulieu et al., 2008}). ^{Guijarro (2014)} destaca la importancia de los servicios meteorológicos en evitar por todos los medios cambios que afecten las mediciones climatológicas, pues ello hace innecesario el uso de los métodos de homogenización.

Las erupciones volcánicas, siendo un fenómeno natural que origina una alteración de los registros climáticos, definen puntos de cambio (change points) que no implican pérdida de homogeneidad; por el contrario, los puntos de quiebre (break points) son puntos de cambio causados por alteraciones no climáticas y por lo tanto implican pérdida de homogeneidad (^{Toreti et al., 2011}).

La homogenización de las series climáticas es una tarea difícil que debe realizarse con cuidado extremo, en especial cuando los datos adicionales o datos históricos (metadata) de la estación climatológica no están disponibles, lo que es un caso común. El objetivo de la homogenización es eliminar, o al menos reducir, las alteraciones no climáticas mientras se preserva la señal climática (^{Toreti et al., 2011}). Se han propuesto muchos métodos para homogeneizar las series climatológicas; una revisión completa que incluye sugerencias para su uso se puede consultar en ^{Peterson et al. (1998)}, así como en ^{Beaulieu, Ouarda y Seidou (2007)}.

Una primera clasificación de los procedimientos y/o pruebas estadísticas que detectan y corrigen la pérdida de homogeneidad establece dos grupos: el primero se basa en la disponibilidad de datos adicionales sobre la evolución histórica de la estación climatológica, conocidos como metadatos, y el segundo en la usencia de tales datos. ^{Rhoades y Salinger (1993)} exponen técnicas para el primer grupo bajo dos enfoques: (1) empleando estaciones circunvecinas; (2) en estaciones aisladas, o de registros antiguos o primigenios. En general, los datos históricos de la estación climatológica son esenciales para validar los cambios detectados por las pruebas y, desafortunadamente, casi nunca están disponibles.

Otra clasificación de las técnicas de homogenización las divide en pruebas absolutas, las cuales emplean en exclusivo la serie o registro bajo análisis, y en pruebas relativas, que emplean registros cercanos o correlacionados para formar la serie de referencia. Ambos enfoques son válidos y útiles, pero también exhiben desventajas. ^{Yozgatligil y Yazici (2016)} comparan siete pruebas de cada tipo y encuentran que los test relativos tienen mejor desempeño; entre ellos, citan en primer lugar la prueba de homogeneidad normal estándar, designada SNHT (Standard Normal Homogeneity Test).

^{Syrakova y Stefanova (2009)} indican que debido a la dificultad para establecer a priori cuáles estaciones son homogéneas, optan por seguir un esquema iterativo de mejoramiento gradual de la homogeneidad, el cual fue propuesto por ^{Hanssen-Bauer y Forland (1994)}, y ha sido aplicado también por ^{Moberg y Alexandersson (1997)}; ^{González-Rouco, Jiménez, Quesada y Valero (2001)}, y ^{Tuomenvirta (2001)}.

En cambio, ^{Wijngaard, Klein-Tank y Können (2003)} aplican y citan dos ventajas de las pruebas absolutas: (1) fácil aplicación en redes poco densas; (2) superiores a las pruebas relativas en redes con cambios simultáneos, ya que éstas no los detectan. ^{Dhorde y Zarenistanak (2013)} utilizan un enfoque mixto en el estudio de homogeneidad de 20 estaciones climatológicas, aplicando las pruebas de Pettitt y SNHT a cada serie, y la prueba SNHT para cambios graduales en la media con serie de referencia.

^{Ducré-Robitaille, Vincent y Boulet (2003)} evaluaron ocho técnicas de homogenización de series de temperatura, empleando series simuladas homogéneas, y conteniendo uno o varios cambios en la media (steps). Encontraron que la mayoría tiene buen desempeño, pero dos de ellas tienden a ser más eficientes: la prueba SNHT y la técnica de la regresión lineal múltiple.

^{Beaulieu et al. (2008)} compararon también ocho pruebas estadísticas para detectar falta de homogeneidad en series de precipitación anual, utilizando varios miles de series sintéticas tanto homogéneas como no homogéneas; encontraron que ningún método es eficiente en todos los tipos de pérdida de homogeneidad, pero tres de ellos tienen un desempeño mucho mejor: el test bivariado, el método de Jaruskova y la prueba SNHT.

En México, y en específico para el estado de Veracruz, ^{Guajardo-Panes, Granados-Ramírez, Sánchez-Cohen, Díaz-Padilla y Barbosa-Moreno (2017)} han aplicado el enfoque absoluto con las pruebas SNHT, Pettitt y Buishand, bajo un esquema de verificación espacial, para seleccionar estaciones climatológicas circunvecinas a la estación que se prueba, dentro de una misma unidad climática, definida con base en el intervalo de variación del elemento del clima que se analiza. Lo anterior equivale a procesar la información climática por subregiones o zonas climáticas.

Concepto de control de calidad

^{González-Rouco et al. (2001)} han definido que los procesos de control de calidad de las series de datos climáticos involucran tres etapas: (1) detección y corrección de valores dispersos (outliers); (2) técnicas de interpolación de datos faltantes, y (3) homogenización de la serie. En realidad, la estimación de datos faltantes ocurre de manera mensual, antes de la homogenización y posterior a ésta, cuando se intenta completar series para tener un periodo común de información. Debido a que todas las series de datos climáticos son extremadamente sensibles a la presencia de valores erróneos y datos dispersos, ^{Eischeid, Baker, Karl y Diaz (1995)} abordan la detección y corrección de datos dispersos según un enfoque objetivo de verificación temporal en el registro analizado y espacial con base en seis técnicas de interpolación de datos mensuales.

Objetivo

Se ha formulado un enfoque práctico o estrategia de verificación de la homogeneidad en registros de precipitación anual, el cual tiene las siguientes características relevantes: (1) se trabaja por regiones o subregiones geográficas, utilizando datos mensuales; (2) se estiman datos mensuales faltantes con una técnica estadística; (3) se detalla cómo se verifican y corrigen los valores mensuales máximos extremos, por truncamiento; (4) se deducen datos anuales faltantes, con una técnica de ponderado regional; (5) se utiliza un enfoque de contraste absoluto, empleando cuatro pruebas estadísticas (Von Neumann, Pettitt, Buishand y SNHT), y (6) se revisa gráficamente cada registro que no fue homogéneo o confiable, y se corrige la pérdida de homogeneidad cuando existe un punto de quiebre. En este estudio se trabajó en la región climática semiárida del Altiplano Potosino de México, procesando 16 registros de periodo amplio, cuyo número de datos anuales varió de 42 a 53.

Corrección de valores dispersos

Los valores dispersos (outliers) son observaciones que distan mucho de la tendencia mostrada por los datos, los cuales corresponden a mediciones erróneas o a eventos meteorológicos extremos. Existen diversos enfoques para contrastar y verificar la variabilidad temporal y/o espacial de la variable climática estudiada, con objeto de identificar los valores dispersos y diagnosticar datos erróneos o valores factibles de ocurrir. Cuando los valores dispersos son sin duda observaciones erróneas, tales datos se eliminan, pero entonces se tiene un problema de datos faltantes. Si son valores factibles, se debe decidir si se corrigen y cómo conviene hacerlo (^{Eischeid et al., 1995}; ^{González-Rouco et al., 2001}).

Debido a las dificultades prácticas que implica la verificación de valores dispersos ―en los registros originales o de campo―, se adoptó un enfoque de reducción de la dispersión en la cola derecha de la distribución de probabilidades, para buscar mejor desempeño de las técnicas homogenización que no son resistentes (^{Lanzante, 1996}). Tal enfoque consiste en acotar los valores dispersos a una magnitud umbral extrema (P _vd ), que es (^{Eischeid et al., 1995}):

donde P ₅₀ es la mediana o el percentil de 50%; P ₇₅ y P ₂₅ son los percentiles de 75 y 25%; y FM es el factor de multiplicación del rango del intercuartil (RIC), con un valor empírico adoptado de 2.75 para los datos de temperatura media y de 4.00 para la precipitación, ambas mensuales. ^{González-Rouco et al. (2001)} utilizan una expresión similar para sus datos de precipitación mensual:

Los valores dispersos superiores a P _vd se sustituyen por tal límite. Este enfoque de truncamiento de los valores dispersos reduce el sesgo que éstos originan y también mantiene información relativa a los valores extremos. ^{Eischeid et al. (1995)} indican que el RIC es utilizado en los procesos de control de calidad debido a que es resistente (^{Lanzante, 1996}) a los valores dispersos. En este estudio se adoptó el criterio definido por la Ecuación (1), por tener mayor justificación teórica y empírica.

Pruebas del enfoque absoluto

Los test que se aplican a la serie estudiada son los del planteamiento de ^{Wijngaard et al. (2003)}, que utiliza cuatro métodos para probar la pérdida de homogeneidad: el test de Von Neumann, el test del rango de Pettitt, el test del rango de Buishand y la prueba SNHT. Estas cuatro pruebas o test suponen que bajo la hipótesis nula, los valores anuales de la serie (Y _i ) de la variable Y que se estudia son independientes y están idénticamente distribuidos. Bajo la hipótesis alternativa ―el primer test―, se asume que la serie no está distribuida de forma aleatoria, y los tres restantes consideran que existe un cambio en la media o punto de quiebre. Lógicamente, la homogeneidad o consistencia de la serie implica que sus datos pertenecen a una sola población que tiene una media invariante con el tiempo (^{Machiwal & Jha, 2012}).

Las pruebas de Pettitt, Buishand y SNHT permiten localizar el año donde es probable ocurriera el punto de quiebre y, en general, son diferentes en su planteamiento o concepción. Se ha encontrado que las dos primeras son más sensitivas a los puntos de quiebre que suceden hacia la mitad de la serie; por el contrario, la prueba SNHT detecta cambios con más facilidad hacia el inicio y final de serie.

Los test de Von Neumann, Buishand y SNHT asumen que los datos Y _i tienen distribución normal; en cambio, la prueba de Pettitt no requiere tal consideración, pues se basa en los rangos y no en los valores de la serie, lo cual también la vuelve menos sensible a los valores dispersos (outliers). La prueba de Von Neumann se considera complementaria a las otras tres, debido a que detecta pérdida de homogeneidad por causas diferentes al cambio en la media (^{Wijngaard et al., 2003}).

Corrección de la pérdida de homogeneidad

Cuando un punto de quiebre se ha identificado en un registro de precipitación anual, los datos previos a tal cambio se corrigen al multiplicarlos por el siguiente factor (^{González-Rouco et al., 2001}):

siendo P-d y P-a las medias de los valores anuales de precipitación después y antes del cambio. ^{Syrakova y Stefanova (2009)} exponen cómo se realiza la corrección de registros de temperatura anual y estacional que se contrastan con pruebas relativas.

Clasificación de las series probadas

La clasificación depende del número de pruebas que rechazaron la hipótesis nula, definiéndose tres categorías, cuyas designaciones según ^{Wijngaard et al. (2003)}, y ^{Dhorde y Zarenistanak (2013)} son las siguientes: útil, dudosa y sospechosa. De acuerdo con ^{Guajardo-Panes et al. (2017)} son las siguientes: confiable, medianamente confiable y no confiable. En este estudio se adopta la designación siguiente: en la clase 1 o “confiable” se agrupan las series que fueron detectadas homogéneas por las cuatro pruebas, o bien, sólo una prueba detectó pérdida de homogeneidad; la clase 2 o “menos confiable” engloba series con dos pruebas que detectan falta de homogeneidad, y la clase 3 o “no confiable” son series con tres o cuatro pruebas que rechazan la hipótesis nula.

Con las series clase 1 se pueden realizar análisis de tendencia y estudios de variabilidad climática; con las clase 2, los resultados de los análisis citados deben ser tomados con reservas, pero son factibles si no existen otras series confiables disponibles; las series clase 3 no se deben utilizar para las estimaciones mencionadas.

Existe una diferencia notable al establecer la clasificación anterior, ya que ^{Wijngaard et al. (2003)} definen un nivel de significancia (α) en las pruebas de 1%; en cambio, en los otros dos trabajos citados se emplea un 5%. En este estudio se aceptó α = 5%, para coincidir con las recomendaciones de la prueba SNHT (^{Alexandersson, 1986}).

Selección de estaciones climatológicas

Este estudio empleó la información climatológica mensual del archivo en Excel proporcionado por la Dirección Local San Luis Potosí de la Comisión Nacional del Agua (Conagua), la cual está agrupada en las tres zonas geográficas del estado: Altiplano Potosino, Zona Media y Región Huasteca. Los registros de datos procesados corresponden a la variable climática precipitación anual en milímetros, integrado cada dato por la suma de los doce valores de lluvia mensual.

El archivo del Altiplano Potosino abarca más de 60 estaciones climatológicas, pero algunas están repetidas por pertenecer a la extinta Secretaría de Recursos Hidráulicos (SRH) y al Servicio Meteorológico Nacional (SMN). Así, se eliminaron las estaciones con menos de 40 años de registro y con lapsos de años faltantes. Con tales restricciones, se obtuvo una muestra de 16 estaciones, cuyo año inicial de registros se indica con AI en la Tabla 2 y Tabla 3 de datos por procesar. Todos los registros abarcan hasta el año 2016, excepto la estación Los Pilares, que fue suspendida en 2009.

Para establecer el mayor periodo común de datos se seleccionó el año 1964 para su inicio y, con ello, 13 registros procesados tienen 53 o 52 valores; tres, La Presa, El Grito y Los Pilares, abarcan 42, 48 y 45 años, respectivamente. En la estación Moctezuma, los siete años faltantes de 1964 a 1970, se tomaron del Boletín Climatológico No. 3 (^{SARH, 1980}), con datos del inicio de las estaciones al año de 1978. Estos valores se han indicado en el inicio de la Tabla 3, con paréntesis de llave. En el boletín citado se verificaron los inicios de las estaciones La Presa y El Grito.

En la Figura 1 se muestra la localización geográfica de las 16 estaciones climatológicas procesadas del Altiplano Potosino, México, destacando que las cuatro primeras forman un subgrupo norte-extremo; las tres siguientes integran un subgrupo norte; las cinco que continúan definen el subgrupo del centro; las cuatro últimas establecen el subgrupo del sur.

Figura 1 Ubicación geográfica de las 16 estaciones climatológicas procesadas del Altiplano Potosino, México. 

Estimación de datos mensuales faltantes

En este estudio se adoptó como valor mensual faltante de los años incompletos la moda mensual respectiva, estimada con base en el ajuste de la distribución Gamma Mixta de dos parámetros (^{Campos-Aranda, 2005}), bajo tres condiciones específicas: (1) se estima sólo en el semestre lluvioso, es decir, de mayo a octubre; (2) se completa un máximo de tres meses en cada año incompleto; y (3) cuando la moda no existe, porque la distribución es de tipo J invertida, se emplea la precipitación mensual con probabilidad de no ser excedida de 25% (P ₂₅), como el caso de octubre en la estación Charcas (ver Tabla 4). En la Tabla 2 y Tabla 3, los valores anuales completados con el procedimiento anterior se han indicado entre paréntesis circulares.

En la Tabla 4 se exponen resultados de los ajustes mensuales ―sólo dos meses de cada registro incompleto― por limitación de espacio. De la Tabla 4 se deduce que sólo cuatro estaciones pluviométricas tuvieron registros completos: Vanegas, La Maroma, Los Pilares y Los Filtros.

Corrección de valores mensuales máximos

El ajuste de la distribución Gamma Mixta permitió estimar los elementos de la Ecuación (1) para calcular el valor disperso máximo mensual (P _vd ), citado en la Tabla 4. También se indican entre paréntesis circulares los valores máximos que deberán ser truncados al valor de P _vd por haberlo excedido. Los respectivos valores anuales corregidos por el truncamiento descrito se destacan entre asteriscos en la Tabla 2 y Tabla 3. Debido al proceso adoptado de deducción de datos mensuales faltantes, sólo se corrigieron valores máximos extremos en el semestre lluvioso de mayo a octubre, pero se buscaron valores erróneos en el resto de los meses para eliminarlos.

Estimación de datos anuales faltantes

Las estaciones pluviométricas Matehuala, Charcas, Palo Blanco y El Mezquite presentaron datos anuales faltantes; la primera, tres; la última, dos; y el resto, un valor. Al observar en la Tabla 2 y Tabla 3 que los valores anuales disponibles en las estaciones circunvecinas a cada registro incompleto eran bastante similares, se optó por deducir cada valor anual faltante con base en los tres valores de las estaciones más cercanas, según un ponderado con el inverso de la distancia en kilómetros entre cada estación vecina y la incompleta. Este procedimiento fue descrito y aplicado por ^{Eischeid et al. (1995}), ^{Mishra y Singh (2008)}, y por ^{Campos-Aranda (2013)}; es adecuado en áreas de topografía plana y para pocos datos faltantes. En la Tabla 2 y Tabla 3, los siete valores anuales estimados con el enfoque descrito se citan entre paréntesis rectangulares.

Cuando se deduzcan varios datos anuales faltantes en estaciones climatológicas de áreas de montaña se debe usar el método de la relación normalizada (^{Paulhus & Kohler, 1952}), aplicando factores de ponderación definidos con base en los coeficientes de correlación entre los valores anuales de la estación incompleta y las circunvecinas, como los definió ^{Young (1992)} y aplicaron ^{Eischeid et al. (1995)}.

En la Figura 2 se citan los conceptos y técnicas numéricas que incluye la estrategia propuesta para la verificación de la homogeneidad en registros de precipitación anual de una determinada área geográfica.

Figura 2 Diagrama de bloques de los conceptos y técnicas estadísticas, que incluye la estrategia sugerida para probar y alcanzar homogeneidad de los registros climatológicos. 

Parámetros estadísticos de los registros

En los diez últimos renglones de la Tabla 2 y Tabla 3 se muestran los indicadores estadísticos relativos a cada registro completado. Con AI se indica el año de inicio de cada serie de precipitación anual, cuyo número de datos, dentro del periodo común establecido, es n. Después de la mediana (Med) del registro se citan los coeficientes de variación (Cv), asimetría (Cs) y curtosis (Ck), así como el de correlación serial de orden uno (r ₁).

Verificación de la normalidad

Se aplicó el test clásico de ^{Shapiro y Wilk (1965)} a los registros de precipitación anual con menos de 50 datos, y el de ^{Shapiro (1998)} en aquellos con más 50 valores. Se obtuvo que los seis registros de las estaciones Santa Ma. del Refugio, Palo Blanco, Reforma, El Mezquite, El Peaje y Villa de Arriaga no proceden de una distribución normal.

A estos seis registros se les aplicó la prueba de Cociente de Geary (^{Machiwal & Jha, 2012}; ^{Campos-Aranda, 2015}) y con ello se ratificó que los cuatro primeros registros no proceden de la distribución Normal. En la Tabla 2 y Tabla 3 se observa que los coeficientes de asimetría y curtosis de cuatro registros son los que más difieren de los valores normales de cero y tres. Tales registros serán procesados usando el logaritmo natural de sus datos para volverlos normales.

Valores críticos de las pruebas estadísticas

En la Tabla 1 se observa que los valores críticos de los estadísticos N, P _k , R/(ns) y T _k de cada prueba son función del tamaño del registro (n) y como varía de 42 a 53, en la Tabla 5 se exponen sus magnitudes respectivas interpoladas.

Resultados de las pruebas estadísticas

En la Tabla 6 se han concentrado las magnitudes máximas que se obtuvieron en cada prueba para su indicador (N, P _k , R/(ns) y T _k ), cuyo contraste contra los valores de la Tabla 5 permitieron indicar si cada registro resultó homogéneo (H) o no homogéneo (NH). Se obtuvo que 11 registros son confiables, es decir, 68.8% de las series de precipitación anual analizadas; tres son menos confiables (Vanegas, Palo Blanco y Villa de Arriaga), y dos son no confiables (Moctezuma y El Mezquite). Para estas dos últimas estaciones, la prueba de Pettitt aporta sus valores máximos en los 16 registros procesados: 236 y 234.

Análisis de los registros no homogéneos

Con respecto a los registros menos confiables de las estaciones Vanegas, Palo Blanco y Villa de Arriaga, es probable que la prueba de Von Neumann detecte la no aleatoriedad debido a la persistencia que muestran tales series, cuyos coeficientes de correlación serial de orden uno (r ₁) fueron 0.241, 0.347 y 0.489; es decir, los más grandes de la Tabla 2 y Tabla 3, junto con el mayor de todos (r ₁ = 0.587) de la estación El Mezquite. A los registros de las tres estaciones citadas se les aplicó la prueba de Anderson (^{WMO, 1971}; ^{Linsley, Kohler, & Paulhus, 1988}), resultando que sus coeficientes r ₁ son significativos.

En la Figura 3 se muestran los gráficos de los estadísticos de las pruebas de Pettitt, Buishand y SNHT, aplicadas al registro de la estación Moctezuma. Se observa que las dos primeras definen valores máximos extremos positivos y negativos, y que la prueba SNHT define un máximo positivo en k = 45 (año 2008). Por lo tanto, no existe un único punto de quiebre según las tres pruebas; pero dos de ellas lo definen en k = 45. En realidad, en 2008 comienza un periodo de siete años extremadamente secos, cuya presencia fue ratificada en tal zona, para una sequía regional que afectó toda la altiplanicie mexicana. Por lo anterior, no procede la aplicación de la Ecuación (10) en k = 45, para corregir el registro previo.

Figura 3 Gráficos de las pruebas estadísticas indicadas, aplicadas al registro de precipitación anual de la estación Moctezuma del Altiplano Potosino, México. 

Por otra parte, en la Figura 4 se exponen los gráficos de las pruebas de Pettitt, Buishand y SNHT, aplicadas al registro de la estación El Mezquite. Las dos primeras definen valores máximos extremos positivos y negativos, y la tercera un único máximo positivo en k = 37 (año 2000). Aceptando tal año como un punto de quiebre, se define un lapso de 2000 a 2016 con 17 valores, y una media de 579.365 mm para después del cambio, y de 36 datos y una media de 327.811 mm para antes del cambio. Con la Ecuación (10) se obtiene un valor de 1.7674 como factor de corrección (f).

Figura 4 Gráficos de las pruebas estadísticas indicadas, aplicadas al registro de precipitación anual de la estación El Mezquite del Altiplano Potosino, México. 

En registro corregido de la estación El Mezquite tiene los siguientes parámetros estadísticos: mínimo = 59.9 mm; máximo = 1344.1 mm; X- = 579.4 mm; mediana = 567.3 mm; Cv = 0.480; Cs = 0.467; Ck = 3.486, y r ₁ = 0.503. Además, de acuerdo con las pruebas de Shapiro y Geary, resultó que puede proceder de una distribución normal. Los nuevos valores máximos de los estadísticos de las pruebas de Von Neumann, Pettitt, Buishand y SNHT fueron N = 1.016; P _k = 182.0; R/(ns) = 2.03, y T _k = 8.02. Por lo cual resulta no homogéneo según la primera y tercera prueba, pasando a ser un registro menos confiable con la corrección realizada.