Introducción
Los lagos de alta montaña, generalmente, se encuentran en áreas protegidas; sus cuencas presentan gradientes pronunciados, temperaturas frías y radiación solar y ultravioleta de alta incidencia [1]-[3]. Estas características hacen que este tipo de cuerpos de agua presenten una gran sensibilidad a los cambios ambientales, atmosféricos y climáticos asociados a eventos de origen natural o la acción antrópica [4], [5].
Lo expuesto lleva a que las entidades encargadas de la preservación natural realicen seguimientos continuos a estos ecosistemas estratégicos, con el fin de conservar un registro de sus variables físicas, químicas y biológicas [6], [7]. Sin embargo, este tipo de monitoreos sigue enfoques convencionales (muestreos in situ) que tienden a ser limitados en términos de cobertura espacial y frecuencia [8], [9].
Con el fin de mejorar la representatividad y comprender la dinámica espaciotemporal de los lagos de alta montaña, se hace necesario utilizar herramientas como la teledetección. Por medio de esta, pueden obtenerse datos sinópticos que logran abarcar grandes áreas y largos periodos, de forma repetitiva y no intrusiva [10], [11].
En ese contexto, se encuentra que el uso de los productos de teledetección satelital aumentaría sustancialmente la comprensión de la dinámica del cuerpo de agua y los procesos biológicos (producción primaria como la clorofila-a), físicos (transparencia, color y sedimentos en suspensión) y químicos (ciclos de nutrientes) que tienen lugar en la columna de agua y dependen del campo de luz subacuático [9], [12]-[14].
Dentro de los métodos para determinar la calidad del agua, destaca por su sencillez (toma de muestra y lectura de resultado) el disco de Secchi (DS). La profundidad medida con este dispositivo está influenciada por los tres componentes ópticos principales: clorofila-a (Chl-a), materia orgánica disuelta coloreada (CDOM) y el total de materia en suspensión (TMS). Al ser un parámetro cuya lectura depende de la luz, ha sido ampliamente acogido en estudios de teledetección en aguas continentales, además, son necesarios solo unos pocos puntos para obtener modelos basados en la transparencia del agua [14]-[18].
En los últimos treinta años, el uso de imágenes satelitales para el análisis de cuerpos de agua en sus componentes ópticamente activos se ha popularizado. Sensores diversos como Sentinel, Landsat 5, 7, 8, MERIS y MODIS, entre otros, proveen imágenes utilizadas luego para estimar la profundidad medida con el disco de Secchi [18]-[20]. De estos sensores, el de mayor potencial para proveer datos satelitales y correlacionarlos con medidas tomadas en campo de las variables físicas, químicas o biológicas es el sensor MODIS.
Este sensor puede generar productos hasta con 36 bandas espectrales, con una resolución espacial de 250 m y una frecuencia temporal diaria. Por ello, es posible, por medio de estos productos crear series temporales útiles para observar la dinámica y evolución del lago [21]-[26]. Otra gran ventaja de este sensor es que genera productos de reflectancia de superficie (MODU9GA) los cuales no necesitan calibración atmosférica y geométrica, al venir ya sus datos generados en la proyección geográfica sinusoidal [12]-[14].
En ese contexto, numerosos modelos empíricos multibanda o de una banda han sido desarrollados para calcular la transparencia del agua, correlacionando las mediciones de campo de la PDS con las imágenes de sensores remotos, en este sentido, podemos encontrar algoritmos que utilizan diferentes operaciones entre las longitudes de onda localizadas en los espectros verde y rojo y azul, para poder correlacionar este valor, ya sea con las mediciones de PDS o con su logaritmo natural, la representatividad de los modelos revisados en la bibliografía es superior a 0,7 en su R2 [18], [27]-[29].
En lo que respecta al sensor MODIS se encuentra que, además de las diferentes combinaciones que pueden implementarse para las bandas ubicadas en los espectros azul, verde y rojo, es posible utilizar un algoritmo basado en la banda centrada en los 645 nm (banda 1). De ello se encuentra que puede obtenerse algoritmos con estimaciones aceptables que oscilan entre 0,6 y 0,92 en su R2 [12]-[14], [27], [28], [30], [31].
La banda centrada en los 858,5 nm (banda 2) tiene una aplicación más práctica en la estimación del total de los sólidos en suspensión en aguas interiores [32]-[34]. Sin embargo, esta misma longitud de onda puede ser utilizada en la derivación de valores de la PDS, teniendo en cuenta que la turbidez del agua está directamente relacionada con las medidas que pueden obtenerse con el disco de Secchi [33], [35]-[37].
El propósito de este estudio fue probar la practicidad de un modelo basado en imágenes de reflectancia de superficie (RS) generadas por el sensor MODIS para monitorear y caracterizar la PDS en el LG, esta situación es motivada por la necesidad de contar con un registro de monitoreo espaciotemporal de este parámetro dentro del cuerpo de agua. A partir de ello, es posible generar información que sirva para verificar la tendencia y dinámica evolutiva del lago. Para alcanzar nuestro objetivo fue necesario:
1) Validar los algoritmos, utilizando el procedimiento Loocv para encontrar el mejor ajuste en las estimaciones PSD-RS, y evaluar la predicción del modelo utilizando medidas estadísticas como R2 y RMSE.
2) Generar visualizaciones temporales (Horizon-Graphs) que permitan observar la dinámica temporal del cuerpo de agua en 2001-2020.
Materiales y métodos
Lugar de estudio
El lago Guamués (LG) se encuentra a 2765 msnm, al suroeste de San Juan de Pasto, capital del departamento de Nariño (figura 1a-b) [38], [39]. Es un lago altoandino, localizado a lo largo de la falla tectónica de Algeciras, fenómeno geológico que, a su vez, se extiende entre los departamentos del Putumayo, Nariño, Caquetá, Cauca, Huila y Meta [40]. El LG es un importante cuerpo de agua que hace parte de la cuenca alta del río Guamués, su localización geográfica presenta las siguientes coordenadas: 0°50", -1°15" latitud Norte, y 77°05", -77°20" longitud Oeste [41].
Este lago es el mayor humedal de los Andes colombianos, con un espejo de agua que cubre 41,74 km2. La longitud de este lago continental es 14,1 km y, en su parte más ancha, cuenta con 6,2 km; su parte más profunda mide 70 m y la temperatura promedio es de 11 °C (figura 1c) [42].
Trabajo de campo
Se hicieron sesenta mediciones de parámetros fisicoquímicos del agua, producción primaria (Chl-a) y la PDS, las mediciones se adelantaron en seis trabajos de campo (diez estaciones de monitoreo para cada trabajo de campo), que se desarrollaron entre los meses de enero y junio de 2013 (figura 1c).
Transparencia de la columna de agua medida con el disco de Secchi
La medición se llevó a cabo con un disco Secchi (DS), un instrumento bicromático de medición de la penetración luminosa, de 30 cm de diámetro, divididos en cuartos de color negro y blanco alternados y peso de plomo unido a una cuerda con guía por cada metro. El procedimiento consiste en introducir el DS en el agua y dejarlo caer hasta que se pierde de vista. Este proceso se repitió dos veces y se registraron los valores de longitud, promediados luego, para obtener la PDS.
Procesamiento de las imágenes de satélite, generación de imágenes y representaciones gráficas
En la figura 2, se ilustra el desarrollo del trabajo de procesamiento digital de imágenes y generación de resultados. Asimismo, se muestran las etapas desde la selección de las imágenes en términos de amplitud, resolución temporal, adquisición y tratamiento de datos, la generación de la visualización de datos HorizonGraph, así como todos los análisis realizadas en el proceso.
Nota: las líneas punteadas indican los procesos necesarios en la PDI, creación de modelo de estimación de PDS, generación de mapas, series del tiempo y control de calidad de las imágenes. Fuente: elaboración propia.
Descarga de imágenes y preprocesamiento de las imágenes. Una vez definido que el producto MOD-09GA es el producto más adecuado para calcular la PDS [27], [28], [43]-[46], se descargaron las imágenes del sitio web https://earthexplorer.usgs.gov/. Se seleccionaron y descargaron 291 imágenes para crear una colección de datos temporales referentes al periodo enero de 2001 a abril de 2020. Las imágenes fueron procesadas utilizando rutinas computacionales desarrolladas en los lenguajes Python y R, las cuales permitieron realizar la reproyección, recorte espectral y espacial del Tile MODIS H10V8, lo cual permitió aislar el cuerpo de agua. A partir de una máscara creada con la aplicación Envi del LG, se utilizaron las bandas de control de calidad QA y se delimitó la superficie acuática del lago, al utilizar solo píxeles clasificados efectivamente como agua.
Selección de imágenes MODIS. Para seleccionar las imágenes que serían utilizadas, las bandas de control de calidad (QA), junto con rutinas computacionales desarrolladas en Python, permitieron eliminar las imágenes del LG que presentaban más de 10 % de errores en sus píxeles. Este proceso permitió obtener la imagen para realizar la calibración y validación (imagen del 24 de abril) y seleccionar también la mejor imagen MODIS por año (20 imágenes).
Antes de la generación del algoritmo predictivo, se realizó la conversión de los números digitales de las imágenes a valores de reflectancia de superficie. Este valor se derivó a partir de la ecuación 1.
Donde: P¡ es el número digital registrado en la banda MODIS y Pf es el número digital final.
Calibración y validación del modelo. La imagen MODIS que mejor calidad radiométrica presentó en el periodo del trabajo de campo (enero a junio de 2013) se obtuvo el 24 de abril de 2013. Una vez realizados los procesos de preprocesamiento y extracción del cuerpo de agua, fueron identificados los píxeles que coinciden con las coordenadas geográficas de los puntos de muestreo. Este vector espectral extraído de la imagen fue utilizado para evaluar iterativamente modelos lineales y polinomiales (segundo y tercer grado) en las bandas 1 y 2 (valores de reflectancia B x y sus inversos B x -1 ) del sensor MODIS.
Con los vectores extraídos de la imagen y los datos de campo de la PDS se realizó el modelo de calibración y validación. Esta fase es fundamental en un proceso de estimación [47], [48], ya que sus resultados permiten evaluar la confiabilidad del algoritmo y las incertidumbres asociadas a los valores calculados. Por tal razón, teniendo en cuenta la limitante en el número de imágenes coincidentes con los trabajos de campo y que, además, el modelo debe ser evaluado con observaciones independientes, fue necesario utilizar el procedimiento Loocv, para analizar el desempeño del modelo de estimación.
En estos procedimientos, excepto uno de los pares de datos de la relación, todos son utilizados para validar el modelo. A continuación, se calculó el error de estimación (PDS estimado frente a PDS observado), a partir del par excluido. El procedimiento se repitió excluyendo, a su vez, todos los N pares disponibles. Adicionalmente, para verificar la validez del modelo, se utilizaron los estimadores estadísticos R2 y RMSE (detalles del procedimiento en [35], [48]-[51]).
Visualización temporal de datos. Para obtener observaciones realistas en las series temporales generadas por sensores remotos, es necesario, por un lado, utilizar observaciones de referencia tomadas en campañas de campo de muestreos in situ (dominio de origen); por otro, asociar los modelos de calibración obtenidos extrapolando esos modelos en la colección de imágenes en el tiempo (dominio de destino) [52]. En esta investigación, el dominio de origen fue la imagen MODIS relacionada con el trabajo de campo de abril de 2013 y el dominio destino las demás imágenes de la serie temporal. Con la aplicación de los modelos predictivos en las imágenes del dominio destino, ahora estas representaron las estimaciones del PDS.
Una vez obtenida esta información, se extrajeron vectores temporales coincidentes con las coordenadas geográficas de los puntos de muestreo; además, utilizando rutinas computacionales del lenguaje de programación Python y R se generó la visualización de datos HorizonGraph (HG) para el periodo 2001-2020.
La HG es una técnica diseñada para optimizar el espacio y la comprensión de series del tiempo multivariadas o con múltiples puntos de observación [53]-[55]. Dentro de la gráfica, se utilizan colores contrastantes (rojo y azul) y cada cambio de color indica que el valor leído en ese punto está por encima o por debajo del promedio. Además, se utilizan tres tonalidades de cada color que muestran la intensidad de la variable observada en ese momento. Este tipo de gráficas permite la presentación simultánea de mayor volumen de datos y la comparación entre series de tiempo multivariadas o univaridas multipunto [54], [56]-[59].
Resultados y discusión
Parámetros de la calidad del agua
Los parámetros de calidad del agua fueron tomados en seis campañas de muestreo, llevadas a cabo entre los meses de enero y junio de 2013. En la tabla 1, se muestran algunas medidas de estadística descriptiva para la PDS medida en campo.
Estación | Promedio | Máximo | Mínimo | Desviación estándar |
---|---|---|---|---|
1 | 4,37 | 5,1 | 3,7 | 0,60 |
2 | 4,23 | 5,6 | 2,7 | 1,02 |
3 | 4,87 | 6,2 | 3,4 | 1,08 |
4 | 4,93 | 7 | 2,9 | 1,33 |
5 | 4,70 | 6,4 | 3,4 | 1,15 |
6 | 4,67 | 7,4 | 3,1 | 1,61 |
7 | 4,65 | 8 | 3,2 | 1,84 |
8 | 4,6 | 7 | 3,2 | 1,53 |
9 | 4,08 | 6,5 | 2,4 | 1,58 |
10 | 3,45 | 4,3 | 2,2 | 0,82 |
Profundidad en metros, medida con el disco de Secchi.
Fuente: elaboración propia.
En general, se observa que el promedio de PDS en los puntos medidos in situ varió entre 3,45 y 4,93 metros. En los datos, se infiere que junio fue el mes que presentó valores mínimos en las lecturas, mientras que los máximos medidos para este parámetro se distribuyeron entre los meses de febrero y marzo. La desviación estándar en los puntos observados osciló entre 0,60 y 1,84 m. Los puntos de muestreo con mayor y menor profundidad fueron el siete y el diez, respectivamente.
Sensibilidad del modelo
En este trabajo, se desarrolló un modelo empírico para la recuperación de la PDS en el lago Guamués, a partir de la Rs captada por la banda 2 del sensor MODIS. Para evaluar el rendimiento del modelo de estimación de la PDS propuesto, se utilizaron siete pares de datos; píxeles banda 2 (858,5 nm) -datos de campo de PDS-; se excluyeron tres puntos (8, 9 y 10), debido a que algunos no fueron medidos en el trabajo de campo (punto 8) y, en otros casos, se presentaron errores en los píxeles extraídos de la imagen MODIS (puntos 9 y 10).
En total, fueron evaluados doce modelos que usan las bandas 1 y 2 (valores de reflectancia B x y sus inversos B x -1 ) del sensor MODIS para estimar los valores de PDS; los estimadores estadísticos RMSE y R2 fueron utilizados para evaluar la sensibilidad de los algoritmos y sus resultados se muestran en la tabla 2.
Fuente: elaboración propia.
A partir de los datos de la tabla 2, se observa que los modelos lineales presentan desempeños modestos en los estimadores estadísticos R 2 , sobre todo, los errores asociados al modelo predictivo que, en algunos casos, llegan a superar los 5,32 m en el RMSE. Por otra parte, los modelos basados en polinomios de tercer grado muestran un R 2 para la banda 2 superior a 0,9. Sin embargo, el RMSE es muy alto, pues, en algunos casos, alcanza los 29,2 m.
Teniendo en cuenta lo anterior, el mejor algoritmo predictivo es el modelo cuadrático basado en la banda 2. Este algoritmo cuenta con un ajuste moderado de R 2 = 0,79 en la calibración y R 2 = 0,74 en la validación, lo que indica que el modelo puede explicar el 74 % de la varianza de la PDS dentro del LG, con errores asociados al modelo predictivo inferiores a los 0,013 m. Por ello, muestra una alta sensibilidad a los cambios en la reflectancia capturada por la imagen MODIS.
Los resultados obtenidos por el modelo utilizado aquí, en muchos casos, son superiores a los que se presentan en otros estudios, como es el caso de la investigación desarrollada en el lago Rawal en Pakistán [60] o de los algoritmos utilizados para estimar PDS en los lagos de Michigan en Estados Unidos [26], donde los coeficientes de determinación (R2) oscilaron entre 0,49 a 0,79 y 0,32 a 0,71, respectivamente.
Una vez probada la validez del algoritmo basado en la banda 2 del sensor MODIS, fue modelado en términos matemáticos (ecuación 2) y representado gráficamente (figura 3). Asimismo, la validación del modelo presentada en la figura 4 muestra la dispersión de los valores de PDS estimados y observados obtenidos a partir del procedimiento LOOCV.
Distribución espacial
Generalmente, la variabilidad en la calidad del agua en un lago no uniforme es fácilmente visible a través de imágenes derivadas de sensores remotos (figura 5). Esta imagen MODIS coincidente con la fecha del trabajo de campo de abril de 2013 se utilizó para producir un mapa que muestra la dinámica espacial de la PDS sobre la superficie del lago Guamués. En esta figura, se revela que, aun siendo la resolución espacial de 250 m para los píxeles del sensor MODIS, fue posible captar la homogeneidad que el lago presenta en gran parte del espejo de agua, así como las zonas que presentan niveles que indican la existencia de un gran flujo de sólidos en suspensión, llevados al lago por algunos cursos de agua de magnitud considerable.
Realizando un análisis más detallado de la imagen MODIS, se tiene que la profundidad del disco Secchi o PDS promedio del lago fue de 3,87 m y los valores más altos se localizaron, como se esperaba, en la parte central del lago. Esto se debe a que esta zona se caracteriza por sus aguas profundas. Otro punto interesante que presenta valores superiores al promedio es la desembocadura del río Guamués. Esta área muestra un valor superior a los cuatro metros para la PDS.
Por último, los valores por debajo del promedio de la PDS, se ubican particularmente en las desembocaduras de las quebradas Flautal, Romerillo y El Motilón. En estos lugares, el algoritmo de estimación calculó medidas inferiores a los 0,5 m.
Visualización de la serie del tiempo en HorizonGraphs
La visualización HorizonGraphs (HG) se utiliza para observar el comportamiento simultáneo de los vectores temporales de los píxeles extraídos de las imágenes MODIS coincidentes con las coordenadas geográficas de las estaciones de muestreo para el parámetro PDS.
La visualización HG permite, de un lado, observar la dinámica temporal de los diez puntos de muestreo en conjunto y, de otro, apreciar con un mejor nivel de detalle los lapsos de tiempo y los puntos que presentan mayores y menores variaciones a lo largo de la serie temporal. La figura 6 utiliza esta técnica de visualización para indicar la variabilidad de las diez estaciones de muestreo para la variable PDS. Para las diez series de tiempo observadas, se tiene que su promedio es 3,32 m y su desviación estándar es de 1,44 m.
Los periodos en que se observa disminución de los valores de la PDS están resaltados en líneas discontinuas rojas. En la serie del tiempo, se observan dos periodos con valores menores del promedio en 2005-2007 (a) y 2017-2020 (b). Estos momentos de bajo PDS se hacen más notorios en la estación de muestreo siete, donde para los periodos indicados se tiene un PDS menor a 0,42 m. Otro punto con alta variabilidad es el seis, pero la oscilación dentro de esta estación de muestreo está en 1,2-3,0 m.
Otros puntos que presentan valores oscilantes entre 1,5 y 3,2 m son los puntos cuatro y cinco. Los restantes puntos de muestreo indican lecturas homogéneas con poca variabilidad, con una PDS promedio de 3,9 m.
Espacialmente, se observa que los puntos más - dinámicos son los comprendidos entre las estaciones de muestreo 4-7.
Conclusiones y recomendaciones
Los productos MOD09GA fueron muy útiles para crear una línea base que explica la variabilidad anual de la transparencia en el lago Guamués, no obstante, la continua oclusión por nubes dentro del área de estudio no permitió verificar la variación interestacional de la PDS dentro del cuerpo de agua en las épocas de menor y mayor precipitación.
Incluso con esta limitante, fue posible construir un modelo bio-óptico con una moderada sensibilidad entre la reflectancia captada por el sensor y las medidas realizadas en campo. De la misma forma, las imágenes MODIS permitieron observar que existen cursos de agua que aportan una gran cantidad de sedimentos de forma continua y en diferentes años, como las quebradas el Flautal, Romerillo y Motilón, cercanas a los asentamientos humanos de Mojondinoy y Santa Teresita (Flautal), Romerillo (Romerillo) y la vereda el Carrizo (Motilón).
El mapa generado a partir de los productos MO-DIS del 23 de abril de 2013 permite observar que, después de realizado un conteo con píxeles menores a un metro, este valor llegó a representar el 15 % del área total del lago. La serie de tiempo visualizada en HG también fue de utilidad para comprender la dinámica espacio temporal del lago, indicando los años 2006, 2007 y 2017-2020, y las zonas con mayor dinámica y variabilidad (puntos 3-6).
Para futuros estudios, se recomienda usar drones para el monitoreo dentro del lago Guamués, lo cual anulará la oclusión ocasionada por las nubes y maximizará la cobertura temporal de imágenes para realizar estudios de mayor precisión.