El aumento progresivo de las temperaturas a escala planetaria y la redistribución de las precipitaciones globales que acompañan al cambio climático inducido y, en general, el incremento general de la variabilidad de su comportamiento, constituyen a día de hoy uno de los principales vectores de cambio en el funcionamiento hidrológico de lagos y humedales (Bergkamp y Orlando, 1999). Este tipo de ecosistemas húmedos es especialmente sensible a los cambios en el régimen termo-pluviométrico traduciéndose rápidamente sus consecuencias en alteraciones de su funcionamiento hidro-geomorfológico (Borja, Camacho y Florín, 2010). El cambio climático inducido, entendido según Duarte (2011) “…como el efecto de la actividad humana sobre la dinámica climática”, constituye un impulsor de cambio directo en el funcionamiento de los ecosistemas húmedos con un alto impacto y una tendencia continua, que irá incrementando su intensidad en las próximas décadas (Burkett y Kusler, 2000; Erwin, 2009; Borja, Camacho y Florín, 2012). En esta línea, los trabajos aplicados a la Península Ibérica (e.g. Ayala-Carcedo, 2002; Sánchez y Míguez-Macho, 2010) pronostican un descenso en las precipitaciones y un aumento de la evapotranspiración como consecuencia del calentamiento global, lo que acarreará importantes alteraciones en el régimen hidrológico de la mayoría de los humedales. Abundando sobre este particular, Álvarez-Cobelas, Catalán y García de Jalón (2005) afirmaron que los cambios en los patrones hídricos serían especialmente visibles en todo tipo de humedales, muchos de los cuales pasarían de ser permanentes a estacionales, y otros incluso podrían llegar a desaparecer (Moreno-Rodríguez, 2005). En el caso de Andalucía las previsiones sobre las consecuencias del cambio climático sobre lagos y humedales no difieren en nada de las realizadas para el conjunto de la Península Ibérica. Incluso podrían ser peores (IPCC, 2007). Las afecciones más significativas se prevén que se produzcan en los humedales de montaña (Sierra Nevada) y en ambientes dependientes de las aguas subterráneas (Álvarez-Cobelas, Catalín y García de Jalón, 2005). En lo que se refiere al ámbito de Doñana, según Fernández y Borja (2006), el rasgo más significativo de los posibles efectos del cambio climático inducido lo constituye la tendencia hacia una progresiva aridificación, puesta de manifiesto en un aumento de las temperaturas y una disminución de las precipitaciones, como asimismo han puesto de manifiesto Vetter, Aguilar y Sousa (2016). El escenario previsto bajo estas circunstancias incluye una tendencia general caracterizada por una caída de los aportes hídricos, una disminución de la humedad del suelo, un aumento del estrés hídrico, una intensificación de los procesos aluviales, una reactivación de procesos dunares, o un abatimiento de los acuíferos, etc. (Borja, Díaz del Olmo y Borja, 2008).

En general, el estudio de las repercusiones del comportamiento del clima sobre los sistemas lacustres y palustres fundamenta sus predicciones en modelizaciones genéricas caracterizadas, principalmente, por asumir escenarios de impacto antrópico simulados y estar formuladas para escalas espaciales de poco detalle. De igual modo, este tipo de predicciones no suele tener en consideración la dinámica hidrológica propia de los acuíferos, ni sus posibles afecciones antrópicas por cambios en los patrones de explotación de los mismos. Estas son cuestiones que requieren una escala de análisis mucho más específica, lo que, por el momento, no tiene cabida en los modelos utilizados tradicionalmente.

La consideración de las afecciones a los acuíferos constituye una cuestión primordial a la hora de establecer predicciones a futuro del funcionamiento hidrogeomorfológico de los sistemas acuáticos vinculados directamente con una alimentación hipogénica. Este particular es especialmente significativo en el caso de acuíferos con problemas de sobrexplotación como puede ser el caso de Doñana, donde la afección antrópica es cada vez más evidente e intensa (Trick, 1998; Guardiola-Albert, Murillo, Martín Machuca, y López Geta, 2004; Guardiola-Albert, García-Bravo, Mediavilla y Martín Machuca 2009; Custodio, Manzano y Montes, 2009; Guardiola-Albert y Jackson, 2011). Entre las más significativas se encuentra la sobreexplotación del acuífero Almonte-Marismas, principal causa del aumento de la estacionalidad de los humedales con alimentación hipogénica, o directamente, de su pérdida. Las extracciones en el entorno de Matalascañas, que se disparan en los meses de verano, la proliferación de pozos ilegales, el crecimiento de cultivos bajo plásticos, etc. (Rebollo, Mediavilla, Ruiz y Díaz, 2008), han hecho que los niveles del acuífero sean cada vez más bajos.

Los modelos de evaluación del comportamiento futuro de niveles piezométricos, basados en la contrastación de series de datos de campo con los modelos climáticos generales existentes, pueden proporcionar información ajustada, según los casos concretos, sobre la evolución a corto y medio plazo de los niveles freáticos y su relación con la alimentación de humedales. El presente trabajo avanza en esta línea, constituyendo su principal objetivo el desarrollo de una propuesta metodológica que permita la modelización de diferentes escenarios a futuro del funcionamiento de los niveles piezométricos, aplicado al caso concreto de la laguna del Charco del Toro, humedal integrante del complejo palustre de Doñana. Para ello se plantean, como objetivos específicos: caracterizar las relaciones entre las variables climáticas y el nivel piezométrico de la laguna; desarrollar modelos matemáticos que permitan interpretar el comportamiento intraanual de su nivel piezométrico a partir de variables ambientales (fundamentalmente climáticas); y elaborar modelos que permitan predecir, para el horizonte de finales del siglo XXI, el comportamiento de sus niveles piezométricos.

La laguna de Charco del Toro se localiza al suroeste del Parque Nacional de Doñana, en la zona conocida tradicionalmente como los cotos o las dunas estabilizadas (Figura 1). Ubicada aproximadamente a un 1,5 km al este de la población turística de Matalascañas (Almonte, Huelva), Charco del Toro ocupa la zona de contacto entre tres de las cinco unidades morfosedimentarias que conforman el Manto Eólico Litoral de El Abalario-Doñana (MELAD) (Borja y Díaz del Olmo, 1994 y 1996; Borja, et al., 1999; Montes, Borja, Bravo y Moreira, 1998; Zazo et al., 1999). Emplazada concretamente sobre el Alto Manto Eólico Húmedo (AMEh), su extremo meridional queda limitada por el Manto Eólico de Dunas Semiactivas (MEDS), mientras que en el sector occidental hace lo propio el Alto Manto Eólico Seco (AMEs) (Figura 2). Desde el punto de vista genético, su origen se asocia al desarrollo del sistema de dunas parabólicas correspondiente al episodio morfosedimentario del MEDS que, dinamizado por vientos de dirección WSW predominante, delimita la laguna en su sector SE. Las dataciones OSL practicadas sobre los tres mantos eólicos anteriormente mencionados han permitido establecer una primera crono-secuencia con fechas que confieren unas edades situadas en el tránsito entre el Holoceno medio y el superior para el sistema dunar que acoge la cubeta lagunar (AMEh), y plenamente históricas para los sistemas de dunas parabólicas del MEDS (siglos VIII a X dC) (Borja, Díaz del Olmo, Borja, Recio, y Lama, 2014).

El funcionamiento hidrogeomorfológico de la laguna de Charco del Toro es de carácter predominantemente hipogénico (Borja, 2011), constatándose que el aporte de aguas procedentes del acuífero regional constituye la principal entrada de agua al humedal durante un ciclo hidrológico referido a la media regional de las precipitaciones. El análisis de los datos proporcionados por el piezómetro profundo ubicado al suroeste del humedal instalado por el IGME (Ref. SGOP49-S1 y nº de inventario 1143310036), y su cotejo con un levantamiento topográfico de detalle llevado a cabo en el ámbito lagunar, indican que, en efecto, es el flujo subterráneo el responsable de la aparición y el mantenimiento de la lámina de agua en el vaso lagunar durante la fase húmeda del ciclo hidrológico (Borja et al., 2008). No obstante, para que este mecanismo entre en funcionamiento es necesaria la acumulación previa de entre 200 y 300 mm de precipitación, dependiendo de si el ciclo hidrológico anterior ha sido húmedo o seco. En cualquier caso, una vez alcanzado este umbral, la respuesta del acuífero a cada evento pluviométrico es prácticamente inmediata (de horas a días) transmitiéndose rápidamente sus efectos en forma de acumulación de agua en la cubeta lagunar. Cabe destacar, sin embargo, que en la fase inicial de acumulación de precipitaciones, es posible reconocer la presencia de flujos subsuperficiales condicionados por la existencia de horizontes edáficos de carácter hidromorfo que, actuando como niveles de alta impermeabilidad dirigen, en ciertos sectores de la laguna, la circulación de los flujos de forma lateral. Después de que esto ocurre, el agua procedente del acuífero inunda la cubeta y anula cualquier tipo de efecto relacionado con los flujos subsuperficiales, erigiéndose en el principal protagonista del funcionamiento hídrico de la misma durante el resto del ciclo hidrológico.

En la actualidad, el modelo de funcionamiento hidrogeomorfológico para un año en la media pluviométrica se aleja notablemente de lo que parece haber sido la dinámica habitual de este humedal. El análisis de la bibliografía disponible muestra que la laguna de Charco del Toro debería exhibir un hidroperiodo permanente fluctuante con una duración de la inundación bastante mayor, tanto en permanencia como en extensión de la lámina de agua, de la que actualmente presenta. En este sentido, desde los años 90 viene afirmándose que esta alteración en su hidroperiodo puede venir motivada por el abatimiento experimentado por el acuífero durante las últimas décadas, debido al aumento de las extracciones de agua subterránea para el abastecimiento del cercano núcleo turístico de Matalascañas (Vela, Rodríguez, y Tenajas, 1991; Custodio y Palancar, 1995; Serrano y Serrano, 1996; Muñoz Reinoso, 2001; Coleto, 2003; Lozano, 2007; Serrano y Zunzunegui, 2008; Custodio et al., 2009). Ciertamente, el estudio detallado de los niveles piezométricos en el entorno lagunar pone en evidencia dicho abatimiento del acuífero (Rebollo et al., 2008), comprobándose igualmente que la cubeta media y alta de la laguna registran un progreso notable del pinar (Pinus pinea) y el jaguarzal (Hallimium hallimifolium), así como un menor acaecimiento de la saturación de los perfiles edáficos, lo que denota una menor frecuencia de la inundación en estos sectores. Todo apunta, así pues, a que estamos asistiendo a una importante alteración en el funcionamiento hidrogeomorfológico del humedal.

Los datos utilizados en el análisis de la dinámica hídrica y la valoración de las fluctuaciones del nivel freático del entorno del Charco del Toro proceden del piezómetro profundo localizado en el sector suroccidental de la laguna. Se trata de una instalación realizada por el IGME (SGOP49 S1), de 14,50 m de profundidad, con la rejilla colocada entre los 11,40 y 14,20 m, en cuyo interior se ha colocado un medidor automático que proporciona información cada 6 horas (Figura 3). Para nuestro estudio se ha contado con la serie de datos correspondiente a los ciclos hidrológicos 2001-02 al 2004-05. Por su parte, la Estación Biológica de Doñana ha proporcionado la información meteorológica con registros mensuales de las variables (1 a 11 de la Tabla 1), sobre temperatura, radiación, precipitación y evapotranspiración. El intervalo temporal que abarcan dichas series es de septiembre de 2001 a diciembre de 2004.

Por último, se ha empleado la serie temporal de niveles del mar medidos en el mareógrafo de Huelva para el mismo intervalo temporal (variable 12 en la Tabla 1), en formato RLR (Fraile, 2011; Fraile y Fernández 2016), proporcionada por el Permanent Service for Mean Sea Level, a partir de los registros gestionados por la red de mareógrafos REDMAR de Puertos del Estado.

Adicionalmente y con la intención de relacionar los niveles piezométricos de la laguna no sólo con los registros de cada una de las variables independientes de la tabla 1, sino con los totales acumulados o promediados en los n meses anteriores, se han generado 44 variables complementarias a partir del cálculo del promedio de los n meses anteriores al del registro, para las variables de temperatura media, radiación media, precipitación mensual y evapotranspiración media (Tabla 2), sumando un total de 55 variables independientes a emplear en el análisis del comportamiento de la variable dependiente nivel piezométrico.

Puesto que uno de los objetivos del trabajo consiste en elaborar modelos matemáticos que permitan predecir el comportamiento futuro de los niveles piezométricos de la laguna de Charco del Toro ante el cambio climático a finales del siglo XXI, ha sido también necesario utilizar las expectativas previstas para las diferentes variables analizadas. Debido al consenso existente acerca de que el cambio climático es un fenómeno de una contrastada variabilidad espacial (IPCC, 2007), se han utilizado modelos climáticos que proporcionan salidas de carácter regional a partir de diferentes escenarios (entre ellos algunos de los empleados por el IPCC):

La obtención de los datos de temperatura y precipitación regionales para el horizonte 2070-2099 se realizó desde el portal de la Junta de Andalucía Cambia Clima SIG-Monitor de Cambio Climático (Álvarez et al., 2014). Los escenarios para los que se obtuvieron resultados locales para estos modelos fueron los escenarios sra1b, sra2, srb1, elaborados por el IPCC para su informe de 2007. Los dos primeros de estos escenarios están elaborados a partir de la consideración de que las emisiones de gases de efecto invernadero continuarán incrementándose durante el siglo XXI a un ritmo superior al del siglo XX en numerosos países, mientras que el escenario srb1 presenta unas características más moderadas (Hollweg et al., 2008).

Desde el punto de vista metodológico, la integración y el análisis de las distintas series de datos antes indicadas se ha desarrollado en dos fases. La primera de ellas consistió en la elaboración de un modelo explicativo del comportamiento del nivel piezométrico en la laguna, a partir de todas las variables analizadas, así como las derivadas de éstas. En una segunda fase se elaboró un modelo predictivo a partir del cual interpretar el comportamiento intraanual futuro del nivel piezométrico para finales del siglo XXI en un año promedio. Aunque las características del modelo matemático obtenido en la primera fase implican que éste debe ser capaz de predecir el comportamiento futuro de la laguna, fue imperativo rediseñar el proceso de modelado a partir de las únicas variables climático–ambientales que estén espacializadas de acuerdo a escenarios de cambio climático futuro, puesto que de las 54 variables independientes utilizadas para la elaboración de un modelo predictivo, sólo se pudieron emplear las 24 correspondientes a precipitación y temperatura para la elaboración del mismo.

En lo que se refiere a la primera fase de elaboración de modelo explicativo se realizó un análisis exploratorio de los datos disponibles a partir de correlaciones entre todas las variables definidas y la variable nivel piezométrico. En la Tabla 3 se muestran los coeficientes de correlación de Pearson obtenidos para los 32 parámetros en las que se identificó una correlación estadísticamente significativa. Este análisis permitió descartar la inclusión de la variable nivel del mar, incluida en los análisis inicialmente por su importante incidencia en otros procesos costeros (Ojeda, Álvarez, Martín y Fraile, 2009; Fraile, Álvarez y Ojeda, 2018) en los análisis subsiguientes, puesto que mostró una correlación inversa con el nivel piezométrico: a mayores niveles del mar se correspondieron niveles piezométricos más profundos. Aunque la correlación se identificó como estadísticamente significativa, es evidente que se trata de un fenómeno de solapamiento temporal en los ciclos anuales de ascenso/descenso del nivel del mar y del nivel piezométrico, sin que existan evidencias físicas que den soporte a dicha relación inversa. En caso de haber existido una correlación positiva, habría sido necesario profundizar más en la naturaleza de dicha relación, para determinar si realmente existiese un fenómeno hidrológico de conexión entre nivel del mar y nivel piezométrico.

Seguidamente se realizó un análisis de regresión lineal múltiple mediante el software estadístico SPSS 23.0 y el matemático R, utilizando como variable dependiente el nivel piezométrico de la laguna, y como posibles variables independientes el resto de las consideradas, tanto las originales como las derivadas, con el propósito de elaborar diferentes modelos que permitieran recrear el comportamiento de los niveles de la laguna durante el periodo analizado. Los métodos elegidos para el desarrollo de los modelos fueron los más comúnmente utilizados en la literatura especializada, esto es, el método de introducción de variables por pasos, stepwise (Bendel y Afifi, 1977), el de la eliminación de variables, los de introducción de variables hacia delante (forward) y hacia atrás (backward) (Whittingham, Stephens, Bradbury, y Freckleton (2006), todos ellos elaborados mediante el software SPSS 23.0, así como modelos logarítmicos múltiples (Aneuryn-Evans y Deaton, 1980) mediante el software R.

Puesto que los coeficientes de determinación obtenidos fueron semejantes entre todos los modelos generados por diferentes métodos, se eligió el de introducción de variables por pasos, por proporcionar diversos modelos de complejidad creciente para cada una de las distintas iteraciones (con una variable, con dos, con tres…), al tiempo que posibilitan una interpretación sencilla del modelo. Al incorporar las variables paso a paso, esta aproximación permite distinguir aquellas de mayor peso en el comportamiento de la variable dependiente, en este caso el nivel piezométrico de la laguna objeto de estudio. En cualquier caso, se corroboró que la relación entre la varianza explicada y el número de variables incluidas en los modelos no era inferior en los casos elaborados mediante este método que frente a los demás métodos existentes. Las sucesivas iteraciones de los modelos propuestos por el método de introducción por pasos aparecen recogidas en la Tabla 4.

En una segunda fase del trabajo, se elaboró un modelo predictivo sobre el comportamiento del nivel piezométrico de la laguna a partir del análisis multivariante tomando como variables independientes tanto las originales como las derivadas de precipitación y como dependiente el nivel piezométrico. El modelo resultante (Tabla 5) se sumó al obtenido en la primera iteración de la anterior fase, en la que se obtuvo un modelo explicativo compuesto por una única variable (precipitación acumulada en los 6 meses anteriores, Tabla 3). El segundo modelo requirió un proceso de reanálisis con variables exclusivamente desarrolladas a partir de los datos de precipitación y temperatura (Tabla 1 de nuevo). Se repitió el mismo tipo de análisis de regresión múltiple tomando como variable dependiente el nivel piezométrico de la laguna, y como independientes las 15 basadas en las medias de los n meses anteriores. De modo semejante al paso anterior, se obtuvo un modelo de regresión múltiple basado, en este caso, en variables exclusivamente calculadas a partir de temperatura y precipitación, y se realizó un análisis de regresión múltiple hacia atrás (backward). Tras 10 iteraciones se obtuvo el modelo de la Tabla 5, basado en 5 variables: precipitación mensual, precipitación acumulada en los 12 meses anteriores, temperatura mensual, temperatura media de los 4 meses anteriores al de estudio y temperatura media de los 8 meses anteriores al de estudio.

Finalmente, para poder estimar el comportamiento futuro intraanual de la laguna se introdujeron los valores calculados para cada una de las variables de los modelos en las predicciones de cambio climático regionalizadas en diferentes escenarios de cambio climático. Los modelos utilizados fueron los bcn2, cncm3, echam5 y egmam, y las salidas se obtuvieron para la celdilla correspondiente al área del Charco del Toro a partir de la aplicación web Cambia Clima SIG – Monitor de Cambio Climático (Álvarez, Camarillo, Limones y Pita, 2014). Se eligieron estos modelos debido a que desde la aplicación empleada permiten obtener predicciones locales sobre temperatura y precipitación, evitando de este modo el empleo de datos globales, que especialmente en el caso de la precipitación pueden resultar poco evidentes en cuanto a su distribución espacial en el futuro.

Del anterior desarrollo metodológico se obtuvieron resultados para los dos modelos seleccionados (denominados simple y múltiple), a partir de las expectativas de precipitación identificadas por los modelos climáticos bcn2 (Figura 4), cncm3 (Figura 5), echam5 (Figura 6) y egmam (Figura 7), y de los escenarios de emisión de gases de efecto invernadero sra1b, sra2 y srb1. De este modo se consiguieron 24 simulaciones del futuro comportamiento intraanual de la laguna Charco del Toro para el período 2079-2100.

Es común a todas las simulaciones realizadas una diferencia de 1 - 1,5 m del nivel piezométrico entre el modelo de regresión lineal simple basado exclusivamente en el modelado de la precipitación de los 6 meses anteriores, y el de regresión múltiple basado en variables de precipitación y temperatura. Todas las predicciones que indican un descenso del nivel piezométrico de acuerdo a este último modelo son independientes de los datos de entrada (modelo climático y/o escenario). Dentro de cada uno de los dos tipos de modelos climáticos, los resultados obtenidos evidencian que en el caso de los bcn2 y echam5 los niveles piezométricos serán significativamente inferiores independientemente del escenario de cambio climático elegido, identificándose en los resultados obtenidos para ambos que, a finales del siglo XXI, el nivel promedio de encharcamiento de la laguna se encontrará por debajo del nivel actual, independientemente el escenario considerado. En los resultados obtenidos mediante la aplicación de los datos del modelo echam, dependiendo de los escenarios elegidos, pueden llegar a identificarse meses en los que el nivel piezométrico del Charco del Toro estará ligeramente por encima del actual, al final del verano y principios del otoño.

En el caso de los resultados obtenidos a partir de los datos del cncm3, se identifica un retraso en el proceso de recarga del nivel piezométrico de la laguna en todos los escenarios climáticos, con un otoño más seco y consecuentemente con niveles más bajos que en la actualidad, un invierno con niveles piezométricos inferiores a los actuales, manteniéndose las diferencias observadas durante el otoño entre el nivel medio actual y el nivel medio predicho por los modelos elaborados, y una primavera en la que el proceso de progresivo abatimiento de nivel piezométrico, que va a asociado al incremento de la ETP durante el transcurso de la primavera, resulta mucho más lento en los modelos elaborados que en los registros analizados.

Dentro de un predominio de resultados en los que el nivel piezométrico de la laguna estimado para finales del siglo XXI se encontrará por debajo del nivel medio actual, el escenario srb1 muestra para este mismo periodo una mayor semejanza a los registros presentes analizados. El escenario sra2 muestra valores muy inferiores a los actuales en los escenarios bcn2 y echam, mientras que el escenario sra1b muestra un comportamiento mucho más variable, mostrando las predicciones más parecidas a los valores presentes en el escenario de egman.

En lo que atañe al proceso metodológico, los resultados obtenidos evidencian una extrema dependencia del modelo estadístico elaborado, frente a una menor importancia de modelos climáticos o escenarios. En este sentido, resulta perentorio elegir entre elaborar predicciones a partir de modelos sencillos, más fáciles de interpretar, o predecir el comportamiento futuro del nivel piezométrico mediante modelos de regresión múltiple, de mayor complejidad pero con una mayor probabilidad de implicar errores, así como de mayor dificultad para su interpretación. No obstante, aunque el modelo de regresión múltiple generado incluye 5 variables, estas se basan únicamente en dos variables principales, precipitación y temperatura. Así pues, según afirman Barreira, Aguilera, y Gómez (2012), uno de los principales puntos débiles de la investigación centrada en evaluar el futuro de determinadas variables dependientes de modelos y/o escenarios radica en la imposibilidad de contrastar los resultados obtenidos con registros reales, tal y como se ha contrastado, por ejemplo, en trabajos sobre el futuro de la erosión costera en determinados ámbitos (Fenster, Dolan y Elder, 1993), inundaciones futuras por subida del nivel del mar (Fraile-Jurado y Ojeda-Zújar, 2013; Fraile-Jurado et al. 2017) o, dentro del campo de estudio de este trabajo, en ambientes dependientes de aguas subterráneas (Lerner, 1986; Kurtz, Hendricks Franssen, Kaiser y Vereecken, 2014). Aunque el procedimiento común para este tipo de aproximaciones es el calibrado del modelo siguiendo procesos estandarizados (Nash, Eamonn y Sutcliffe 1970) basados en la valoración de los residuos de predicciones realizadas para algunos de los registros disponibles, la corta duración de la serie empleada ha dificultado esta aproximación. No obstante, y considerando que los resultados de este trabajo se derivan de la aplicación de dos modelos de cambio del nivel piezométrico basados en procesos de modelado de complejidad creciente, cabe destacar que es posible identificar en la bibliografía existente trabajos que han modelado elementos semejantes a partir de una única variable (Bouraoui, Vachaud, Li, Le Treut, y Chen, 1999) o de varias (Dudula y Randhir, 2016). Aunque existe un claro consenso acerca de que un menor número de variables es deseable para cualquier proceso de modelado (Rodríguez, 2005), es evidente que los modelos propuestos en este trabajo, con una y cinco variables (derivando estas cinco únicamente de dos variables originales) presentan un número aceptable de ellas. El modelo predictivo emplea únicamente variables para las cuales existen salidas a escalas de escaso detalle como las necesarias para este trabajo. A pesar de que en el último lustro se han dado pasos hacia el modelado de otras variables que en esta investigación han presentado correlaciones estadísticamente significativas con el nivel piezométrico, como los trabajos sobre modelización del viento de Sydeman et al. (2014), frecuentemente los resultados publicados no permiten una homogeneización escalar con respecto a las salidas de los modelos de cambio en precipitación y temperatura, dada la elevada precisión de estas últimas, o bien aún se encuentran ante incertidumbres de difícil gestión en la elaboración de los escenarios de emisiones, tal y como sucede con el potencial incremento de la radiación solar (Bais et al., 2015). Por su parte, Hosseinzadehtalaei, Tabari y Willems (2016) subrayan los problemas de utilizar como referencia una variable como la ETP, al ser esencialmente un indicador integrado por otras variables, lo que, de cara a pronosticar cambios en el futuro, implica la asunción de amplias incertidumbres independientemente del escenario desarrollado. La existencia de una correlación inversa estadísticamente significativa entre la variable nivel del mar y el nivel piezométrico pone de manifiesto el problema de las falsas correlaciones, es decir, el riesgo que se asume en ocasiones en el análisis temporal al valorar como verdaderas correlaciones estadísticamente significativas que, sin embargo, no implican la existencia de relaciones naturales entre dichos elementos, tal y como sostienen Stevens et al. (2017).

Por su parte, en lo que se refiere a los resultados, la consideración de variables derivadas de las expectativas de incremento térmico asociadas al cambio climático suponen una aportación relevante, puesto que la mayoría de la bibliografía consultada que ha abordado este tipo de aproximaciones lo ha hecho desde el modelado en los cambios esperables en la precipitación, tanto en la evolución de cuencas fluviales (Shrestha, Cochrane, Caruso, Arias y Piman, 2016), como de lagos (Shrestha, Dibike y Prowse, 2012). No obstante, existe un importante precedente en el trabajo de Dudula y Randhir (2016), donde, de cara al modelado de cuencas, se consideran no sólo los cambios en la precipitación, sino también los registrados por las temperaturas. El hecho de que las variables que mostraron una mayor correlación con respecto al comportamiento del nivel piezométrico fueron las derivadas de la precipitación acumulada en períodos de 5, 6 y 7 meses, coincide con las aportaciones hidrogeomorfológicas realizadas por Borja (2011) y Borja, Díaz del Olmo y Borja (2009). No obstante, tanto los elevados valores obtenidos con el coeficiente de correlación de Pearson para las variables derivadas de la radiación incidente, de la ETP y de la temperatura acumulada, como el hecho de que una de las variables incluidas en el modelo explicativo sea el promedio de la ETP durante los 10 meses anteriores al registro analizado, ponen de manifiesto la importancia de variables frecuentemente descartadas en la literatura científica, tal y como se refleja en los trabajos en los que se analizan únicamente variables vinculadas a la precipitación de manera directa (Bouraoui et al., 1999) o indirecta mediante la consideración de diferentes retrasos en el período de recarga (Pinault, Amraoui, y Golaz, 2005; Tomás, Domenech, Mira, Cuenca y Delgado, 2007; Khan, Gabriel y Rana, 2008; Caballero y Ladouche, 2015). Finalmente, los resultados obtenidos en el cálculo de las expectativas futuras del comportamiento del nivel piezométrico en el modelo predictivo elaborado mediante regresión lineal múltiple evidencian el riesgo al que está expuesto el humedal objeto de estudio. Aunque los resultados obtenidos se han calculado a partir de cuatro modelos climáticos y unos escenarios de emisiones de gases de efecto invernadero concretos, de acuerdo con la práctica totalidad del resto de informes científicos sobre cambio climático en la zona mediterránea, las asunciones realizadas por estos modelos son realistas, tanto en lo que se refiere al descenso de las precipitaciones (Hendrix y Salenyah, 2012), como en términos de incremento de la temperatura (Giorgi y Lionello, 2008). Los resultados obtenidos se adecúan también a las valoraciones realizadas por Klove et al. (2014), por las cuales el impacto del cambio climático en humedales de pequeño tamaño será significativamente más elevado que en los de mayores dimensiones.

El método empleado permite identificar los principales factores presentes en el comportamiento reciente del nivel piezométrico de la laguna de Charco del Toro mediante la elaboración de modelos explicativos de regresión lineal múltiple a partir de diferentes variables ambientales. Los resultados sugieren que la interpretación realizada con anterioridad a este análisis estadístico es esencialmente válida. Según la cual, el funcionamiento hidrogeomorfológico de la laguna de Charco del Toro depende directamente del nivel del acuífero en su entorno inmediato y, éste, a su vez, de la recarga directa por precipitaciones. No obstante, se han identificado otras variables con una fuerte correlación con el nivel piezométrico, como la temperatura, la ETP y la radiación solar. El papel de estas variables, no siempre consideradas en la literatura, no está aun suficientemente contrastado, por lo que serán necesarios análisis más profundos mediante series temporales más largas, tanto en este ámbito natural de estudio como en otros semejantes. Por otra parte, pese a que la laguna estudiada está localizada muy cerca de la costa y a escasa altitud, el análisis estadístico realizado ha permitido descartar la influencia de los cambios del nivel del mar en el nivel piezométrico de la laguna objeto de estudio, al mostrar una correlación inversa con este.

Los modelos de regresión lineal elaborados permiten realizar predicciones sobre el futuro de la laguna mediante el uso de los resultados locales de modelos climáticos a partir de cualquier escenario de cambio climático. No obstante, al estar limitadas las salidas regionales de estos modelos a las variables de precipitación, temperatura, y las que se derivan de éstas, la complejidad de los modelos predictivos obtenidos es menor que la de los modelos explicativos. Las diferencias observadas entre los modelos que únicamente incorporan las variaciones futuras en la precipitación y los que emplean variables no sólo de precipitación sino de temperatura, son de un rango tan amplio que exigirán futuros procesos de análisis de ambos tipos de modelos. En cualquier caso, el futuro del Charco del Toro estará marcado, ante cualquier escenario de cambio climático, por un aumento progresivo de la estacionalidad de su hidroperiodo, sin que se pueda descartar una posible desecación a largo plazo. De acuerdo con los modelos predictivos más complejos, esta desecación será total, mientras que según los modelos predictivos sencillos, los efectos irán encaminados tanto a una mayor estacionalidad como a un retraso de los periodos de ascenso y descenso del nivel piezométrico, y, por tanto, a mantener un cierto embalsamiento de agua en la cubeta lagunar.

Somos conscientes, por último, de que el escaso desarrollo temporal de la serie de datos con la que se ha trabajado puede restar alcance a los resultados. Aun así, este ejercicio de análisis explicativo y, sobre todo, predictivo, posee un alto valor metodológico, ya que constituye un tratamiento novedoso, al aplicar métodos de análisis multivariante a la interpretación actual y futura del funcionamiento de los niveles piezométricos a escala de humedal, lo que puede constituir un punto de partida de cara a la definición de estrategias de adaptación a las consecuencias del cambio climático.