Geostatistical Analysis

An exploratory data analysis was initially performed on the values obtained for organic matter, beginning with a univariate description. Measures of location, dispersion and shape were calculated, as well as the histogram and the normality curve. Structural analysis (variogram) of the data was performed to investigate if the values showed a spatial structure that would allow the use of the kriging interpolation technique, considered the best unbiased linear estimator (Cressie, 1990CRESSIE, N.: “The origins of kriging”, Mathematical geology, 22: 239-252, 1990, ISSN: 0882-8121, Publisher: Springer, DOI: https://doi.org/10.1007/978-3-319-78999-6.). Interpolation was performed by kriging to obtain maps of organic matter content at unmeasured points. Surfer 8 software (Golden Software, Inc.) was used.

The experimental semivariogram γ (h) was obtained from expression (1) $\gamma \left(h\right)=\frac{1}{2Np\left(h\right)}\sum _{i=1}^{Np\left(h\right)}{\left[Z\left(x_i\right)-Z\left(x_i+h\right)\right]}^2$ (Journel & Huijbregts, 1978JOURNEL, A.G.; HUIJBREGTS, C.J.: Mining geostatistics, Ed. Academic Press, London, England, 600 p., 1978.).

Where Np(h) is the number of pairs of observations separated by distance h, Z(xi) is the value of the variable at site xi and Z(xi + h) is the value of the variable at a site located at a distance h from site xi.

The adjustment of the experimental semivariogram to theoretical models was carried out, obtaining the one with the best adjustment according to the methodology proposed by Legrá-Lobaina & Atanes-Beatón (2010)LEGRÁ-LOBAINA, A.A.; ATANES-BEATÓN, D.M.: “Variogramas adaptativos: un método práctico para aumentar la utilidad del error de estimación por kriging”, Minería y Geología, 26(4): 53-78, 2010, ISSN: 1993-8012, Publisher: Instituto Superior Minero Metalúrgico de Moa’Dr Antonio Nuñez Jiménez’.. In this, the semivariogram is generated from the adaptation of its scope, plateau, nugget effect and its model, which allows obtaining local corrections of the variability of the magnitude under study.

Once the theoretical model was established, the values of the nugget effect (Co) were found, which is the lowest value of the semivariance and the maximum semivariance (C₀ + C₁). On the other hand, C₁ is the difference between the maximum semivariance and the value of the pip effect. In order to quickly obtain quantitative information on the spatial dependence of the OM variable, the Degree of Spatial Dependence (GDE) proposed by Cambardella et al. (1994)CAMBARDELLA, C.A.; MOORMAN, T.; NOVAK, J.; PARKIN, T.; KARLEN, D.; TURCO, R.; KONOPKA, A.: “Field‐scale variability of soil properties in central Iowa soils”, Soil science society of America journal, 58(5): 1501-1511, 1994, ISSN: 0361-5995, Publisher: Wiley Online Library, DOI: https://doi.org/10.2136/sssaj1994.03615995005800050033x., was also calculated, which is defined by expression (2) $GDE\left(\%\right)=\left(\frac{C_o}{C_0+C_1}\right)*100$ .

Seidel & Oliveira (2014)SEIDEL, E.J.; OLIVEIRA, M.S.: “Novo índice geoestatístico para a mensuração da dependência espacial”, Revista Brasileira de Ciência do Solo, 38: 699-705, 2014, ISSN: 1806-9657, Publisher: SciELO Brasil, DOI: https://doi.org/10.1590/S0100-06832014000300002 propose the following categories for the GDE: strong spatial dependence (GDE > 75 %), moderate spatial dependence (25 <GDE ≤ 75 %) and weak spatial dependence (GDE ≤ 25 %). To also take into account the effect of the model used to fit the experimental variogram as well as all the other characteristics of the semivariogram, the spatial dependence index of the model (IDE) proposed by Seidel & Oliveira (2014)SEIDEL, E.J.; OLIVEIRA, M.S.: “Novo índice geoestatístico para a mensuração da dependência espacial”, Revista Brasileira de Ciência do Solo, 38: 699-705, 2014, ISSN: 1806-9657, Publisher: SciELO Brasil, DOI: https://doi.org/10.1590/S0100-06832014000300002 and Seidel & Oliveira (2016)SEIDEL, E.J.; OLIVEIRA, M.S.: “A classification for a geostatistical index of spatial dependence”, Revista Brasileira de Ciência do Solo, 40, 2016, ISSN: 1806-9657, Publisher: SciELO Brasil., which is given by the following expressions for the spherical, exponential and Gaussian models, respectively.

Where the practical range is a. MD is the maximum distance. The coefficients that appear at the beginning of each model according to Seidel & Oliveira (2014)SEIDEL, E.J.; OLIVEIRA, M.S.: “Novo índice geoestatístico para a mensuração da dependência espacial”, Revista Brasileira de Ciência do Solo, 38: 699-705, 2014, ISSN: 1806-9657, Publisher: SciELO Brasil, DOI: https://doi.org/10.1590/S0100-06832014000300002, are known as the Model Factor (FM) and express the strength of the special dependency that a given model can achieve, so the higher its value, the greater the strength of the spatial dependency of the model. Seidel & Oliveira (2016)SEIDEL, E.J.; OLIVEIRA, M.S.: “A classification for a geostatistical index of spatial dependence”, Revista Brasileira de Ciência do Solo, 40, 2016, ISSN: 1806-9657, Publisher: SciELO Brasil. propose, according to the IDE values, the classification presented in Table 1.

Satellite Image Processing

An image of April 26, 2022, belonging to the Landsat 9 OLI/TIRS 2 satellite (LC09_L2SP_011046_20220426_20220428_02_T1) of the United States Geological Survey on path 011 row 046 was used and was projected on the WGS 84 System UTM Zone 18 North in the QGIS 3.10 “A Coruña” software and spectral index of soil and vegetation were determined (Table 2), after performing the atmospheric correction to eliminate the effect of clouds on the image.

For the extraction of the digital values of the image, the layer of sampling points was used in the ArcGIS 10.5 software and an Excel database was created with said information for each calculated spectral index. In the Statgraphics Plus software, the exploratory and linear regression analysis was carried out between the values of OM and the spectral index of soil and vegetation.

Exploratory Analysis of the Variables under Study

Table 3 shows the statistics of the variables analyzed. The average OM content was found to be 3.81 % with minimum and maximum values of 1.65 % and 6.75 %, respectively. That could be associated with the fact that the area has remained fallow for three years and the possible presence of grazing animals at some points, which leads to the incorporation of OM into the soil by the decomposition of their excreta.

The median showed a trend of 3.74% with a standard deviation of 1.25%, with a standard error in its determination of 0.12% in the permissible ranges in which the unit of measurement of this property oscillates. The coefficient of variation indicated that the values of OM vary moderately for 32.80% according to Wilding (1985)WILDING, L.: “Spatial variability: its documentation, accomodation and implication to soil surveys”, En: Soil spatial variability, Las Vegas NV, 30 November-1 December 1984, Netherlands, pp. 166-194, 1985, ISBN: 90-220-0891-6.. Alexakis et al. (2019)ALEXAKIS, D.D.; TAPOGLOU, E.; VOZINAKI, A.E.K.; TSANIS, I.K.: “Integrated use of satellite remote sensing, artificial neural networks, field spectroscopy, and GIS in estimating crucial soil parameters in terms of soil erosion”, Remote Sensing, 11(9): 1106, 2019, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs11091106. refer that the coefficient of variation reflects the distribution of each soil property and can have characteristic spatial patterns for each experimental area. Ayoubi et al. (2011)AYOUBI, S.; SHAHRI, A.; KARCHEGANI, P.A.; SAHRAWAT, K.L.: Application of artificial neural network (ANN) to predict soil organic matter using remote sensing data in two ecosystems, Ed. InTech Rijeka, Croatia, vol. Biomass and Remote Sensing of Biomass, 181-196 p., 2011, ISBN: 978-953-307-490-0. obtained a moderate variation with a coefficient of variation of 32 % and 34 % in the sites where they sampled the OM content and refer that this variation depends on the accumulation of water on the soil cover.

The asymmetry and kurtosis, both in the OM content of the soil and in the determined spectral index, is in the range of -1 to 1 which indicates that the values do not follow a normal distribution (López-Granados et al., 2005LÓPEZ-GRANADOS, F.; JURADO-EXPÓSITO, M.; PEÑA-BARRAGÁN, J.M.; GARCÍA-TORRES, L.: “Using geostatistical and remote sensing approaches for mapping soil properties”, European Journal of Agronomy, 23(3): 279-289, 2005, ISSN: 1161-0301, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.eja.2004.12.003.).

The average value of NDVI was found to be 0.26, which ranges from -1 to 1 and agrees with what was stated by Rawashdeh (2012)RAWASHDEH, A.S.B.: “Assessment of change detection method based on normalized vegetation index in environmental studies”, International Journal of Applied Science and Engineering, 10(2): 89-97, 2012, ISSN: 1727-7841, Publisher: Chaoyang University of Technology. that for this index, values from 0 to 0.5, there is scarce vegetation and coincides with the current state of the study area. The SAVI vegetation index reported an average of 0.52, closely related to the average index obtained from NDVI and the corresponding vegetation status classification.

Joko-Prasetyo et al. (2020)JOKO-PRASETYO, S.Y.; DWI, K.; CHRISMAWATI-PASELENG, M.; WIDIYANTO, D.C.; WINARKO, E.: “Satellite imagery and machine learning for aridity disaster classification using vegetation indices”, Bulletin of Electrical Engineering and Informatics, 9(3): 1149-1158, 2020, ISSN: 2302-9285, DOI: https://doi.org/10.11591/eei.v9i3.1916. use the NDVI and SAVI as indicators of the state of aridity in Indonesia, obtaining values similar to those found in this research and pointing out that they are areas with low vegetation cover where photosynthetic activity is decreasing and, therefore, the values fluctuate between 0.1 and 0.5.

In line with the type of soil in the study area (Vertisol), characterized by having a high content of monmorillonite-type clays (Hernández et al., 2015HERNÁNDEZ, J.A.; PÉREZ, J.J.M.; MESA, N.A.; BOSCH, I.D.; RIVERO, L.; CAMACHO, E.: Nueva versión de la clasificación genética de los suelos de Cuba., Ed. AGRINFOR, La Habana, Cuba, ISBN: 959-246-022-1, Barcaz L L ed., vol. I, La Habana, Cuba, 64 p., 2015, ISBN: 959-246-022-1., 2019HERNÁNDEZ-JIMÉNEZ, A.; PÉREZ-JIMÉNEZ, J.M.; BOSCH-INFANTE, D.; SPECK-CASTRO, N.: “La clasificación de suelos de Cuba: énfasis en la versión de 2015”, Cultivos Tropicales, 40(1), 2019, ISSN: 0258-5936, Publisher: Ediciones INCA.), the CI yielded mean values greater than 1 with 1.32. The results obtained coincide with what was stated by Boettinger et al. (2008)BOETTINGER, J.; RAMSEY, R.; BODILY, J.; COLE, N.; KIENAST-BROWN, S.; NIELD, S.; SAUNDERS, A.; STUM, A.: Landsat spectral data for digital soil mapping, Ed. Digital soil mapping with limited data Australia: Springer Scienc, A.E. Hartemink, A.B. McBratney, M.L. Mendonca-Santo ed., Australia, publisher: Springer, 2008, ISBN: 978-1-4020-8591-8. who report that the multispectral images of the Landsat satellite in its near-infrared bands can be used to identify the parent material of the soil.

Geostatistical Analysis

Figure 2 illustrates the experimental and theoretical semivariogram of the values obtained from OM in the study area, which had a better fit to an exponential model, agreeing with studies previously carried out by Reza et al. (2016)REZA, S.; NAYAK, D.; CHATTOPADHYAY, T.; MUKHOPADHYAY, S.; SINGH, S.; SRINIVASAN, R.: “Spatial distribution of soil physical properties of alluvial soils: a geostatistical approach”, Archives of agronomy and soil science, 62(7): 972-981, 2016, ISSN: 0365-0340, Publisher: Taylor & Francis, DOI: https://doi.org/10.1080/03650340.2015.1107678.; Bogunovic et al. (2017)BOGUNOVIC, I.; KISIC, I.; MESIC, M.; PERCIN, A.; CAÑIZARES, Z.J.; BILANDŽIJA, D.; JONJIC, A.; PEREIRA, P.: “Reducing sampling intensity in order to investigate spatial variability of soil pH, organic matter and available phosphorus using co-kriging techniques. A case study of acid soils in Eastern Croatia”, Archives of Agronomy and Soil Science, 63(13): 1852-1863, 2017, ISSN: 0365-0340, ublisher: Taylor & Francis, DOI: https://doi.org/10.1080/03650340.2017.1311013.; Durdevic et al. (2019)DURDEVIC, B.; JUG, I.; JUG, D.; BOGUNOVIC, I.; VUKADINOVIC, V.; STIPESEVIC, B.; BROZOVIC, B.: “Spatial variability of soil organic matter content in Eastern Croatia assessed using different interpolation methods”, International Agrophysics, 33(1), 2019, ISSN: 0236-8722, Publisher: Polska Akademia Nauk. Instytut Agrofizyki PAN, DOI: https://doi.org/10.31545/intagr/104372.. They report that most of the soil properties, when performing a structural analysis, have a better fit to an exponential model.

Jian-Bing et al. (2006)JIAN-BING, W.; DU-NING, X.; XING-YI, Z.; XIU-ZHEN, L.; XIAO-YU, L.: “Spatial variability of soil organic carbon in relation to environmental factors of a typical small watershed in the black soil region, northeast China”, Environmental monitoring and assessment, 121: 597-613, 2006, ISSN: 0167-6369, Publisher: Springer, DOI: https://doi.org/10.1007/s10661-005-9158-5. and Rawashdeh (2012)RAWASHDEH, A.S.B.: “Assessment of change detection method based on normalized vegetation index in environmental studies”, International Journal of Applied Science and Engineering, 10(2): 89-97, 2012, ISSN: 1727-7841, Publisher: Chaoyang University of Technology. when studying the spatial variability of chemical properties of the soil, obtained similar results. They refer that OM, pH, electrical conductivity, assimilable potassium and total carbonate had a better fit to an exponential model. In Figure 2, according to the obtained range of 600 m, it can be established that in samples taken at distances less than this one, their values will be spatially related, while those taken at greater distances are not related, due to the fact that the semivariance is made equal to the sample variance (Kerry & Oliver, 2007KERRY, R.; OLIVER, M.: “Comparing sampling needs for variograms of soil properties computed by the method of moments and residual maximum likelihood”, Geoderma, 140(4): 383-396, 2007, ISSN: 0016-7061, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.geoderma.2007.04.019.).

The degree of spatial dependence (GDE) was 43.75%, which corresponds to the classification proposed previously by Cambardella et al. (1994)CAMBARDELLA, C.A.; MOORMAN, T.; NOVAK, J.; PARKIN, T.; KARLEN, D.; TURCO, R.; KONOPKA, A.: “Field‐scale variability of soil properties in central Iowa soils”, Soil science society of America journal, 58(5): 1501-1511, 1994, ISSN: 0361-5995, Publisher: Wiley Online Library, DOI: https://doi.org/10.2136/sssaj1994.03615995005800050033x. as moderate spatial dependence (25 <GDE ≤ 75%) and according to the classification of the effect of the adjusted exponential model, an SDI value of 16.36 % implies a strong spatial dependence (Table 1) as it is greater than 13 %. These values reflect that the spatial dependence is controlled by intrinsic and extrinsic factors influenced by inappropriate agricultural practices in the soil (Liu et al., 2014LIU, Z.; ZHOU, W.; SHEN, J.; HE, P.; LEI, Q.; LIANG, G.: “A simple assessment on spatial variability of rice yield and selected soil chemical properties of paddy fields in South China”, Geoderma, 235: 39-47, 2014, ISSN: 0016-7061, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.geoderma.2014.06.027.). The existing spatial structure allows the use of kriging as an interpolation technique that will allow the preparation of the OM map by estimating values at unmeasured points.

Figure 3 shows the distribution of the OM content obtained by kriging, where the highest values of OM are found to the north with darker tones in an irregular transept in the study area while the lowest values are in the center with lighter shades. It is possible that these higher OM values in this area are associated with the presence of peasant cattle as the area is fallow.

Relationships between Spectral Index and Organic Matter Content

Table 4 presents the statistics of the linear regression analysis between the OM content of the soil and the spectral index of soil and vegetation. There is a high correlation of 0.98 between the OM content of the soil, the NDVI and the SAVI, which may be because both indexes use the red and infrared bands of the sensor for their determination.

The CI had a correlation of 0.94 with respect to the OM content of the soil. The determination coefficient yielded values close to 100%, the NDVI model vs. O.M (%) 95.61%, SAVI vs. O.M (%) 95.66% and CLAY INDEX vs. O.M (%) of 88.92 %. Therefore, it can be affirmed that the spectral index used can predict the OM content, with an error in its forecast in all cases within the permissible ranges in which the determined variables are measured (Ayoubi et al., 2011AYOUBI, S.; SHAHRI, A.; KARCHEGANI, P.A.; SAHRAWAT, K.L.: Application of artificial neural network (ANN) to predict soil organic matter using remote sensing data in two ecosystems, Ed. InTech Rijeka, Croatia, vol. Biomass and Remote Sensing of Biomass, 181-196 p., 2011, ISBN: 978-953-307-490-0.).

Correlation coefficient close to 1.0; of 0.86 and 0.90 were found by Wang et al. (2018)WANG, X.; ZHANG, F.; JOHNSON, V.C.: “New methods for improving the remote sensing estimation of soil organic matter content (SOMC) in the Ebinur Lake Wetland National Nature Reserve (ELWNNR) in northwest China”, Remote Sensing of Environment, 218: 104-118, 2018, ISSN: 0034-4257, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.rse.2018.09.020. with the use of Landsat 8 OLI/TIRS images with the processing of the red and infrared bands of the sensor, they obtained a strong determination of 92 % in the relationship between the OM content of the soil and the combination of sensor bands.

Xu et al. (2023)XU, X.; DU, C.; MA, F.; QIU, Z.; ZHOU, J.: “A Framework for High-Resolution Mapping of Soil Organic Matter (SOM) by the Integration of Fourier Mid-Infrared Attenuation Total Reflectance Spectroscopy (FTIR-ATR), Sentinel-2 Images, and DEM Derivatives”, Remote Sensing, 15(4): 1072, 2023, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs15041072. report that the use of remote sensing establishes a strong relationship between the data captured by the sensor and the OM of the soil with a strong linear relationship between the NDVI variables and index obtained by a digital elevation model (DEM), while no positive relationship was found with the Sentinel 2 spectral bands. Previous studies have shown that the relationship of soil OM content through remote sensing cannot be seen as a methodology that is generalized to various environments, but rather which is unique for each study site, and depends on the type of sensor used, the characteristics of the soil, the relief and the climate (Lamichhane et al., 2019LAMICHHANE, S.; KUMAR, L.; WILSON, B.: “Digital soil mapping algorithms and covariates for soil organic carbon mapping and their implications: A review”, Geoderma, 352: 395-413, 2019, ISSN: 0016-7061, V, DOI: https://doi.org/10.1016/j.geoderma.2019.05.031.).

On the other hand, Prudnikova & Savin (2021)PRUDNIKOVA, E.; SAVIN, I.: “Some peculiarities of arable soil organic matter detection using optical remote sensing data”, Remote Sensing, 13(12): 2313, 2021, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs13122313. found a negative relationship between the CI and the OM content of the soil with a coefficient of determination for the rainy season of 81 % and of 84 % in the dry season, respectively when using Sentinel 2 to estimate the OM content of an Albic Luvisols in Russia.

There are references to other investigations in which estimates of the organic carbon content of the soil are made (Sodango et al., 2021SODANGO, T.H.; SHA, J.; LI, X.; NOSZCZYK, T.; SHANG, J.; ANESEYEE, A.B.; CHAFIK, Z.: “Modeling the spatial dynamics of soil organic carbon using remotely-sensed predictors in Fuzhou city, China”, Remote Sensing, 13(9): 1682, 2021, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs1309168.). This chemical element derived from the OM content of the soil (Rasel et al., 2017RASEL, S.; GROEN, T.A.; HUSSIN, Y.A.; DITI, I.J.: “Proxies for soil organic carbon derived from remote sensing”, International journal of applied earth observation and geoinformation, 59: 157-166, 2017, ISSN: 1569-8432, Publisher: Elsevier, DOI: http://dx.doi.org/10.1016/j.jag.2017.03.004.) is estimated from spectral vegetation index, in which, when using the NDVI with mean values of 0.49 (dense vegetation) a correlation of 0.74 is reached (Sankar-Gouri. et al., 2019SANKAR-GOURI., G.; KUMAR-SHIT, P.; REZA-POURGHASEMI, H.: “Soil organic carbon mapping using remote sensing techniques and multivariate regression model”, Geocarto International, 34(2): 215-226, 2019, ISSN: 1010-6049, Publisher: Taylor & Francis, DOI: https://doi.org/10.1080/10106049.2017.1381179.).

Análisis geoestadístico

A los valores obtenidos de materia orgánica se les realizó inicialmente un análisis exploratorio de datos comenzando por una descripción univariada. Fueron calculadas medidas de localización, de dispersión y de forma así como el histograma y la curva de normalidad. Se realizó el análisis estructural (variograma) de los datos con vistas a investigar si los valores mostraban una estructura espacial que permitiera la utilización de la técnica de interpolación kriging considerada el mejor estimador lineal insesgado (Cressie, 1990CRESSIE, N.: “The origins of kriging”, Mathematical geology, 22: 239-252, 1990, ISSN: 0882-8121, Publisher: Springer, DOI: https://doi.org/10.1007/978-3-319-78999-6.). Se realizó la interpolación por kriging para obtener los mapas de los contenidos de materia orgánica en puntos no medidos. Se utilizó el software Surfer 8 (Golden Software, Inc.).

El semivariograma experimental γ(h) fue obtenido a partir de la expresión (1) $\gamma \left(h\right)=\frac{1}{2Np\left(h\right)}\sum _{i=1}^{Np\left(h\right)}{\left[Z\left(x_i\right)-Z\left(x_i+h\right)\right]}^2$ (Journel y Huijbregts, 1978JOURNEL, A.G.; HUIJBREGTS, C.J.: Mining geostatistics, Ed. Academic Press, London, England, 600 p., 1978.).

En la que Np(h) es el número de pares de observaciones separadas la distancia h, Z(xi) es el valor de la variable en el sitio xi y Z(xi + h) es el valor de la variable en un sitio ubicado a una distancia hdel sitio xi.

Se realizó el ajuste del semivariograma experimental a modelos teóricos, obteniendo el de mejor ajuste según la metodología propuesta por Legrá-Lobaina y Atanes-Beatón (2010)LEGRÁ-LOBAINA, A.A.; ATANES-BEATÓN, D.M.: “Variogramas adaptativos: un método práctico para aumentar la utilidad del error de estimación por kriging”, Minería y Geología, 26(4): 53-78, 2010, ISSN: 1993-8012, Publisher: Instituto Superior Minero Metalúrgico de Moa’Dr Antonio Nuñez Jiménez’., donde el semivariograma se genera a partir de la adaptación de su alcance, meseta, efecto pepita y a su modelo, lo que permite obtener correcciones locales de la variabilidad de la magnitud que se estudia.

Una vez establecido el modelo teórico se hallaron los valores del efecto pepita (Co) que es el menor valor de la semivarianza y de la máxima semivarianza (C₀ + C₁). Por otra parte C₁ es la diferencia entre la máxima semivarianza y el valor del efecto pepita. Para tener de forma rápida información cuantitativa de la dependencia espacial de la variable MO también se calculó el Grado de Dependencia Espacial (GDE) propuesto por Cambardella et al. (1994CAMBARDELLA, C.A.; MOORMAN, T.; NOVAK, J.; PARKIN, T.; KARLEN, D.; TURCO, R.; KONOPKA, A.: “Field‐scale variability of soil properties in central Iowa soils”, Soil science society of America journal, 58(5): 1501-1511, 1994, ISSN: 0361-5995, Publisher: Wiley Online Library, DOI: https://doi.org/10.2136/sssaj1994.03615995005800050033x.), que se define por la expresión (2) $GDE\left(\%\right)=\left(\frac{C_o}{C_0+C_1}\right)*100$ .

Seidel y Oliveira (2014)SEIDEL, E.J.; OLIVEIRA, M.S.: “Novo índice geoestatístico para a mensuração da dependência espacial”, Revista Brasileira de Ciência do Solo, 38: 699-705, 2014, ISSN: 1806-9657, Publisher: SciELO Brasil, DOI: https://doi.org/10.1590/S0100-06832014000300002 proponen para el GDE las siguientes categorías: fuerte dependencia espacial (GDE > 75%), moderada dependencia espacial (25 <GDE ≤ 75%) y débil dependencia espacial (GDE ≤ 25%). Para tener en cuenta también el efecto del modelo empleado para el ajuste del variograma experimental así como todos las demás características del semivariograma, se utilizó el índice de dependencia espacial del modelo (IDE) propuesto por Seidel y Oliveira (2014)SEIDEL, E.J.; OLIVEIRA, M.S.: “Novo índice geoestatístico para a mensuração da dependência espacial”, Revista Brasileira de Ciência do Solo, 38: 699-705, 2014, ISSN: 1806-9657, Publisher: SciELO Brasil, DOI: https://doi.org/10.1590/S0100-06832014000300002 y Seidel y Oliveira (2016)SEIDEL, E.J.; OLIVEIRA, M.S.: “A classification for a geostatistical index of spatial dependence”, Revista Brasileira de Ciência do Solo, 40, 2016, ISSN: 1806-9657, Publisher: SciELO Brasil., que está dado por las siguientes expresiones para los modelos esférico, exponencial y gausiano, respectivamente.

En el que a es el alcance práctico y MD es la distancia máxima. Los coeficientes que aparecen al inicio de cada modelo según Seidel y Oliveira (2014)SEIDEL, E.J.; OLIVEIRA, M.S.: “Novo índice geoestatístico para a mensuração da dependência espacial”, Revista Brasileira de Ciência do Solo, 38: 699-705, 2014, ISSN: 1806-9657, Publisher: SciELO Brasil, DOI: https://doi.org/10.1590/S0100-06832014000300002, son conocidos como Factor del Modelo (FM) y expresan la fortaleza de la dependencia espacial que un modelo dado puede alcanzar, por lo que a mayor valor del mismo, mayor será la fortaleza de la dependencia espacial del modelo.

Seidel y Oliveira (2016)SEIDEL, E.J.; OLIVEIRA, M.S.: “A classification for a geostatistical index of spatial dependence”, Revista Brasileira de Ciência do Solo, 40, 2016, ISSN: 1806-9657, Publisher: SciELO Brasil.proponen según los valores del IDE la clasificación que se presenta en la Tabla 1.

Procesamiento de la imagen satelital

Se utilizó una imagen del día 26 de abril de 2022, perteneciente al satélite Landsat 9 OLI/TIRS 2 ( LC09_L2SP_011046_20220426_20220428_02_T1 ) del Servicio Geológico de los Estados Unidos en el path 011 row 046 y fue proyectada en el Sistema WGS 84 UTM Zona 18 Norte en el software QGIS 3.10 “A Coruña”y se determinaron índices espectrales de suelo y vegetación (Tabla 2), luego de realizar la corrección atmosférica para eliminar el efecto de las nubes sobre la imagen.

Para la extracción de los valores digitales de la imagen se utilizó la capa de puntos de muestreo en el software ArcGIS 10.5 y se confeccionó una base de datos en Excel con dicha información por cada índice espectral calculado. En el software Statgraphics Plus se realizó el análisis exploratorio y de regresión lineal entre los valores de materia orgánica y de los índices espectrales de suelo y vegetación.

Análisis exploratorio de las variables bajo estudio

La Tabla 3 muestra los estadígrafos de las variables analizadas. El contenido medio de materia orgánica se encontró en un 3,81% con valores mínimos y máximos de 1,65% y 6,75% respectivamente lo cual podría asociarse a que el área ha permanecido tres años en barbecho y también pudo existir en algún momento la presencia de animales pastando lo que conlleva a la incorporación de MO al suelo por la descomposición de la excreta de los mismos.

La mediana mostró una tendencia de un 3,74% con una desviación estándar de 1,25%, con un error estándar en su determinación de 0,12% en los rangos permisibles en los que oscila la unidad de medida de esta propiedad. El coeficiente de variación indicó que los valores de materia orgánica varían de forma moderada para un 32,80% según Wilding (1985)WILDING, L.: “Spatial variability: its documentation, accomodation and implication to soil surveys”, En: Soil spatial variability, Las Vegas NV, 30 November-1 December 1984, Netherlands, pp. 166-194, 1985, ISBN: 90-220-0891-6.; Alexakis et al. (2019)ALEXAKIS, D.D.; TAPOGLOU, E.; VOZINAKI, A.E.K.; TSANIS, I.K.: “Integrated use of satellite remote sensing, artificial neural networks, field spectroscopy, and GIS in estimating crucial soil parameters in terms of soil erosion”, Remote Sensing, 11(9): 1106, 2019, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs11091106. refieren que el coeficiente de variación refleja la distribución de cada propiedad del suelo y puede llegar a tener patrones espaciales característicos para cada área experimental. Ayoubi et al. (2011)AYOUBI, S.; SHAHRI, A.; KARCHEGANI, P.A.; SAHRAWAT, K.L.: Application of artificial neural network (ANN) to predict soil organic matter using remote sensing data in two ecosystems, Ed. InTech Rijeka, Croatia, vol. Biomass and Remote Sensing of Biomass, 181-196 p., 2011, ISBN: 978-953-307-490-0. obtuvieron una moderada variación con un coeficiente de variación de 32,00% y 34,00% en los sitios donde muestrearon el contenido de materia orgánica y abordan que esta variación depende de la acumulación de agua sobre la cobertura del suelo.

La asimetría y kurtosis tanto en el contenido de materia orgánica del suelo como en los índices espectrales determinados se encuentra en el rango de -1 a 1 para los cuales indica que los valores no se ajustan a una distribución normal (López-Granados et al., 2005LÓPEZ-GRANADOS, F.; JURADO-EXPÓSITO, M.; PEÑA-BARRAGÁN, J.M.; GARCÍA-TORRES, L.: “Using geostatistical and remote sensing approaches for mapping soil properties”, European Journal of Agronomy, 23(3): 279-289, 2005, ISSN: 1161-0301, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.eja.2004.12.003.).

El valor medio de NDVI se encontrarón de 0,26 el cual oscilan entre -1 a 1 lo que concuerda con lo planteado por Rawashdeh (2012)RAWASHDEH, A.S.B.: “Assessment of change detection method based on normalized vegetation index in environmental studies”, International Journal of Applied Science and Engineering, 10(2): 89-97, 2012, ISSN: 1727-7841, Publisher: Chaoyang University of Technology. que para este índice, valores desde 0 hasta 0,5 se presenta una escasa vegetación que concuerda con estado actual del área de estudio. El índice de vegetación SAVI reportó una media de 0,52 en estrecha relación con el índice medio obtenido de NDVI y la clasificación del estado de vegetación correspondiente.

Joko-Prasetyo et al. (2020)JOKO-PRASETYO, S.Y.; DWI, K.; CHRISMAWATI-PASELENG, M.; WIDIYANTO, D.C.; WINARKO, E.: “Satellite imagery and machine learning for aridity disaster classification using vegetation indices”, Bulletin of Electrical Engineering and Informatics, 9(3): 1149-1158, 2020, ISSN: 2302-9285, DOI: https://doi.org/10.11591/eei.v9i3.1916. utilizan el NDVI y el SAVI como indicadores del estado de aridez en Indonesia obteniendo valores similares a los encontrados en esta investigación y señalan que son áreas con baja cobertura vegetal donde la actividad fotosintética se encuentra en decrecimiento, y por ello los valores fluctúan entre 0,1 hasta 0,5.

En consonancia al tipo de suelo del área de estudio (Vertisol), caracterizado por tener un alto contenido de arcillas del tipo monmorillonitas Hernández et al. (2015HERNÁNDEZ, J.A.; PÉREZ, J.J.M.; MESA, N.A.; BOSCH, I.D.; RIVERO, L.; CAMACHO, E.: Nueva versión de la clasificación genética de los suelos de Cuba., Ed. AGRINFOR, La Habana, Cuba, ISBN: 959-246-022-1, Barcaz L L ed., vol. I, La Habana, Cuba, 64 p., 2015, ISBN: 959-246-022-1., 2019)HERNÁNDEZ-JIMÉNEZ, A.; PÉREZ-JIMÉNEZ, J.M.; BOSCH-INFANTE, D.; SPECK-CASTRO, N.: “La clasificación de suelos de Cuba: énfasis en la versión de 2015”, Cultivos Tropicales, 40(1), 2019, ISSN: 0258-5936, Publisher: Ediciones INCA., el CI arrojó valores medios superiores a 1 con 1,32. Los resultados obtenidos coinciden con lo planteado por Boettinger et al. (2008)BOETTINGER, J.; RAMSEY, R.; BODILY, J.; COLE, N.; KIENAST-BROWN, S.; NIELD, S.; SAUNDERS, A.; STUM, A.: Landsat spectral data for digital soil mapping, Ed. Digital soil mapping with limited data Australia: Springer Scienc, A.E. Hartemink, A.B. McBratney, M.L. Mendonca-Santo ed., Australia, publisher: Springer, 2008, ISBN: 978-1-4020-8591-8. donde refieren que las imágenes multiespectrales del satélite Landsat en sus bandas del infrarrojo cercano puede ser utilizado para la identificación del material parental del suelo.

Análisis Geoestadístico

En la Figura 2 se ilustra el semivariograma experimental y teórico de los valores obtenidos de materia orgánica en el área de estudio, el cual tuvo un mejor ajuste a un modelo exponencial concordando con estudios realizados anteriormente por Reza et al. (2016)REZA, S.; NAYAK, D.; CHATTOPADHYAY, T.; MUKHOPADHYAY, S.; SINGH, S.; SRINIVASAN, R.: “Spatial distribution of soil physical properties of alluvial soils: a geostatistical approach”, Archives of agronomy and soil science, 62(7): 972-981, 2016, ISSN: 0365-0340, Publisher: Taylor & Francis, DOI: https://doi.org/10.1080/03650340.2015.1107678.; Bogunovic et al. (2017)BOGUNOVIC, I.; KISIC, I.; MESIC, M.; PERCIN, A.; CAÑIZARES, Z.J.; BILANDŽIJA, D.; JONJIC, A.; PEREIRA, P.: “Reducing sampling intensity in order to investigate spatial variability of soil pH, organic matter and available phosphorus using co-kriging techniques. A case study of acid soils in Eastern Croatia”, Archives of Agronomy and Soil Science, 63(13): 1852-1863, 2017, ISSN: 0365-0340, ublisher: Taylor & Francis, DOI: https://doi.org/10.1080/03650340.2017.1311013.; Durdevic et al. (2019)DURDEVIC, B.; JUG, I.; JUG, D.; BOGUNOVIC, I.; VUKADINOVIC, V.; STIPESEVIC, B.; BROZOVIC, B.: “Spatial variability of soil organic matter content in Eastern Croatia assessed using different interpolation methods”, International Agrophysics, 33(1), 2019, ISSN: 0236-8722, Publisher: Polska Akademia Nauk. Instytut Agrofizyki PAN, DOI: https://doi.org/10.31545/intagr/104372. donde refieren que la mayoría de las propiedades del suelo al realizarle un análisis estructural tienen un mejor ajuste a un modelo exponencial. Jian-Bing et al. (2006)JIAN-BING, W.; DU-NING, X.; XING-YI, Z.; XIU-ZHEN, L.; XIAO-YU, L.: “Spatial variability of soil organic carbon in relation to environmental factors of a typical small watershed in the black soil region, northeast China”, Environmental monitoring and assessment, 121: 597-613, 2006, ISSN: 0167-6369, Publisher: Springer, DOI: https://doi.org/10.1007/s10661-005-9158-5. y Rawashdeh (2012)RAWASHDEH, A.S.B.: “Assessment of change detection method based on normalized vegetation index in environmental studies”, International Journal of Applied Science and Engineering, 10(2): 89-97, 2012, ISSN: 1727-7841, Publisher: Chaoyang University of Technology. al estudiar la variabilidad espacial de propiedades químicas del suelo obtuvieron resultados similares, donde la materia orgánica, pH, conductividad eléctrica, potasio asimilable y carbonato total tuvieron mejor ajuste a un modelo exponencial,

En la Figura 2, de acuerdo al rango obtenido de 600 m, podremos establecer que muestras tomadas a distancias menores de esta distancia, sus valores estarán relacionados espacialmente, mientras que las que se toman a distancias mayores no están relacionadas, debido a que la semivarianza se hace igual a la varianza de la muestra (Kerry y Oliver, 2007KERRY, R.; OLIVER, M.: “Comparing sampling needs for variograms of soil properties computed by the method of moments and residual maximum likelihood”, Geoderma, 140(4): 383-396, 2007, ISSN: 0016-7061, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.geoderma.2007.04.019.).

El grado de dependencia espacial (GDE) fue de 43.75% lo cual se corresponde de acuerdo a la clasificación planteada anteriormente por Cambardella et al. (1994)CAMBARDELLA, C.A.; MOORMAN, T.; NOVAK, J.; PARKIN, T.; KARLEN, D.; TURCO, R.; KONOPKA, A.: “Field‐scale variability of soil properties in central Iowa soils”, Soil science society of America journal, 58(5): 1501-1511, 1994, ISSN: 0361-5995, Publisher: Wiley Online Library, DOI: https://doi.org/10.2136/sssaj1994.03615995005800050033x. como moderada dependencia espacial (25 <GDE ≤ 75%) y según la clasificación del efecto del modelo exponencial ajustado (Tabla 1) un valor de IDE de 16.36% que, al ser mayor que 13 % implica una fuerte dependencia espacial. Estos valores reflejan que la dependencia espacial es controlada por factores intrínsecos y extrínsecos influenciado por prácticas agrícolas inadecuadas en el suelo (Liu et al., 2014LIU, Z.; ZHOU, W.; SHEN, J.; HE, P.; LEI, Q.; LIANG, G.: “A simple assessment on spatial variability of rice yield and selected soil chemical properties of paddy fields in South China”, Geoderma, 235: 39-47, 2014, ISSN: 0016-7061, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.geoderma.2014.06.027.). La estructura espacial existente permite emplear como técnica de interpolación el kriging que permitirá la confección del mapa de MO estimando valores en puntos no medidos.

La Figura 3 muestra la distribución del contenido de materia orgánica obtenido mediante kriging, donde los mayores valores de materia orgánica se encuentran al norte con tonalidades más oscuras en un transecto irregular en el área de estudio mientras que los valores más bajos están en el centro con tonalidades más claras. Es posible que estos mayores valores de MO en esta zona estén asociados a la presencia de ganado de campesinos al estar el área en barbecho.

Relaciones entre los índices espectrales y el contenido de materia orgánica

La Tabla 4 presenta los estadígrafos del análisis de regresión lineal entre el contenido de materia orgánica del suelo y los índices espectrales de suelo y vegetación. Existe una alta correlación de 0,98 entre el contenido de materia orgánica del suelo, el NDVI y el SAVI lo cual puede estar dado en que ambos índices utilizan las bandas Roja e Infrarroja del sensor para su determinación.

El CI tuvo una correlación de 0.94 con respecto al contenido de materia orgánica del suelo. El coeficiente de determinación arrojó valores cercanos al 100%, el modelo NDVI vs. M.O (%) un 95,61%, SAVI vs. M.O (%) 95,66% y CLAY INDEX vs. M.O (%) de 88,92% por lo que se pude afirmar que los índices espectrales utilizados pueden predecir el contenido de materia orgánica, con un error en su pronóstico en todos los casos dentro de los rangos permisibles en los que son medidos las variables determinadas (Ayoubi et al., 2011AYOUBI, S.; SHAHRI, A.; KARCHEGANI, P.A.; SAHRAWAT, K.L.: Application of artificial neural network (ANN) to predict soil organic matter using remote sensing data in two ecosystems, Ed. InTech Rijeka, Croatia, vol. Biomass and Remote Sensing of Biomass, 181-196 p., 2011, ISBN: 978-953-307-490-0.).

Coeficiente de correlación cercanos a 1,0; de 0,86 y 0,90 fueron encontrados por Wang et al. (2018)WANG, X.; ZHANG, F.; JOHNSON, V.C.: “New methods for improving the remote sensing estimation of soil organic matter content (SOMC) in the Ebinur Lake Wetland National Nature Reserve (ELWNNR) in northwest China”, Remote Sensing of Environment, 218: 104-118, 2018, ISSN: 0034-4257, Publisher: Elsevier, DOI: https://doi.org/10.1016/j.rse.2018.09.020. con el uso de imágenes del Landsat 8 OLI/TIRS con el procesamiento de las bandas roja e infrarroja del sensor obtuvieron una fuerte determinación de un 92,00% en la relación que existe entre el contenido de materia orgánica del suelo y la combinación de bandas del sensor.

Xu et al. (2023)XU, X.; DU, C.; MA, F.; QIU, Z.; ZHOU, J.: “A Framework for High-Resolution Mapping of Soil Organic Matter (SOM) by the Integration of Fourier Mid-Infrared Attenuation Total Reflectance Spectroscopy (FTIR-ATR), Sentinel-2 Images, and DEM Derivatives”, Remote Sensing, 15(4): 1072, 2023, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs15041072. refieren que el uso de sensoramiento remoto es un indicador que establece una relación fuerte entre los datos captados por el sensor y la materia orgánica del suelo con una fuerte relación lineal entre las variables NDVI e índices obtenidos por un modelo digital de elevación (DEM) mientras que no se encontró relación positiva con las bandas espectrales del Sentinel 2. Estudios previos han demostrado que la relación del contenido de materia orgánica del suelo a través de teledetección no se pueden ver como una metodología que se generalice a diversos entornos, sino que es única para cada sitio de estudio, y depende del tipo de sensor que se utilice, las características del suelo, el relieve y el clima (Lamichhane et al., 2019LAMICHHANE, S.; KUMAR, L.; WILSON, B.: “Digital soil mapping algorithms and covariates for soil organic carbon mapping and their implications: A review”, Geoderma, 352: 395-413, 2019, ISSN: 0016-7061, V, DOI: https://doi.org/10.1016/j.geoderma.2019.05.031.).

Por otra parte Prudnikova y Savin (2021)PRUDNIKOVA, E.; SAVIN, I.: “Some peculiarities of arable soil organic matter detection using optical remote sensing data”, Remote Sensing, 13(12): 2313, 2021, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs13122313. encontraron una relación negativa entre el CI y el contenido de materia orgánica del suelo con un coeficiente de determinación para la época de lluvia de 81,00% y en la seca de 84,00% respectivamente al utilizar el Sentinel 2 para estimar el contenido de materia orgánica de un AlbicLuvisols en Rusia.

Existen referencias de otras investigaciones en las que se realizan estimaciones del contenido de carbono orgánico del suelo Sodango et al. (2021)SODANGO, T.H.; SHA, J.; LI, X.; NOSZCZYK, T.; SHANG, J.; ANESEYEE, A.B.; CHAFIK, Z.: “Modeling the spatial dynamics of soil organic carbon using remotely-sensed predictors in Fuzhou city, China”, Remote Sensing, 13(9): 1682, 2021, ISSN: 2072-4292, Publisher: MDPI, DOI: https://doi.org/10.3390/rs1309168., elemento químico que se deriva del contenido de materia orgánica del suelo según Rasel et al. (2017)RASEL, S.; GROEN, T.A.; HUSSIN, Y.A.; DITI, I.J.: “Proxies for soil organic carbon derived from remote sensing”, International journal of applied earth observation and geoinformation, 59: 157-166, 2017, ISSN: 1569-8432, Publisher: Elsevier, DOI: http://dx.doi.org/10.1016/j.jag.2017.03.004.; a partir de índices espectrales de vegetación, en los que al utilizar el NDVI con valores medios de 0,49 (vegetación densa) alcanzaron una correlación de 0,74 (Sankar-Gouri. et al., 2019SANKAR-GOURI., G.; KUMAR-SHIT, P.; REZA-POURGHASEMI, H.: “Soil organic carbon mapping using remote sensing techniques and multivariate regression model”, Geocarto International, 34(2): 215-226, 2019, ISSN: 1010-6049, Publisher: Taylor & Francis, DOI: https://doi.org/10.1080/10106049.2017.1381179.).