30 de septiembre de 2026
Soporte a la decisión en remediación de áreas contaminadas con algoritmos de machine learning
Soporte a la Decisión en Remediación de Áreas Contaminadas con Algoritmos de Machine Learning
Deivid Lucas dos Santos Migueleti; Aimée Shirozono
DOI: 10.22167/2675-6528-202602835
Artículo derivado de Trabajo de Fin de Curso (TCC), con contenido basado en el trabajo original del estudiante y adaptado al formato editorial de la Revista E&S con el apoyo de la herramienta ResumeAI, una solución de inteligencia artificial desarrollada por el Instituto Pecege para la síntesis y organización textual.
Resumen
La remediación de áreas contaminadas por hidrocarburos en acuíferos sujetos a oscilaciones estacionales del nivel freático demanda instrumentos cuantitativos para la planificación y optimización, cuya eficiencia se ve limitada por el atrapamiento de derivados del petróleo. El trabajo demostró la aplicabilidad de modelos estadísticos de aprendizaje automático en la predicción de series temporales de la carga hidráulica en un área contaminada por queroseno de aviación en el polo petroquímico de Paulínia-SP, como evidencia para el soporte a la decisión en remediación. Se utilizaron registros mensuales de carga hidráulica de pozos de monitoreo y precipitación pluviométrica acumulada mensual de un período de 10 años. El análisis del comportamiento de la serie temporal, la ingeniería de variables y el análisis de correlación auxiliaron en la estimación de modelos SARIMA y SARIMAX. Se seleccionaron y ajustaron modelos con capacidad predictiva adecuada basándose en criterios de parsimonia y diagnóstico de residuos. El modelado de series temporales con variables exógenas mejoró el rendimiento de los modelos en comparación con los modelos SARIMA. El modelo con precipitación acumulada en 3 meses presentó un error medio absoluto (MAE) de 0,1041 metros para la previsión de 12 meses, considerado satisfactorio. Esto constituyó una prueba de concepto relevante para apoyar la planificación de la remediación y la gestión a largo plazo de áreas contaminadas por hidrocarburos.
Palabras clave: Agua subterránea; Contaminación; LNAPL; Precipitación; Suelo.
1. Introducción
El agua subterránea representa el 99% de las reservas de agua dulce del planeta, siendo un recurso estratégico fundamental para la adaptación climática de las sociedades. En este contexto, las políticas públicas orientadas a la protección y recuperación de acuíferos subterráneos son de extrema relevancia (UNESCO, 2022).
La contaminación del suelo y del agua subterránea, sin embargo, puede comprometer seriamente el uso de este recurso vital. En el Estado de São Paulo, el Sistema de Información de Áreas Contaminadas y Rehabilitadas (SIACR), mantenido por la Compañía Ambiental del Estado de São Paulo (CETESB), registraba, en octubre de 2025, 7.253 áreas contaminadas, de las cuales 4.846, aproximadamente el 66,8%, estaban afectadas por combustibles automotrices derivados del petróleo (CETESB, 2025).
Los derivados de petróleo presentes en el medio subterráneo en su forma pura se conocen como LNAPL, acrónimo de Light Non Aqueous Phase Liquid, o fase líquida no acuosa inmiscible. Estas mezclas complejas pueden contener sustancias tóxicas, como benceno y naftaleno, que contaminan el agua subterránea y representan riesgos para la salud humana (Alden et al. 2024).
El proceso de rehabilitación de un área contaminada implica etapas de intervención para mitigar o reducir los riesgos de exposición a sustancias químicas de interés (CETESB, 2017). La información sobre la distribución de los contaminantes y su interacción con el medio físico y los receptores se consolida en el Modelo Conceptual Ambiental (MCA) (Alden et al. 2024).
Las medidas de remediación utilizan tecnologías para tratar el medio afectado, eliminando los contaminantes o restringiendo su movilidad (contención) (CETESB, 2017). El Decreto Estatal de São Paulo número 59.263/2013 prioriza las medidas que promueven la eliminación y reducción de la masa de los contaminantes (ALESP, 2013). Las técnicas de remediación más comunes para áreas contaminadas por derivados de petróleo se basan en la eliminación hidráulica del LNAPL del acuífero, incluyendo extracción multifásica, bombeo y tratamiento, recuperación de fase libre y barrera hidráulica (ITRC, 2023; CETESB, 2017).
A pesar de su amplia aplicación, la eficiencia de las técnicas de remoción hidráulica está limitada por diversos factores y debe optimizarse a medida que se actualiza el MCA del área (ITRC, 2023). Uno de los fenómenos que restringe la recuperación de LNAPL es el atrapamiento, donde el LNAPL queda atrapado en los poros del suelo por debajo del nivel freático debido a las diferencias en las fuerzas capilares y la tensión interfacial entre el aceite y el agua (Isler et al. 2018; Teramoto et al., 2020). Este fenómeno fue observado por Pede (2009) en un área del polo petroquímico de Paulínia.
En regiones con regímenes pluviométricos bien definidos, el atrapamiento del LNAPL puede intensificarse por la oscilación estacional del nivel freático en respuesta a la precipitación (Teramoto et al., 2020). El cambio climático y los eventos extremos, como sequías prolongadas y lluvias torrenciales, añaden complejidad a la gestión de áreas contaminadas (UNESCO, 2022). Las precipitaciones muy por debajo del promedio intensifican el atrapamiento, dificultando la remediación por métodos hidráulicos.
Aunque existen desde hace tiempo modelos matemáticos y computacionales para simular y predecir el nivel del agua de acuíferos, su aplicación es limitada debido a la complejidad y la necesidad de recursos especializados. Investigaciones anteriores utilizaron Redes Neuronales Artificiales, Máquinas de Vectores de Soporte y modelos estadísticos como ARMA y ARIMA para predecir niveles de agua, algunos con capacidad de anticipar hasta cuatro meses (Tao et al., 2022; Nayak et al. 2006). Más recientemente, modelos basados en la arquitectura Transformer han demostrado una capacidad predictiva superior (Sun et al. 2024).
La necesidad de modelos de machine learning robustos y más accesibles, capaces de predecir el nivel del agua, es evidente. Dichos modelos, aliados al conocimiento sobre la distribución vertical del LNAPL, pueden proporcionar una línea de evidencia crucial para la toma de decisiones y la negociación con organismos ambientales sobre la discontinuidad de los esfuerzos de remediación en áreas contaminadas donde ha ocurrido el atrapamiento del LNAPL. El objetivo de este trabajo fue demostrar la aplicabilidad de al menos un modelo de aprendizaje automático en la predicción de series temporales del nivel del agua en un área contaminada por LNAPL como subsidio a la toma de decisiones sobre la continuidad de la remediación mediante técnicas de remoción hidráulica.
2. Material y Métodos
El presente estudio se caracterizó como una investigación de naturaleza cuantitativa, con un enfoque metodológico centrado en el análisis de series temporales. El objetivo fue demostrar la aplicabilidad de modelos de aprendizaje automático en la predicción del nivel de agua en un área contaminada por LNAPL, proporcionando un subsidio para la toma de decisiones sobre la continuidad de la remediación mediante técnicas hidráulicas. La investigación se llevó a cabo en un área contaminada por queroseno de aviación en el polo petroquímico de Paulínia, São Paulo, donde el derrame fue detenido en 2006 y el fenómeno de atrapamiento del LNAPL fue previamente observado por Pede (2009).
La unidad de análisis consistió en registros mensuales de carga hidráulica de 118 pozos de monitoreo (enero de 2013 a septiembre de 2025) y 180 registros de precipitación pluviométrica acumulada mensual (enero de 2011 a diciembre de 2025). Los datos fueron recolectados de pozos de monitoreo y de una estación meteorológica de superficie en la misma planta industrial. La utilización de estos datos, propiedad de la empresa, fue autorizada mediante término de anuencia. Los datos brutos fueron inicialmente procesados en Microsoft PowerBI Desktop para limpieza y tratamiento, y la carga hidráulica se calculó restando el nivel de agua de la cota altimétrica de la boca del tubo.
Los pasos subsiguientes de tratamiento, relleno de información faltante y preparación de los datos se llevaron a cabo en el lenguaje Python, en el entorno Google Colaboratory (Google, 2026), utilizando bibliotecas como Pandas, NumPy, Matplotlib, Seaborn, Scipy, Scikit-learn, Statsmodels y pmdarima. Para rellenar huecos en las series temporales, se empleó la función `interpolate` de Pandas. Para la estimación de los modelos, se utilizó el período de enero de 2016 a septiembre de 2025, excluyendo el período de la crisis hídrica de 2014-2015. A partir de las variables originales, se calcularon seis variables de carga hidráulica desfasadas (uno a seis meses) y cinco variables de precipitación basadas en el Índice de Precipitación Estandarizada (SPI) para escalas de uno, dos, tres, seis y doce meses, utilizando el paquete SPEI (Vonk, 2025), según McKee et al. (1993). Variables de precipitación acumulada se crearon con la función `shift` de Pandas.
Se realizó una etapa de clusterización exploratoria de los registros de carga hidráulica, empleando el algoritmo DBScan de la biblioteca Scikit-Learn (Ester et al., 1996), tras la estandarización de los valores con StandardScaler. Los parámetros `eps` y `min_samples` del DBScan se ajustaron basándose en un gráfico de k-vecinos. El análisis exploratorio de los datos incluyó análisis descriptivos, gráficos de box-plot, gráficos de línea y barras, y una matriz de correlación de Pearson, visualizada a través de un mapa de calor generado con el paquete Seaborn.
La descomposición de la serie temporal de la carga hidráulica se realizó utilizando la función `seasonal_decompose` de la biblioteca Statsmodels (Hyndman y Athanasopoulos, 2021). La autocorrelación y la autocorrelación parcial (FAC y FACP) se evaluaron con las funciones `plot_acf` y `plot_pacf` de la misma biblioteca. La estacionariedad de la serie temporal se verificó mediante las pruebas de Dickey-Fuller Aumentado (ADF) y Kwiatkowski-Phillips-Schmidt-Shin (KPSS), utilizando las funciones `adfuller` y `kpss` del módulo `tsa` de la biblioteca Statsmodels, al nivel de significancia del cinco por ciento (Morettin y Toloi, 2018). La serie de la carga hidráulica necesitó dos diferenciaciones para alcanzar la estacionariedad.
Para la estimación de los modelos, se utilizó la función `auto_arima` del paquete pmdarima, que minimiza los criterios de información de Akaike (AIC) y Bayesiano (BIC) (Morettin y Toloi, 2018). Se estimaron modelos Autorregresivos Integrados de Medias Móviles (ARIMA) y sus extensiones estacionales (SARIMA), así como modelos SARIMAX, que incorporan variables exógenas explicativas (Nelson, 1973). El ajuste de los modelos seleccionados consideró solo las variables exógenas con coeficientes estadísticamente significativos. El rendimiento se evaluó mediante predicciones de la carga hidráulica para doce meses, utilizando las funciones `get_forecast`, `predicted_mean` y `conf_int`. Las métricas de error utilizadas fueron el Error Medio Absoluto (MAE) y la Raíz del Error Cuadrático Medio (RMSE) (Hyndman y Athanasopoulos, 2021).
El diagnóstico de los residuos se realizó para verificar su representatividad, incluyendo la generación de gráficos de residuos estandarizados, histogramas, Q-Q plots y correlogramas. La verificación formal de las características de los residuos incluyó la prueba de Ljung-Box para autocorrelación, y las pruebas de Jarque-Bera, Shapiro-Wilk y Heteroskedasticidad para la adherencia a la normalidad y homoscedasticidad (Fávero y Belfiore, 2024).
3. Resultados y Discusión
El análisis de los datos del área contaminada por queroseno de aviación en Paulínia-SP reveló patrones hidrológicos y pluviométricos que influyen directamente en el comportamiento de la carga hidráulica del acuífero, proporcionando subsidios esenciales para la aplicación de modelos predictivos. Inicialmente, la exploración de los registros de carga hidráulica de los pozos de monitoreo permitió identificar agrupaciones espaciales, aunque la información detallada de los clústeres no fue incorporada en los análisis subsecuentes para enfocar en el comportamiento general del nivel del agua en el área. La distribución espacial de los pozos de monitoreo, como se observó, presentó una agrupación que seguía el sentido del flujo del agua subterránea, de noreste a suroeste, con algunos pozos clasificados como “outliers” debido a mediciones discrepantes.
La serie temporal de la carga hidráulica, compuesta por 106 registros de nivel de agua tras una interpolación del 9,40% de los datos faltantes, mostró una mediana de 589,092 metros y una media de 589,093 metros. Se aplicó interpolación lineal para rellenar huecos, una limitación metodológica reconocida, pero cuyo impacto se consideró limitado por la naturaleza de la respuesta del acuífero a la precipitación y por la corta duración de las fallas de medición. La distribución mensual de la carga hidráulica reveló un patrón estacional claro, con picos en los meses de marzo y abril y mínimos en octubre y noviembre, lo que indica la influencia estacional sobre el nivel del agua.
La precipitación acumulada mensual, con 117 registros, presentó una mediana de 75,000 milímetros y una media de 99,150 milímetros, con una leve asimetría a la derecha. La distribución mensual de la precipitación confirmó el patrón estacional esperado para la región, caracterizado por veranos lluviosos e inviernos secos. Esta estacionalidad es un factor crucial para entender la dinámica de la recarga del acuífero y, consecuentemente, las oscilaciones de la carga hidráulica, que son fundamentales para la planificación de la remediación de áreas contaminadas por LNAPL.
Correlación entre variables
La evaluación de la correlación entre las variables originales y construidas, incluida la carga hidráulica, la precipitación mensual, las variables rezagadas de la carga hidráulica y los índices de precipitación estandarizada (SPI), reveló relaciones significativas. Se observó una correlación del 41% entre los valores de la carga hidráulica y el SPI de 6 meses (spi6), del 41% con la precipitación acumulada en 2 meses (pluv_acum2) y del 54% con la precipitación acumulada en 3 meses (pluv_acum3). Estos resultados son consistentes con la comprensión de que la respuesta del acuífero a la precipitación no es inmediata, sino rezagada, debido a los procesos de infiltración y recarga de agua en el suelo.
El análisis de las series temporales de la carga hidráulica junto con la precipitación pluviométrica acumulada mensual, el SPI de 6 meses y la precipitación acumulada en 3 meses, demostró visualmente que la carga hidráulica tiende a aumentar en períodos de mayor precipitación, confirmando la relación de dependencia. Esta observación es vital para la modelización, ya que la inclusión de variables exógenas que capturen este desfase temporal puede mejorar significativamente la capacidad predictiva de los modelos, alineándose con la literatura que destaca la importancia de la precipitación en la dinámica de acuíferos (Teramoto et al., 2020).
Análisis de la descomposición de la variable respuesta
La descomposición de la serie temporal de la carga hidráulica, realizada mediante el método aditivo, permitió identificar claramente la presencia de componentes de tendencia, estacionalidad y un residuo. La tendencia observada en la carga hidráulica a lo largo del tiempo refleja cambios a largo plazo en el sistema hidrogeológico, mientras que la estacionalidad capta los patrones anuales de oscilación del nivel del agua, fuertemente influenciados por el régimen pluviométrico de la región. La identificación de estos componentes es un paso fundamental para la estimación de modelos de series temporales, ya que permite que los modelos capturen y expliquen estos patrones de comportamiento (Hyndman y Athanasopoulos, 2021).
Los gráficos de la Función de Autocorrelación (FAC) y de la Función de Autocorrelación Parcial (FACP) de la variable respuesta confirmaron la existencia de tendencia y estacionalidad en los datos. Los picos en la FACP indicaron que la serie temporal de la carga hidráulica puede ser modelada como un proceso autorregresivo de orden 2 o 3. Esta información fue crucial para la elección de los parámetros iniciales de la función `auto_arima` del paquete `pmdarima`, orientando la estructura de los modelos SARIMA y SARIMAX a ser estimados.
Pruebas de estacionariedad de series temporales
La evaluación de la estacionariedad de la serie temporal de la carga hidráulica se realizó mediante las pruebas de Dickey-Fuller Aumentado (ADF) y Kwiatkowski-Phillips-Schmidt-Shin (KPSS), al nivel de significancia del 5%. La prueba ADF inicialmente indicó la no estacionariedad de la serie original, con un p-valor de 0,286026, que es superior al nivel de significancia. Después de una diferenciación (d=1), el p-valor se redujo a 0,129114, aún indicando no estacionariedad. Solo después de dos diferenciaciones (d=2), el p-valor alcanzó 7,300095e-13, confirmando la estacionariedad de la serie.
La prueba KPSS, que asume la estacionariedad como hipótesis nula, también corroboró la necesidad de diferenciación. Con un p-valor inicial de 0,100000, la prueba no rechazó la hipótesis nula de estacionariedad, pero la combinación con la prueba ADF y el análisis gráfico de la serie diferenciada dos veces, que mostró un comportamiento estable en torno a la media, confirmó la adecuación del parámetro d=2 para la modelización. La estacionariedad es un requisito previo para la aplicación de modelos ARIMA y SARIMA, garantizando que la media y la varianza de la serie sean constantes a lo largo del tiempo (Morettin y Toloi, 2018).
Estimación y ajuste de los modelos
La estimación de los modelos SARIMA y SARIMAX se realizó utilizando la función `auto_arima` del paquete `pmdarima`, buscando minimizar los criterios de información de Akaike (AIC) y Bayesiano (BIC), que evalúan la parsimonia y la calidad del ajuste del modelo. Inicialmente, se estimaron modelos SARIMA, sin variables exógenas, para comparación. El modelo SARIMA(1,2,0)(2,0,1)[12] presentó un AIC de -15,377, mientras que el SARIMA(3,2,0)(2,0,0)[12] tuvo un AIC de -18,023, lo que indica un mejor ajuste.
La introducción de variables exógenas en los modelos SARIMAX demostró una mejora significativa en el rendimiento. El modelo SARIMAX(0,2,3)(2,0,0)[12], con todas las variables exógenas de la Tabla 1, alcanzó un AIC de -38,882, una caída sustancial en comparación con los modelos SARIMA puros. Sin embargo, este modelo presentó autocorrelación en los residuos, como lo indica la prueba de Ljung-Box (p<0,05), lo que sugiere que no toda la información fue capturada adecuadamente.
El análisis de los coeficientes de los modelos SARIMAX con todas las variables exógenas reveló que algunas variables SPI y `pluv_acum2in`, así como ciertos términos autorregresivos y de medias móviles, no fueron estadísticamente significativos. Esta observación puede atribuirse a la multicolinealidad entre las variables predictoras, un fenómeno común en series temporales (Fávero y Belfiore, 2024). Los modelos con variables exógenas como `pluv_acum3` (precipitación acumulada en los 3 meses anteriores) y `spi1` (índice de precipitación estandarizado de 1 mes) presentaron los menores errores, lo que es coherente con la respuesta retardada del acuífero a la precipitación.
El modelo SARIMAX(2,2,3)(2,0,1)[12] con la variable `pluv_acum3` destacó con un error medio absoluto (MAE) de 0,1041 metros para una previsión de 12 meses, considerado satisfactorio para el objetivo del estudio. Otro modelo notable fue el SARIMAX(2,2,2)(2,0,0)[12] con `spi1`, que obtuvo un AIC de -52,144 y BIC de -31,066, indicando un ajuste robusto. La selección final de los modelos para ajuste y diagnóstico de residuos se basó en los criterios AIC, BIC y la estadística Ljung-Box, priorizando la parsimonia y la significancia estadística de los coeficientes.
Diagnóstico de residuos
El diagnóstico de los residuos de los modelos ajustados se realizó para verificar su representatividad, según lo recomendado por Morettin y Toloi (2018). La evaluación de los gráficos de dispersión de los residuos estandarizados de todos los modelos ajustados reveló una distribución aleatoria alrededor de la media cero a lo largo del tiempo, lo que indica homocedasticidad. Esta característica fue confirmada formalmente por la estadística Heteroskedasticity, que, según la documentación de la biblioteca `Statsmodels`, es análoga a la prueba de Goldfeld-Quandt, cuya hipótesis nula es la no heterocedasticidad de los residuos.
Sin embargo, los gráficos Cuantil-Cuantil indicaron que la distribución de los residuos en todos los modelos no se adhirió completamente a la normalidad, posiblemente debido a la presencia de valores extremos. Esta falta de adherencia fue corroborada por las pruebas estadísticas de Jarque-Bera y Shapiro-Wilk, que presentaron p-valores inferiores a 0,05 para todos los modelos ajustados. Aunque la no normalidad de los residuos no invalida los modelos, puede tener implicaciones en la determinación de los intervalos de confianza (Fávero y Belfiore, 2024). La transformación de la variable dependiente, como la transformación de Box-Cox, podría mitigar este problema, pero no fue explorada en el presente trabajo.
En resumen, los resultados demostraron la aplicabilidad de modelos estadísticos de aprendizaje automático en la predicción de la carga hidráulica en un área contaminada, destacando los modelos SARIMAX que incorporan variables exógenas de precipitación. La capacidad de predecir el nivel del agua con un error absoluto medio de 0,1041 metros, utilizando la precipitación acumulada en 3 meses, representa una prueba de concepto robusta. Esta evidencia cuantitativa es crucial para apoyar la planificación de la remediación y la gestión a largo plazo de áreas contaminadas por hidrocarburos, especialmente en contextos de atrapamiento de LNAPL y oscilaciones estacionales del nivel del agua.
4. Conclusión
El presente trabajo demostró la aplicabilidad de modelos de aprendizaje automático en la predicción de series temporales de la carga hidráulica en un área contaminada por queroseno de aviación, proporcionando un subsidio para la toma de decisiones sobre la continuidad de la remediación por técnicas hidráulicas. Se verificó que la carga hidráulica del acuífero exhibe patrones estacionales claros y una respuesta desfasada a la precipitación pluviométrica. La estimación de modelos SARIMAX, que incorporaron variables exógenas de precipitación, mejoró significativamente la capacidad predictiva en comparación con modelos SARIMA puros. Un modelo SARIMAX, utilizando la precipitación acumulada en tres meses, alcanzó un error medio absoluto de 0,1041 metros para una previsión de doce meses, lo que se considera satisfactorio. Este resultado estableció una prueba de concepto robusta, ofreciendo una línea de evidencia cuantitativa crucial para la planificación de la remediación y la gestión a largo plazo de áreas contaminadas por hidrocarburos, especialmente en escenarios de atrapamiento de LNAPL y oscilaciones estacionales del nivel del agua.
A pesar de la relevancia de los hallazgos, el estudio presentó limitaciones, como la interpolación lineal del 9,40% de los datos faltantes, cuyo impacto en la sensibilidad de los modelos no fue evaluado, y la no adherencia total de los residuos a la normalidad, lo que puede influir en la determinación de los intervalos de confianza. También se observó multicolinealidad entre algunas variables exógenas de SPI, lo que resultó en coeficientes no significativos para ciertas variables. Para estudios futuros, se sugiere profundizar la investigación sobre la multicolinealidad y la utilización de los clusters identificados por el algoritmo DBScan en la estimación de los modelos. Se recomienda además probar enfoques como el Método de Predicción Global y explorar modelos más avanzados, como LSTM, NeuralProphet y aquellos basados en Transformers, para capturar no linealidades y mejorar la capacidad predictiva en pozos piezométricos específicos.
Referencias Bibliográficas
Alden, D.F.; García-Rincón, J.; Rivett, M.O.; Wealthall, G.P.; Thomson, N.R. 2024. Complexities of Petroleum Hydrocarbon Contaminated Sites. p. 1-19. In: García-Rincón, J.; Gatsios, E.; Lenhard, R.J.; Atekwana, E.A.; Naidu, R. Advances in the Characterisation and Remediation of Sites Contaminated with Petroleum Hydrocarbons. Environmental Contamination Remediation and Management. Springer, Cham.
Companhia Ambiental do Estado de São Paulo [CET
Companhia Ambiental do Estado de São Paulo [CETESB], 2025. Sistema de Áreas Contaminadas e Reabilitadas. Relatório de Áreas Contaminadas e Reabilitadas no Estado de São Paulo. Disponível em: < https://sigam.ambiente.sp.gov.br/sigam3/Default.aspx?idPagina=17653>. Acesso em: 10 out. 2025.
Organização das Nações Unidas para a Educação, a Ciência e a Cultura [UNESCO]. 2022. UN World Water Development Report 2022: Groundwater: Making the invisible visible. Disponível em: <https://www.unwater.org/publications/un-world-water-development-report-2022>. Acesso em: 12 out. 2025.
Artículo originario del Trabajo de Conclusión de Curso de la Especialización en Data Science y Analytics del MBA USP/Esalq
Para saber más sobre el curso, haga clic aquí y acceda a la plataforma MBX Academy