1. INTRODUCCIÓN

El pH es una medida de la acidez o alcalinidad de un suelo. Esta propiedad controla las reacciones químicas que determinan si los nutrientes van a estar o no disponibles para su absorción y, por lo tanto, afecta directamente el crecimiento y desarrollo de las especies de plantas que se encuentren sobre el suelo (Rivera et al, 2018). Además, se ha demostrado que cambios en el pH por el uso de fertilizantes y enmiendas puede alterar la composición de las poblaciones microbianas en el suelo (Liu, et al., 2020), lo que ocasionaría problemas en procesos como la fijación de nitrógeno. Por otro lado, conforme el suelo se hace ácido, el pH se incrementa en un factor de 10 por cada punto de caída. Esto significa que los efectos de un suelo bajo en pH son muy notables y rápidamente se disminuye el crecimiento de las raíces, la disponibilidad de fósforo y molibdeno, así como la capacidad de utilización del nitrógeno, fósforo, potasio y de los micronutrientes (Aguilar, et al., 2018). Igualmente, se ha evidenciado que el principal efecto de un pH inadecuado en el suelo es la reducción del crecimiento radical, es decir, no se logra un enraizamiento adecuado lo que conlleva a caída de plantas, poca o nula captación de nutrientes de baja movilidad en el suelo, incluso cuando la nutrición de la planta sea adecuada (Barbaro, et al., 2019). Por lo tanto, es importante desarrollar técnicas que permitan determinar de manera rápida, precisa y económica la variabilidad espacial del pH, lo cual sin duda ayuda a cuantificar el impacto de esta variación sobre la producción y las posibles pautas de manejo requeridas para optimizar los rendimientos.

Generalmente se han usado muestreos discretos para comprender las propiedades del suelo, sin embargo, entre las problemáticas que existen en la interpretación de tales análisis se encuentra su propia variabilidad, la cual está influenciada directamente por la dinámica nutricional, la topografía del terreno, la génesis del suelo y su manejo, entre otros. Debido a esto se han desarrollado herramientas que incorporan los datos reales para predecir las propiedades del suelo y permitir una representación continua de todo el campo (Radocaj, et al., 2020), entre estas se encuentran los Sistemas de Información Geográfica (SIG). Una de las utilidades más importantes que poseen los SIG es el uso de interpolaciones de diferentes tipos de variables. En el campo de la agronomía, el uso de esta herramienta permite analizar la variabilidad de diferentes características sobre el paisaje lo que contribuye a mejorar las prácticas de manejo (Henríquez, 2013). El uso de interpoladores posibilita analizar distribuciones y variaciones espacio-temporales de observaciones, basadas en el análisis de dependencia espacial entre los valores de una variable medidos en sectores vecinos, con el fin de generar estimaciones en otras posiciones dentro de la región estudiada (Pose et al., 2001). Los métodos de interpolación se dividen en deterministicos y geoestadísticos. Los primeros utilizan los valores muestreados y funciones matemáticas para la interpolación; mientras que, los métodos geoestadísticos también aplican funciones estadísticas para cuantificar la incertidumbre de las predicciones (Oliver, 2010). De acuerdo con lo anterior, comparar los métodos de interpolación permite selecionar el más óptimo y preciso para una zona de estudio en particular, lo cual tiene un impacto directo en la creación de planes de gestión del uso y manejo del suelo (Garbanzo, et al., 2017).

Finalmente, el propósito de este trabajo es comparar tres métodos de interpolación para el pH del suelo en el departamento de Cesar, de manera que sea posible utilizar el método con menor incertidumbre para futuras investigaciones y para mejorar el manejo de los suelos en tal zona y así lograr incrementar la productividad de los cultivos.

2. DESCRIPCIÓN DE LA ZONA DE ESTUDIO

El área de estudio es el departamento de Cesar, ubicado al noreste del país, en la llanura del Caribe; entre los 7°41’16’’ de latitud norte y los 72°53’27’’ y 74°8’28’’ de longitud oeste. Posee una extensión de 22.905 km^2, lo que representa 2% del territorio nacional y, según el DANE (2005), tiene una población de 903.279 habitantes. Las actividades agropecuarias ocupan el principal renglón económico del departamento, pues de éstas deriva 47% de sus ingresos (Arias & Acevedo, 2009), ya que el departamento cuenta con buena calidad de suelos y tiene un alto potencial de adecuación de las tierras mediante el riego (URosario, 2018). Por lo tanto, resulta importante conocer la variabilidad espacial de las propiedades del suelo, más específicamente del pH. Por otra parte, el departamento tiene una doble connotación: es Caribe y es fronterizo. La primera hace referencia a su identificación con los 8 departamentos del Caribe, aunque no tenga límites con el mar, ancestralmente ha compartido cultura, costumbres y comportamientos sociales y económicos con el Caribe. El carácter fronterizo se debe a que 9 de sus 25 municipios tienen límites con la fronte de Venezuela, lo cual hace que participe de las dinámicas entre Colombia y Venezuela. Debido a lo mencionado anteriormente, el departamento de Cesar resulta importante por su diversidad social, étnica, cultura y natural (DNP, 2011).

Debido a su topografía predominantemente plana, el comportamiento anual de la precipitación en el departamento presenta una reducida variabilida espacial. El régimen de lluvias es de tipo bimodal en la mayor parte del departamento y, en cuanto a la temperatura, Cesar constituye una de las zonas más calientes del país, con promedios superiores a 28°C. El clima cálido semiárido se encuentra principalmente en el norte del departamento, mientras que el clima cálido semihúmedo ocupa sectores del centro y del sur (IDEAM, 2020).

3. DESCRIPCIÓN DE DATOS Y MÉTODOS

Los datos corresponden al pH en agua de los primeros 5 cm del suelo y fueron obtenidos de SoilGrids, un sistema de cartografía automatizada del suelo basado en métodos de predicción espacial de última generación.

Los métodos de interpolación seleccionados para la evaluación del pH en el departamento de Cesar fueron Polígonos de Thiessen, Distancia Inversa Ponderada (IDW) y Kriging.

El método de los Polígonos de Thiessen consiste en delimitar áreas de influencia (unidades discretas) a partir de un conjunto de puntos, por lo que el tamaño y la configuración de los polígonos depende de la distribución de los puntos originales. Una limitante de este método es que no se puede estimar el error asociado, pues el valor para cada polígono se obtiene a partir de un solo punto. Inicialmente se forman triángulos entre los puntos muestreados más cercanos y se unen con segmentos rectos sin que se corten entre sí. A partir de allí se trazan líneas perpendiculares a todos los lados del triángulo, las que al unirse en un punto común dentro de cada triángulo conforman una serie de polígonos que delimitan el área de influencia de cada punto (Linsley, 1982). La fórmula general puede expresarse como:

Donde P es el valor estimado, Pi es el valor muestreado, Ai es el área del polígono correspondiente al sitio i, A es el área de la zona de estudio y n es el número de datos.

El método de distancia inversa ponderada (IDW), es un modelo determinístico que se basa en el supuesto de que el valor de una propiedad del suelo en un sitio no muestreado es el valor medio ponderado por la distancia de los valores de esa propiedad del suelo que se encuentran alrededor, dentro de un radio establecido (Li, 2011).En este método se le asigna una ponderación más alta a los puntos más cercanos a la posición por estimar, por lo que, a medida que la distancia se incrementa la ponderación asociada con el valor muestreado decrece (Toro y Melo, 2009)

Como fórmula general se tiene:

donde Zi son los valores muestreados. Di es la distancia de cada Zi al punto por estimar, P es la potencia, n es el número de ubicaciones con valores conocidos dentro de un radio establecido y Zp es el valor estimado.

Por otro lado, el método de Kriging es el método de interpolación geoestadística más utilizado, el cual incorpora un semivariograma como herramienta para el modelado de la dependencia espacial entre los valores muestreados y la distancia entre ellos (Radojac, et al., 2020). Para emplear este método se supone que la variable tiene distribución normal, la estructura global es constante y la covarianza entre los valores de los residuales depende únicamente de la distancia (Castro, et al., 2017). Para la construcción del semivariograma se usa la fórmula:

donde γ(h) es la semivarianzam, N es el número de parejas de muestras, Z (xi) es el punto muestreado y Z(xi+h) es el punto por estimar, el cual está a una distancia h de Z(xi).

4. PRESENTACIÓN DE RESULTADOS

Inicialmente se debe acceder a los datos de SoilGrids usando el protocolo webDAV, además, se usan las funciones **gdalUTils* para descargar los datos de ISRIC. ISRIC es el Centro Internacional de Información y Referencia de Suelos, su ojetivo es aumentar la conciencia y comprensión de los suelos en las principales problemáticas mundiales. Por otro lado, SoilGrids es un sistema para el mapeo de suelos digital y global que hace uso de la información del perfil del suelo y datos de covariables para modelar la distribución espacial de las propiedades del suelo en todo el mundo.

pH
## class      : RasterLayer 
## dimensions : 58034, 159246, 9241682364  (nrow, ncol, ncell)
## resolution : 250, 250  (x, y)
## extent     : -19949750, 19861750, -6147500, 8361000  (xmin, xmax, ymin, ymax)
## crs        : +proj=igh +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
## source     : https://files.isric.org/soilgrids/latest/data/phh2o/phh2o_0-5cm_mean.vrt 
## names      : phh2o_0.5cm_mean 
## values     : -32768, 32767  (min, max)

La salida anterior nos muestra el tamaño de esta capa, la cual tiene 58.034 filas con 159.246 columnas, por lo que hay más de 9.204 millones de pixeles.

Ahora, es posible filtrar los datos y utilizar únicamente los de Cesar

Cesar_igh
hist(Cesar_pH_wgs84)

El histograma obtenido muestra claramente que los suelos que predominan en el departamento de Cesar son ligeramente ácidos con un pH de 5 a 5.5, lo cual podría ser un impedimento para algunas plantas así como para la absorción de micronutrientes como el hierro, manganeso, boro, cobre y zinc, lo que conllevaría a posibles deficiencias en los cultivos.

plot(pH.mask, main= "pH de los suelos del departamento de Cesar")
plot(Aoi, add=TRUE)

El mapa corrobora el análisis realizado a partir del histograma, pero además permite observar la variación espacial del pH, obteniendo así que en la parte sur del departamento predominan los suelos con pH de 5 a 5.5, en menor proporción de 4.5 a 5; mientras que, en la zona norte del Cesar se encuentran principalmente suelos con un pH de 5.5 a 6 y en menor proporció de 6.0 a 6.5. Además, añadir las divisiones municipales permite un mayor análisis de la variable puesto que facilita entender su distribución a lo largo del territorio departamental y permite un análisis más exacto del comportamiento del pH.

Por otro lado, varíar la escala en la que se representan los valores da una idea más detallada de la distribución del pH ya que se tienen en cuenta variaciones más pequeñas que, como se mencionó anteriormente, pueden afectar significativamente el crecimiento y desarrollo de las plantas.

4.1 INTERPOLACIÓN CON POLÍGONOS DE THIESSEN

Para comenzar con la interpolación es necesario convertir el raster a puntos y disminuir la cantidad de pixeles

El mapa muestra la distribución espacial del pH en los suelos de Cesar luego de haber disminuido la cantidad de datos.

## class       : SpatialPointsDataFrame 
## features    : 232 
## extent      : 999430.4, 1128974, 1346222, 1683996  (xmin, xmax, ymin, ymax)
## crs         : +proj=tmerc +lat_0=4.59620041666667 +lon_0=-74.0775079166667 +k=1 +x_0=1000000 +y_0=1000000 +ellps=GRS80 +units=m +no_defs 
## variables   : 3
## names       :      MUNIC, CODIGO,               pH 
## min values  :  AGUACHICA,  20001, 5.12337388243505 
## max values  : VALLEDUPAR,  20787, 6.12426346668117

La salida anterior permite observar que el pH mínimo en el departamento se encuentra en el municipio de Aguachica y es de 5.12, mientras que el máximo valor está en Valledupar con 6.12.

El siguiente mapa ha convertido el raster en puntos, observando así los lugares donde se realizaron las mediciones

tm_shape(Cesar2) + tm_polygons() +
  tm_shape(pH2) +
  tm_dots(col="pH", palette = "BrBG", title="Valores de pH en los primeros 5 cm de suelo", size=0.1) +
  tm_text("pH", just="center", xmod=.6, size = 0) +
  tm_legend(legend.outside=TRUE)
## Warning in sp::proj4string(obj): CRS object has comment, which is lost in output

##
th.z     <- over(th, pH2, fn=mean)
th.spdf  <-  SpatialPolygonsDataFrame(th, th.z)
th.clp   <- raster::intersect(Cesar2,th.spdf)
###

El mapa permite observar la superficie generada mediante el método de interpolación de polígonos de Thiessen.

tm_shape(th.clp) + 
  tm_polygons(col="pH", palette="BrBG",
              title="Interpolación del pH del suelo usando Polígonos de Thiessen") +
  tm_legend(legend.outside=TRUE)
## Warning in sp::proj4string(obj): CRS object has comment, which is lost in output

###

4.2 INTERPOLACIÓN CON IDW

El mapa permite observar la superficie continua de los datos de pH en el departamento de Cesar a partir de la interpolación con IDW

##
tm_shape(r.m) + 
  tm_raster(n=10,palette = "BrBG", auto.palette.mapping = FALSE,
            title="Interpolación con IDW del pH") + 
  tm_shape(pH2) + tm_dots(size=0.09) +
  tm_legend(legend.outside=TRUE)
## Warning: The argument auto.palette.mapping is deprecated. Please use midpoint
## for numeric data and stretch.palette for categorical data to control the palette
## mapping.

###

La siguiente gráfica presenta la relación entre los datos de pH y los datos estimados usando la interpolación

OP <- par(pty="s", mar=c(4,3,0,0))
  plot(IDW.out ~ P$pH, asp=1, xlab="Observed", ylab="Predicted", pch=16,
       col=rgb(0,0,0,0.5))
  abline(lm(IDW.out ~ P$pH), col="red", lw=2,lty=2)
  abline(0,1)

###

Y también es posible calcular el error medio cuadrático:

## [1] 0.1105045

Finalmente, se muestra el mapa de validación cruzada, el cual mide la incertidumbre de la interpolación, con un intervalo de confianza del 95% y una potencia de 2 (idp=2.0)

4.2.1 Ajuste polinomial de primer grado

También es posible realizar un ajuste con un polinomio de primer grado

4.2.2 Ajuste polinomial de segundo grado

4.3 INTERPOLACIÓN CON KRIGING

Primero se deben obtener los parámetros del variograma y luego con estos es posible hacer la representación gráfica

A continuación se muestra la superficie generada con el método de Kriging

Y a continuación el mapa con un intervalo de confianza del 95% generado a partir del mapa de varianza

5. ANÁLISIS DE RESULTADOS

Los tres métodos de interpolación permiten un acercamiento diferente a los datos reales debido a que algunos se asemejan más a los datos muestreados. Inicialmente, el método de polígonos presenta a grandes rasgos la distribución de la variable, por lo que podría servir como primer acercamiento a estos; sin embargo, debido a que el muestreo se realiza en forma de cuadrícula, los polígonos que se forman son cuadrados que no representan detalladamente la variación espacial del pH. Otra desventaja de este método es que no es posible estimar el error asociado, debido a que cada punto se obtiene a partir de un único punto.

Por otra parte, el método de interpolación IDW tiene una baja incertidumbre, ya que se obtuvo un error medio cuadrático de 0.11, es decir, los datos estimados pueden estar 0.11 unidades por encima o por debajo del valor real, lo cual no representaría un error significativo y podría hacer considerar este método como el más exacto. Asimismo, al analizar el mapa del intervalo de confianza de 95% se obtuvo que los valores del error no son mayores a 0.015, sino que en la mayoría del departamento se encuentran entre 0.005 a 0.010.

El método de Kriging, aunque se esperaba que generara una superficie continua muy similar a los datos inciales no fue así, esto pudo deberse a la disminución de la cantidad de datos. Por el contrario, se observan fácilmente 8 áreas de diferentes valores de pH y, por tanto, de diferentes colores, sin tener en cuenta la variación dentro de cada área que se observa en el primer mapa de los datos originales. Lo anterior se corrobora al analizar el mapa del intervalo de confianza del 95%, en el que se muestra que los errores obtenidos con este método de interpolación se encuentran entre 0.10 y 0.15 a lo largo de todo el territorio, incrementándose en las zonas periféricas.

Este último método se esperaba que fuera el más preciso debido a que realiza los cálculos en función de la variabilidad espacial o autocorrelación espacial, asegurando la mínima varianza, y no solamente tiene en cuenta los valores vecinos al momento de realizar la estimación, tal como ocurre en IDW. Sin embargo, como ya se mencionó, esto pudo verse afectado por la disminución de la cantidad de datos.

6. CONCLUSIONES

Los valores de pH más altos se encuentran en la parte norte del departamento, lo que corresponde a los municipios de Valledupar, Pueblo Bello, el Copey y Bosconia; mientras que, los de pH más bajo corresponden a San Alberto, San Martín y Aguachica. Los valores altos de pH en Valledupar pueden deberse a prácticas de encalado, puesto que este municipio tiene el mayor rendimiento en producción de hortalizas.

La interpolación tiene aplicaciones en diferentes ámbitos académicos. En la agronomía resulta de gran importancia debido a que permite generar superficies continuas a partir de datos discretos, lo cual contribuye a un mejor entendimiento de la variabilidad espacial de la propiedad o medida de interes. Lo anterior contribuye a mejorar las practicas de manejo implementadas en los sitemas productivos, lo que a su vez, podría reflejarse como un aumento de la producción y disminución de los gastos.

Los métodos de polígonos y de IDW generaron una superficie continua más similiar a la real, mientras que, el método de Kriging tuvo inconvenientes al interpolar valores medios de pH (5.5-5.7). Además, debido a que, al realizar la interpolación con IDW es posible calcular el error medio cuadrático, el cual fue muy bajo, podría concluirse que este fue el método más exacto para estimar los valores de pH en los suelos de Cesar.

Finalmente, para mejorar las interpolaciones realizadas con cada método se debería usar un área de estudio más grande, de manera que el error no aumente en las zonas periféricas. Además, para el caso del método de polígonos, se debería usar un muestreo irregular que permita una estimación de los valores de pH más exacta. Igualmente, se podrían implementar técnicas de interpolación como Kriging Universal, el cual asume que la estructura global no es constante a lo largo del espacio, mientras que Kriging ordinario sí lo asume.

7. REFERENCIAS

Aguilar, E., Aldana, M., García, M.(2018). El pH del suelo. Conceptos Fundamentales. West Analítica y Servicios S.A de C.V https://westanalitica.com.mx/wp-content/uploads/2018/08/EL-pH-DEL-SUELO.-CONCEPTOS-FUNDAMENTALES.pdf

Arias, A., Acevedo, T. (2009). Monografía Político Electoral Departamento de Cesar. 1977 a 2007.

Barbaro, L., Karlanian, M., Mata, D. (2019). Importancia del pH y la Conductividad Eléctrica (CE) en los sutratos para plantas. Ministerio de Agricultura, Ganadería y Pesca de Argentina. https://inta.gob.ar/sites/default/files/script-tmp-inta_-_importancia_del_ph_y_la_conductividad_elctrica.pdf

Castro, M., García, D., & Jiménez, A. (2017). Comparación de técnicas de interpolación espacial de propiedades del suelo en el piedemonte llanero colombiano. Tecnura, 21(53), 78-95.

DANE. (2005). Estadísticas del departamento de Cesar.

DNP. (2011). Visión de Desarrollo Territorial Departamental. Visión CESAR Caribe 2032: Un departamento en crecimiento generando bienestar. https://colaboracion.dnp.gov.co/CDT/Desarrollo%20Territorial/VISION%20CESAR.pdf

Garbanzo, G., Alemán, B., Alvarado, A., & Henríquez, C. (2017). Validación de modelos geoestadísticos y convencionales en la determinación de la variación espacial de la fertilidad de suelos del Pacífico Sur de Costa Rica. Investigaciones Geográficas, Boletín del Instituto de Geografía, 2017(93), 20-41. https://www-sciencedirect-com.ezproxy.unal.edu.co/science/article/pii/S0188461117300481#bib0310

Henríquez, C., Méndez, J. C., & Masís, R. (2013). Interpolación de variables de fertilidad de suelo mediante el Análisis kriging y su validación¹. Agronomía costarricense, 37(2), 71-82. https://www.scielo.sa.cr/scielo.php?script=sci_arttext&pid=S0377-94242013000200006

IDEAM. 2020. Cesar. Atlas IDEAM. http://atlas.ideam.gov.co/basefiles/cesar_texto.pdf

Li, J., y Heap, A. (2011). A Review of Comparative Studies Of Spatial Interpolation Methods in Environmental Sciences: Performance and Impact Factors. Ecological Informatics, 6(3-4), 228-241

Linsley, K. (1982); “Hidrología para Ingenieros”, Editorial McGraw-hill, Bogotá (Colombia).

Liu, T., Wu, X., Li, H., Alharbi, H., Wang, J., Dang, P.& Yan, W. (2020). Soil organic matter, nitrogen and pH driven change in bacterial community following forest conversion. Forest Ecology and Management, 477, 118473. https://www-sciencedirect-com.ezproxy.unal.edu.co/science/article/pii/S0378112720312421

Oliver, M. A. (2010). The variogram and kriging. In Handbook of Applied Spatial Analysis (pp. 319-352). Springer, Berlin, Heidelberg.

Pose, M., González, A., & Valcárcel, M. (2001). Dependencia espacial de datos topográficos a escala de ladera y pequeña cuenca agrícola.

Radočaj, D., Jurišić, M., Zebec, V., & Plaščak, I. (2020). Delineation of Soil Texture Suitability Zones for Soybean Cultivation: A Case Study in Continental Croatia. Agronomy, 10(6), 823.

Rivera, E., Sánchez, M., & Domínguez, H. (2018). pH como factor de crecimiento en plantas. Revista de Iniciación Científica, 4, 101-105. https://core.ac.uk/download/pdf/234019718.pdf

Toro, G., y Melo, C. (2009). Aplicación de métodos de interpolación geoestadísticos para la predicción de niveles digitales de una imagen satelital con líneas perdidas y efecto sal y pimienta. Revista Tecnura, 12(24), 55-67.

URosario. (2018). CESAR. https://repository.urosario.edu.co/bitstream/handle/10336/8681/HurtadoEscobar-PaulaAndrea-2014.pdf?sequence=18&isAllowed=y