INTRODUCCIÓN

La evaluación de la calidad del suelo, junto con su fertilidad, es sumamente importante para tomardecisiones que permitanoptimizar la producción de alimentos y así contribuir a la capacidad de adaptación al cambio climático. Los métodos de laboratorio convencionales permiten caracterizar la física, química y biología de los suelos. Sin embargo, la exploración en uso de métodos novedosos permite, eventualmente, contar con métodos más rápidos, precisos y accesibles económicamente para evaluar la calidad del suelo. Actualmente se utilizan muestras de suelo y métodos de laboratorio convencionales como la Espectroscopia de Absorción Atómica o Espectroscopia de emisión por plasma de acoplamiento inductivo y colometríapara cuantificar el contenido de nutrientes en suelo. Esos métodos convencionales involucran una serie de pasos e insumos necesarios para poder obtener resultados.

Sin embargo, la espetroradioscopia visible e infrarroja ofrece una técnica con menos pasos e insumos, por lo cual tiene potencial para ser menos costosa y más rápida en comparación con los métodos de análisis convencionales. Con frecuenciaesos métodos convencionales requieren preparación de muestras que son propensas a la contaminación y errores durante el muestreo, la manipulación, la extracción, etc.,aspectos que instan a la aplicación de métodos más rápidos y accesibles para evaluar la calidad del suelo.

La espectroradioscopia visible e infrarroja (Vis-NIR-SWIR) ha surgido como opción para caracterizar eficazmente varias propiedades del suelo como carbono orgánico, pH, textura, capacidad de intercambio catiónico, contenido de CaCO₃, entre otros (^{Angelopoulouet al. 2020}, ^{Curcioet al. 2013}, ^{Demattêet al. 2016}, ^{Liu et al. 2016}, ^{Miloš y Bensa 2017}, ^{Silva et al. 2019}, ^{Steinberg et al. 2016},^{Viscarra Rossel et al. 2016}).

Las características de absorción en el espectro visible están dominadas por las excitaciones electrónicas moleculares y las de la gama, NIR contienen una combinación de vibraciones moleculares de armónicos (^{Adeline et al. 2017}). Sin embargo, caracterizar el efecto de cada una de las propiedades del suelo en los espectros de Vis-NIR, es una tarea difícil debido a la naturaleza compleja de la matriz del suelo con múltiples superposiciones de características espectrales y colinealidades espectrales fuertes entre las propiedades del suelo (Adeline et al. 2017). Las firmas espectrales codifican información sobre la composición inherente del suelo, que comprende la composición mineral, el contenido de nutrientes, los compuestos orgánicos y el agua.

El objetivo principal de este proyecto fueconstruir bibliotecas espectrales para los suelos tropicales de Costa Rica, y determinar las bandas hiperespectrales óptimas en el rango espectral visible (400nm - 700nm), infrarrojo cercano (NIR, 701-1000 nm) e infrarrojo de onda corta (SWIR, 1001-2500 nm)para desarrollar modelos predictivos de propiedades de suelo mediante firmas hiperespectrales.

^{Rinnanet al. (2009}) y ^{Hadouxet al. (2014}) dan una buena reseña de diversos métodos matemáticos para predecir las propiedades del suelo a partir de los espectros Vis-NIR del suelo, como la regresión lineal múltiple por pasos (^{Shibusawaet al. 2001}), la regresión por componentes principales (^{Chang et al. 2001}), el análisis de regresión múltiple (^{Ben-Dor y Banin 1995}), lossplines de regresión adaptativa multivariable (^{Shepherd y Walsh 2002}) y la regresión por support vector machine (^{Stevens et al. 2010}). El método de regresión de mínimos cuadrados parciales (PartialLeast Square Regresión - PLSR) es la técnica estadística multivariada más común para la calibración espectral y la predicción de las propiedades del suelo (^{Viscarra Rossel et al. 2006}). En este análisisse evaluó la capacidad de modelos PLSR para la predicción de 29 propiedades del suelo de muestras pre-procesadas en laboratorio a partir de firmas hiperespectrales Vis-NIR-SWIR generadas por un espectróradiometro (ASD).

MATERIALES Y MÉTODOS

Compilación de la librería hiperespectral

Los reportes de laboratorio para las propiedades de suelo (249 en total)se encontraban en diferentes formatos Excel,según el tipo de análisis. En otro orden, las 5500 firmas hiperespectrales generadas por el equipo FieldSpec 4,y a través del programa RS3, fueron guardadas en el formato binario (archivo “*.asd”),propietario de la compañía ASD. Se escribió un programa en MATLAB para extraer los resultados de los análisis de suelo de los archivos Excel yla depuración de los datos erróneos (problema de sintaxis, formato o a fuera de rangos) yluego se procedió a guardar cada resultado de forma sistemática en una base de datos. El programa también favoreció el cálculo de diferentes propiedades de suelo comocapacidad de intercambio catiónico efectivo (CICe), saturación de bases (SB), saturación de acidez (SA), relaciones catiónicas (Ca/Mg, Ca/K, Mg/K, (Ca+Mg)/K) y relación de absorción de sodio (RAS). Cada muestra tuvo un código único, el cual permitió integrar las firmas hiperespectralesal asociar el promedio de las 4 firmas sacadas a 90º de orientación por cada muestra y al seguir este código único. El programa de MATLABtambién permitió eliminar errores provenientes del cambio de sensores en el equipo FieldSpec 4 en las longitudes de onda de 1000 y 1800 nm (Figura 1).

Figura 1.  Efecto del algoritmo de LOWESS para eliminar errores provenientes del cambio de sensor en el equipo FieldSpec 4 a la longitud de onda de 1800 nm. 

En estas longitudes de onda, las firmas presentan un cambio abrupto que fue suavizado con algoritmo de regresión local (también conocido por sus siglas en inglés, LOWESS - LocallyWeightedScatterplotSmoothing). Adicionalmente, el código de MATLAB logra eliminar las reflectancias de los primeros 50nm (o sea de 350 a 400nm), debido a que la firmas contienen mucho ruido para estas longitudes de onda. Finalmente, el programa compila la librería hiperespectral en formato MATLAB (archivo “*.MAT”) y Excel (Archivo “*.xlsx”) para el análisis PLSR.

Análisis PLSR y desarrollo de Modelos Predictivos

PLSR es una técnica de modelado popular utilizada en quimiometría y se utiliza comúnmente para análisis espectral cuantitativo; contribuye a la formulación de modelos predictivos cuando hay diversidad de variables predictoras que son altamente colineales. La técnica está estrechamente relacionada conregresión de componentes principales (PCR). El método PLSR se basa en los mínimos cuadrados parciales (PLS), yaque permite reducir el espacio dimensional al maximizar la covarianza entre las puntuaciones de los espectros y cada propiedad del suelo. Este método realiza proyecciones rotativas iterativas para estimar matrices de variables latentes y puntuaciones asociadas tanto al espectro como a la variable de respuesta (propiedad de suelo).La regresión de mínimos cuadrados parciales es la asociación de una reducción de PLS con una regresión lineal multivariada clásica, que explica la correlación entre los espectros Vis-NIR y la propiedad del suelo. Los modelos predictivos tienen la forma de:

S = C ₄₀₀ R ₄₀₀ + C ₄₀₁ R ₄₀₁ + ...+ C _i R _i + m

Donde:

S = a la estimación de la propiedad de suelo.

C = al coeficiente de regresión PLSR para cada longitud de onda.

R= a la reflectancia generada por el espectroradiómetro (de 400 a 2500 nm).

m = alaintersección.

Antes del análisis estadístico, los espectros de reflectancia (R) se centraron en el promedio más/menos 2 desviaciones estándares para eliminar firmas atípicas (Figura 2). Luego, para cada propiedad del suelo, la base de datos de cada configuración espectral se dividió en un conjunto de calibración de 80% de las muestras y un conjunto de validación de 20% de las muestras,para predecir propiedades de suelo, con un juego de datos independientes a la calibración del modelo. Los algoritmos de PLSR fueron desarrollados en MATLAB.

Figura 2 Eliminación de las firmas hiperespectrales atípicas. a) Firmas antes y b) después de aplicar el filtro.  

Las líneaspunteadas representan el promedio de las firmas y las discontinuadas representan más o menos 2 desviaciones estándares del promedio.

RESULTADOS Y DISCUSIÓN

Se predijeron 29 propiedades del suelo a partir de un conjunto de datos de 1375 muestras de suelo pre-procesadas en laboratorio recolectadas en Costa Rica; su pH, acidez intercambiable (AE), K, Ca, Mg, P, Fe, Cu, Zn, Mn, Na, Si, B, S, C, N, Materia Orgánica (MO), textura (porcentaje de arena, limo arcilla), densidad aparente (DA), capacidad de intercambio catiónico (CICe), saturación de bases (SB), saturación de acidez (SA), relaciones catiónicas (Ca/Mg, Ca/K, Mg/K, (Ca+Mg)/K) y relación de absorción de sodio (RAS).Estas propiedades fueron seleccionadas por sus diferentes comportamientos espectroscópicos espectrales y predichos porPLSR.

La Figura 3 presenta los resultados de las diferentes métricas (R2, RMSE, PRESS, AIC, Wold’s R, AIC Normalizado y BIC) para la variable Ca como ejemplo. Dichas métricas fueron utilizadas para evaluar la calidad de las predicciones PLSR y determinar la cantidad óptima de componentes en los modelos para cada una de las propiedades de suelo. Los mejores modelos presentaron valores de R2,próximos a 1(uno) tanto con la dinámica de datos de calibración (80%), como con los datos de validación (20%), que llegaron a minimizar los valores de PRESS AIC y BIC. Las líneas verticales punteadas muestran la cantidad óptima de componentes que se emplearon en los modelos de predicción para cada propiedad de suelo.

Figura 3 Resultados de las métricas para evaluar la calidad de las predicciones PLSR y determinar la cantidad óptima de componentes en los modelos para calcio.  

La línea punteada corresponde al número de componentes óptimos para cada modelo.

Se determinó a qué longitud de onda se encontraban los coeficientes PLSR más sobresalientes para la descripción y la predicción de cada propiedad de suelo;posteriormente se eligieron los 10 picos con el valor absoluto más alto. La Figura 4 muestra resultados para pH, Mg, y C. Las flechas indican los picos que tienen mayor peso para describir el parámetro de interés. Es interesante notar que ciertas longitudes de ondas muestran rasgos similares entre losparámetrosque fueron analizados.

Figura 4.  Coeficientes PLSR versus longitud de onda para a) pH del agua, b) Magnesio y c) Carbono.  

Los coeficientes de regresión PLS pueden utilizarse para seleccionar predictores (longitud de onda) más relevantes de acuerdo con la magnitud de sus valores.

Las longitudes de mayor importancia en los modelos PLSR están indicadas por las flechas.

La Tabla 1 presenta las 10 longitudes de onda de mayor relevancia para predecir cada propiedad de suelo. No fue posible reconocer un patrón claro, pero si se superpusieron los 10 picos más significativos para el pH, los macronutrientes y micronutrientes y CICe (Figura 5). Se puede observar que las bandas se agrupan en 4 zonas del espectro del FieldSpec 4. Más precisamente, se puede observar que las bandas en el espectro visible tienen un papel importante en las predicciones PLSR así que 3 zonas de las longitudes de onda corta infrarroja (SWIR de 1001 a 2500 nm).

Tabla 1 Longitudes de onda más importantes para las predicciones PLSR para cada propiedad de suelo. 

Figura 5 Longitudes de onda de mayor importancia para los modelos PLSR. 

a) todos los elementos confundidos, b) para solo pH del agua y Acidez Extractable (AE), c) para fósforo y potasio y d) para magnesio y calcio.

Como en el estudio de ^{Adelineet al.(2017}), se puede agrupar las bandas espectrales más sobresalientes para predecir propiedades de suelo en 4 rangos espectrales. La primera zona se ubica de 400 a 750 nm;según Adelineet al.(2017) esta zona responde a la presencia de óxidos de hierro en el suelo. La segunda zona incluye longitudes de onda de 1350 a 1450 nm;según ^{Poppielet al. (2018}); ese rango responde a la presencia de minerales arcillosos tipo 2:1 (ej. caolinita) y agua adsorbida en la superficie de las partículas de suelo. La zona 3 de 1850 a 2000 nm estaba correlacionado con propiedades asociadas al contenido de agua (Adelineet al.2017). La cuarta zona incluye bandas espectrales cercanas a 2200 nm y está relacionado con la presencia de CaCO₃, minerales de arcilla en el suelo y materia orgánica.

La Figura 6expone losgráficos de dispersión de valor predicho frente al valor observado para algunas propiedades de suelo. La línea punteada negra en diagonal representa el eje de predicción perfecta. La punteada color gris obedece ala de ajuste del mejor modelo PLSR para esta variable. Entre más cerca se encuentren las 2 líneas punteadas, mejor ajuste predictivo. Las 2 continuas exponenuna desviación estándar,aspecto que permite identificar los datos atípicos. El valor de R² sin paréntesis (esquina izquierda superior de cada gráfico)muestrael coeficiente de determinación para todos los datos disponibles para esa propiedad de suelo (“n”es la cantidad de muestras usadas para desarrollar el modelo PLSR), y el número entre paréntesis es el coeficiente de determinación para las predicciones con los datos apartados para la validación que correspondió al 20% de todos los datos. Por su parte el RMSE mostró el error cuadrático medio de predicción en porcentaje.

Figura 6 Gráficos de dispersión describiendo el ajuste de las predicciones con los resultados de laboratorio. 

pH del agua, b) calcio, c) magnesio, d) fosforo, e) hierro, f) cobre, g) nitrógeno y h) capacidad de intercambio catiónico efectivo para muestras de suelo pre-procesadas tomadas en Costa Rica mediante el modelo PLS

La Tabla 2aporta un resumen de la calidad de los modelos para la mayoría de las propiedades de suelo. El filtro aplicado a las muestras consistió en centrar el promedio a 2 desviaciones estándar, para eliminar firmas irregulares. Algunos elementos y/o propiedades de suelo tenían menos análisis de laboratorio, por lo que se realizó la modelización con menos datos,que corresponden a la segunda columna de la Tabla 2.

No se eliminaron datos de calibración para tratar de ajustar los modelos. Se generaron buenos modelos PLSR (R²>0,8, RMSE<10%) para Ca, Mg, Fe, C, N y CICe.

El número óptimo de componentes para el modelo PLSR se escogió basado en R², esto debido a que es el más, comúnmente, utilizado en literatura, lo que permite realizar comparaciones entre los resultados obtenidos y otras investigaciones. Los resultados para los modelos de parámetros físicos de suelo como el porcentaje de arena, limo y arcilla, y densidad aparenteindicaron una pobre capacidad de predictiva (R²<0,5, RMSE entre15-20%). Los resultados para carbono coincidieron con ^{Reyna et al. (2017}). Este análisisfue realizado con un número de 70 muestras, que aunque limitado sereportó que la modelización PLSR es útil para obtener información sobre el contenido de carbono. ^{Kawamuraet al. (2019}) obtuvieron por su parte, predicciones para fósforo (R²=0,78) al combinar un análisis PLS con un algoritmo genético con solamente 103 muestras. En el caso de esta investigación, se contócon una base de datos mucho más sólida con más de 1300 muestras.

La manipulación digital de cada propiedad de suelo,se identificó como limitante, ya que el tiempo promedio para evaluar los primeros 40 componentes de los modelos PLSR, con un procesador con 4cores (Intel Core i7-4712HQ) y 16,0 GB de RAM, fue de 17 horas por cada propiedad de suelo.En ese sentido, ^{Adeline et al. (2017}) investigaron la posibilidad de reducir la resolución de las bandas hiperespectrales y obtuvieron buenos resultados. Ese proceso es conocido como degradación espectral y se presenta como una opción para llegar a reducir el tiempo de utilizado.

Con los modelos que se generaron, se pudieron predecir las variables de suelo analizado a partir de su firma hiperespectral con certidumbre, como fue en el caso de calcio, magnesio, hierro, carbono, nitrógeno y CICE, que reportaron R² de validación, superior al 0,8.Se escribió un programa en MATLAB que usa los coeficientes PLSR generados para predecir cada propiedad de suelo con su firma espectral en menos de 5 minutos. Esta técnica facilitó una alternativa para aproximar rápidamente las deficiencias nutricionales que puede tener un suelo degradado a partir de muestras con un pre-procesamiento previo como secado y molido. A futuro, se consideraincorporar técnicas de machine learningpara cuantificar las propiedades de los suelos de Costa Rica y considerar así un mejor ajuste del modelo predictivo.

CONCLUSIONES

Los resultados obtenidos mostraron que la espectroscopía Vis-NIR-SWIR puede predecir de forma relativamente rápida en comparación a las técnicas convencionales de análisis de suelo con extracción en Mehlich III u Olsen y con espectroscopia de absorción atómica o espectroscopia de emisión por plasma de acoplamiento inductivo,una amplia gama de propiedades del suelo, independientemente de las diferentes ubicaciones geográficas en Costa Rica. A partir de los coeficientes PLSR encontrados en este estudio, se ha trabajado en un programa basado en MATLAB que permite predecir el valor de cada propiedad del suelo a partir de su firma hiperespectral. Se logró una estimación precisa del contenido de diferentes componentes (Ca, Mg, Fe, C, N y CICe) con un R² superior a 0,8 y un error cuadrático medio (RMSE) inferior a 10%.

Se identificó que los análisis espectroscópicos combinados con la regresión PLS pueden proporcionar una herramienta muy útil para la agricultura de precisión en los suelos tropicales de Costa Rica.A futuro, la capacidad predictiva se complementará con técnicas de pre-procesamiento, con el análisis de regresión PLS y otros algoritmos de aprendizaje automático, que permitangenerar las bases para investigaciones en la predicción de las variables analizadas en campo.