Next Article in Journal
Rocky Area Inhabiting Daddy Long-Legs Spiders, Pholcus Walckenaer, 1805 (Araneae: Pholcidae) in Mountainous Mixed Forests from South Korea
Next Article in Special Issue
Influence of the Type and Use of Soil on the Distribution of Organic Carbon and Other Soil Properties in a Sustainable and Resilient Agropolitan System
Previous Article in Journal
Rain-Driven Failure Risk on Forest Roads around Catchment Landforms in Mountainous Areas of Japan
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Forests
Volume 14
Issue 3
10.3390/f14030539
forests-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Mapping Forest Cover and Estimating Soil Organic Matter by GIS-Data and an Empirical Model at the Subnational Level in Mexico

by
Itzel Arroyo
1,
Víctor Tamaríz-Flores
2 and
Rosalía Castelán
2,*
1
Posgrado en Ciencias Ambientales, Instituto de Ciencias, Benemérita Universidad Autónoma de Puebla, Puebla 72580, Mexico
2
Centro de Investigación en Ciencias Agrícolas, Benemérita Universidad Autónoma de Puebla, Puebla 72592, Mexico
*
Author to whom correspondence should be addressed.
Forests 2023, 14(3), 539; https://doi.org/10.3390/f14030539
Submission received: 31 January 2023 / Revised: 3 March 2023 / Accepted: 6 March 2023 / Published: 9 March 2023
(This article belongs to the Special Issue Models for Soil Organic Matter Estimation on Natural Ecosystems, Agricultural Soils and Urban Soils)
Browse Figures
Versions Notes

Abstract

:
Forests play an essential role in climate change as they are the terrestrial ecosystems that store the highest C content in their soils and biomass. Despite this, the lack of information at the subnational level hinders their proper management and conservation. This study aimed to identify the extension and distribution of forests and to develop an empirical model for the spatial prediction of soil organic matter (SOM) in Ixtacamaxtitlan, Puebla, Mexico, based on environmental variables generated through Geographical Information Systems. A supervised classification in Landsat 8 images was used to define the forest cover, and environmental variables related to topography, climate and vegetation were generated. Finally, a Multiple Linear Regression model validated with the leave-one-out cross-validation method was used to examine the relationships between the covariates and the SOM and estimate its content in forest. The results show that the forest cover extension is 41%, with an overall accuracy of 97.7%. The model shows a good fit (R2cv = 0.69, RMSEcv = 1.53). The mean of SOM was 5.2%, and upper values were consistent with higher altitude, precipitation and cooler temperature. Estimating SOM content in forest areas is essential in developing planning strategies at the subnational level to mitigate the harmful effects of climate change.

1. Introduction

The increase in global temperature has caused several modifications in climate patterns, with direct effects on natural and human systems. Deforestation is considered one of the main drivers of the alteration in global atmospheric composition and the production of greenhouse gases (GHGs) [1,2]. In this sense, CO2 is the most abundant GHG, accounting for approximately 74% of anthropogenic GHGs, of which between 17 and 25% is produced by deforestation activity [2,3,4].
Deforestation promotes a significant increase in atmospheric CO2 due to the burning of vegetation, the decomposition of biomass and the exposure of the soil to climatic factors that facilitate erosion and oxidation of organic compounds [2,5,6]. Despite their importance in the global C cycle, the constant degradation to which forest ecosystems have been exposed has led to a significant loss of C previously stored in soils and a modification in carbon fluxes caused by decreased photosynthesis due to the absence of vegetation [7,8,9].
Forest ecosystems play an essential role in climate change, as they contain more C per unit area than any other type of land cover/use [10]. In particular, temperate forests are estimated to store about 70 and 90% of total terrestrial C in their soils and biomass [11,12,13]. Therefore, it is necessary to have up-to-date and accurate information on forest resource distribution, extent, and health that promote the generation of strategies to manage and conserve forest areas and maintain carbon fluxes between soil and plant. However, most of the available and updated information on forest resources is usually at coarse spatial resolutions [14]. Although in Mexico there are recent data published by the National Institute of Statistics and Geography (INEGI, Mexico) on Land Use and Vegetation [15], these are presented at a smaller spatial scale (1:250,000) that may not elucidate the realities at subnational scales. In addition, according to our review, at the national level, there is no published and accessible information on soil organic carbon (SOC) or soil organic matter (SOM), which globally limits international commitments on climate change and the production of accurate reports on C content in forests [16]. However, more importantly, at the subnational level, the lack of information limits decision-making, forest resource management, and sustainable development, mainly affecting rural communities that are the most vulnerable to climate change.
One of the most frequent causes of the lack of information on SOC/SOM at large scales is that soil sampling is complex due to the steep topography of temperate forests, and its determination in the laboratory can be time-consuming and costly [17,18]. Fortunately, remote sensing, in combination with the use of Geographic Information Systems (GIS), is considered a powerful resource that allows the identification and characterization of diverse biophysical factors, their analysis, and monitoring at a low cost. These tools and predictive modeling are considered efficient alternative methods to generate spatiotemporal information and describe soil and vegetation patterns, properties, and processes [19,20,21].
There are different methods for spatial prediction and modeling. The most used are interpolation methods (i.e., kriging, cokriging) [18,22,23,24] and regression methods (i.e., simple or multiple linear regression, generalized linear model, geographically weighted regression) [19,20,21,23,25,26]. These methods are used according to the data’s nature and the sampling type. However, it is considered that prediction methods incorporating environmental variables, such as regression methods, perform better than interpolation models [19,21]. Furthermore, since SOM content is influenced by topographic, climatic, and biotic factors such as vegetation type and productivity [23,24,27,28], it is possible to generate spatial information to estimate SOM in a given region and at a major scale. Despite the importance of organic matter as a key indicator of soil quality [26,29], limited articles worldwide focus on SOM estimation and mapping [24,26,30,31,32] and hardly any focus on temperate mountain forests [33]. For its part, SOC has been the object of several research projects in Mexico [11,16,21,34,35] and Latin America [25,36,37,38,39]. The method generally used to estimate SOC stock is through soil type maps in combination with soil profile sampling at different depths or horizons; this generates SOC maps at the local, regional, or national level. Although the research has meant a breakthrough in the estimation of SOC, amounts are estimated at depths of one meter that do not necessarily reflect the soil depths in the varied land covers, especially in mountainous temperate forest ecosystems. Likewise, the mapping and modeling of SOC generated in both natural and productive systems is generalized without considering, in many cases, that the factors that influence the SOC stock are distinct for the different land uses/vegetation covers [26,40]. In this regard, publications on SOC in temperate mountain forests are also scarce [41,42].
The municipality of Ixtacamaxtitlan, located in the Sierra Norte of Puebla, is endowed with contrasts between vegetation, relief and climate that extend along a mountainous topography with temperate coniferous forests. Information on these forests is usually obtained from national forest inventories that do not present details at the municipal or local scale or from state inventories that were last updated in 2013. In this sense, there is a lack of recent information on the distribution and extent of forest areas or their SOM content. In this framework, this study aimed to identify the forest areas in the municipality of Ixtacamaxtitlan, Puebla, Mexico, and to develop an empirical model for the spatial prediction of SOM based on data generated through GIS and set for the study site.

2. Materials and Methods

2.1. Study Site and Soil Sampling

Ixtacamaxtitlan municipality covers an area of approximately 562 km2. It is located in the north of the state of Puebla, Mexico, between coordinates 19°27′ and 19°45′ N and 97°41′ and 98°03′ W (Figure 1). The study site has a temperate and semi-cold sub-humid climate with rainfall in summer. The temperature varies between 5 and 22 °C throughout the year, and the mean annual precipitation is typically 800 mm [15]. The altitude varies between 2000 and 3500 masl with a volcanic mountain range topography of steep slopes. The geology is dominated by volcanic rocks represented by pyroclastic and andesitic deposits. The predominant soils are Andosols [43]. The main land cover/uses are agriculture, grassland, desert scrub, and forest. The forest area is mainly composed of coniferous forests (pine, fir, juniper), which in some areas are associated with broadleaf trees (pine-oak forests) forming small patches of mixed forests. The most representative species are Pinus patula, Pinus montezumae, Abies religiosa, Juniperus spp. and, to a lesser extent, Quercus spp. [44].
The samples used in this study were taken at a depth of 0 to 30 cm with a soil auger, totaling 22 composite soil samples collected in two coniferous forest areas in the municipality of Ixtacamaxtitlan, Puebla, Mexico. The forest areas were selected to represent the whole area with different relief forms. The sampling represented an area equivalent to 6% of the forest areas, and was performed in 2020. The CO content in soil samples was determined through the wet chemical oxidation procedure of the Walkley and Black method [45]. Then, the result was multiplied by the Van Bemmelen factor to obtain the SOM content in percent [46].

2.2. Satellite Data Collection and Pre-Processing

Two Landsat 8 multispectral space-borne images level-1 precision and terrain-corrected product (L1TP) with a spatial resolution of 30 m were downloaded freely from https://earthexplorer.usgs.gov/ (accessed on 25 November 2022). The images were dated 10 February and 7 October 2020. The study area was in a single scene (path 025, raw 046). The images correspond to days without clouds in the study area, so atmospheric corrections were unnecessary [47]. Likewise, the products are geometrically corrected and orthorectified [14,48].
A Digital Elevation Model (DEM) taken by the ALOS/PALSAR sensor with Beam Mode FBD and a Radiometrically Terrain-Corrected (RTC) processing level was downloaded freely through the Alaska Satellite Facility portal (https://search.asf.alaska.edu/#/ accessed on 14 October 2022). The DEM was composed of an image in raster format corresponding to one scene (path 183, frame 380, 17 November 2010) with a spatial resolution of 12.5 m.
Climatic data (temperature and precipitation) were obtained through the Network of Weather Stations from the National Weather Service (SMN, Mexico) [49]. Fifteen stations located within the municipality and at 15 km or less outside were considered. Therefore, the climatic data are the long-term average of those stations’ 1993 to 2018 period.
All data were geographically referenced to the UTM Zone 14 N coordinate system based on the WGS 84 ellipsoid. The cartographic processing was performed in ArcGIS 10.4 software.

2.3. Forest Cover and Accuracy Assessment

In the Landsat 8 satellite image dated 7 October 2020, the RGB bands 5, 4, and 3 were composited to obtain a False Colour Composite (FCC) image since this combination highlights the vegetation of the region studied, assigning a dark red color to intense vegetation such as forests [50].
Forest cover in the study area was undertaken through the supervised maximum likelihood algorithm method in the ArcGIS 10.4 program [51]. The spectral signature was defined by placing 50 training sites in each land cover class through visual interpretation [52]. As reference training sites, we used Land Use and Vegetation data and high spatial resolution imagery from Google Earth Pro and the World Imagery folder of ArcGIS 10.4 [15,51,53]. The assigned land cover classes in the study area were five: forest, secondary forest, cropland, grassland, and settlement. The secondary forest was included as a class to avoid spectral confusion due to the edge effect with desert scrub.
After obtaining the classified image, a reclassification was performed to define the areas with forest/non-forest cover (Table 1). In addition, polygons of fewer than two hectares were eliminated to reduce noise on the map.
A confusion matrix was performed using 1000 randomized sampling points based on Land Use and Vegetation maps [15] and photo interpretation in Google Earth Pro software [53] to validate and assess the accuracy of the forest cover distribution and extension. The number of validation points used for each new class (Forest/non-forest) was assigned according to their approximate extension in the territory (40/60%, respectively) [54].
The confusion matrix identifies the number of pixels classified in a given class, and it is reported in terms of overall accuracy (OA), user’s accuracy (UA), producer’s accuracy (PA), and kappa coefficient. Overall accuracy refers to the ratio of correctly classified pixels to those evaluated. UA and PA are individually assessed for each classification result. UA refers to the ratio between the total number of correctly classified points of a class and the total number of reference validation points, in this case, 1000. PA is the ratio of correctly classified reference points to a class’s total number of validation points [7]. Finally, the kappa coefficient evaluates the classification performance and the agreement with reality in the field.
Additionally, two specific measures for forest classification were quantified. These are the true positive rate (TPR) and the false positive rate (FPR) [8]. The TPR shows the rate of correctly classified forests, while the FPR shows the proportion of samples that were classified as forests but were not.

2.4. Environmental Factors

Information on environmental variables, including topographic (elevation and slope percentage), climatic (temperature and precipitation), and vegetation (Normalized Difference Vegetation Index), were used. In addition, composite covariables were used to generate more information in the model [55]. These variables were based on the factors used to calculate soil loss according to the Universal Soil Loss Equation (USLE) [56]. Specifically, the variables were the rainfall erosivity factor (R-factor) and the slope length and grade factor (LS-Factor).
Landsat 8 images were processed to generate the Normalized Difference Vegetation Index (NDVI) [57], which is used to indicate the vegetation’s distribution, density, and health [58]. The index is dimensionless with a range value from −1 to +1 and is obtained using bands 4 and 5, the visible red band (RED), and the near-infrared band (NIR), respectively. NDVI is a dynamic indicator that presents seasonal variability. To reduce this variability, the average NDVI was obtained from the two Landsat 8 satellite images, corresponding to the end of the dry season (February) and the end of the rainy season (October) [59]. The calculation was carried out according to the formula:
NDVI = N I R − R e d N I R + R e d
Topographic variables such as slope (%), altitude and LS-factor were generated through DEM. The slope in percentage units was produced through the ArcGIS 10.4 program with the slope tool [60,61].
The topographic factor length-slope (LS) was modeled using the approach of Foster et al. [62] and Desmet and Govers [63] to obtain the slope length factor, and the slope grade factor (S) was modeled according to McCool et al. [64]. This methodology is used for LS-factor modeling in GIS as it can identify complex topography and can be applied to landscape scale modeling [65]. LS-factor is calculated as follows:
F = s i n β / 0.0896 3 s i n β 0.8 + 0.56
m = F 1 + F
L i , j = ( A i , j + D 2 ) m + 1 − A i , j m + 1 x m ∗ D m + 2 ∗ 22.13 m  
S i , j = 10.8   s i n θ i , j + 0.03 θ i , j < 0.09 16.8   s i n θ i , j − 0.5 θ i , j ≥ 0.09
where F is the slope angle factor, β is the slope angle in radians, m = slope length exponent, L(i,j) is the slope length in the pixel (i,j), A(i,j) is the flow accumulation in the pixel (i,j), D is the pixel side (30 m), x is the shape coefficient (1), S(i,j) is the steepness slope (i,j) and θ(i,j) is the slope angle (i,j) in grades.
The climatic variables of precipitation, temperature and R-factor were generated using the kriging interpolation method. This interpolation method is characterized by its high accuracy for modeling environmental variables in areas with pronounced relief dynamics and high spatial variability [22,66].
The rainfall erosivity factor (R-factor) was mapped using the Modified Fournier Index (MFI) by the equation proposed by Arnoldus et al. [67]. This equation is easily applicable in places without continuous precipitation records [65]. Finally, the MFI was replaced to obtain the R-factor [68]:
M F I = ∑ i = 1 12 P i 2 P t
R = 30.4 M F I + 28.3
where MFI is the modified Fournier index, Pi is the mean monthly precipitation, Pt is the mean annual precipitation, and R is the rainfall erosivity (R-factor).
All layers were resampled, adjusted, and spliced to a spatial resolution of 30 m using as a basis the position and cell size of the satellite images obtained. The DEM was resampled using the bilinear method [7]. Subsequently, with the definition of the areas with forest cover, the digital maps of the environmental variables were cut out.

2.5. Model Construction

The main objective of this study was to develop a valid empirical model to estimate the SOM content in the forest areas of Ixtacamaxtitlan, Puebla, Mexico. Multiple Linear Regression (MLR) is a statistical technique whose objective is to determine to what extent the independent variables can explain the response or dependent variable; in turn, it is used to predict the response of interest through the covariates [20]. For this purpose, data on the independent variables generated above were extracted at the soil sampling points.
The Shapiro-Wilk test was performed to check the normality of SOM (p > 0.05) [69]. Subsequently, a Pearson correlation analysis was performed to identify the relationships between the variables. The stepwise technique was used to select the most significant variables for the SOM estimation. This technique adds or eliminates variables according to their significance level at each step and chooses the model with the lowest Akaike Information Criterion (AIC) [70,71]. With the model defined, the parameters of the MLR, the coefficient of determination (R2) and the root-mean-square error (RMSE) were estimated.
The assumptions of the MLR were checked. The normality of residuals was verified by the Shapiro-Wilk test, whose power has proven to be better than other normality tests [72]. The Breusch-Pagan, Variance Inflation Factor (VIF) and Cook’s distance tests were also performed to evaluate the heteroskedasticity, multicollinearity and the existence of potential outliers, respectively. The leave-one-out cross-validation (LOOCV) method was used to validate the model. LOOCV is an effective technique to evaluate the model’s accuracy with limited samples [73]. This process calibrates and validates as many times as the number of samples (k = 22), leaving out one sample each time and using the rest to predict the value of the omitted sample [74]. Model performance is evaluated by the coefficient of determination (R2cv), and the root mean square error (RMSEcv) [75], which are considered the accurate metrics of the prediction model. The ratio of prediction to deviation (RPD), which identifies the accuracy and reliability of the results, was also evaluated. The RPD is the ratio between the standard deviation of the observed SOM and the RMSEcv.
All statistical analyses were performed in the statistical software R Studio Version 2022.07.2 using additional packages. Finally, the coefficients of the predictive model were applied to the area of interest using map algebra in the ArcGIS program [51].

3. Results

3.1. Forest Cover Accuracy Assessment

The accuracy assessment of forest/non-forest cover performed using Landsat 8 L1TP satellite images with the supervised classification method showed a high level of agreement with the distribution of forest on the ground (Table 2). The overall accuracy was 97.7%, with a kappa coefficient of 0.952. UA and PA are greater than 95% for both classifications (forest/non-forest). UA for forest cover is higher by almost 2% than UA for non-forest cover, reaching nearly 99%. PA shows a lower percentage for forest cover (95.64) than non-forest cover, whose accuracy is slightly higher than 99%. The rate of correctly classified forests (TPR) is 98.8%. According to the results of the present study, the forest area covers 228.68 km2, representing 40.7% of Ixtacamaxtitlan.

3.2. Forest Cover and Environmental Variables

The spatial distribution of forest cover in the study area and the environmental variables analyzed in this forest cover are shown in Figure 2. Forests are distributed throughout the territory, with a higher proportion to the south of Ixtacamaxtitlan. In contrast, areas of fragmented forests are observed mainly in the central part of the study area (Figure 2a).
The descriptive statistics of the environmental variables in the forest areas are presented in Table 3. NDVI in forest cover shows a value < 0.45, with a mean of 0.28 ± 0.05. In the spatial distribution of NDVI (Figure 2b), values < 2.3 in the center and northwest can be observed, associated with the location where more significant fragmentation is observed, while the highest values are kept in the south of the studied area.
The altitude in forest zones is between 2056 and 3506 masl. In its spatial distribution (Figure 2c), it can be observed that the elevations lower than the average (2555 masl) are located in the central strip from west to east of the territory. Associated with the above, the mean slope is 39.85 ± 19.44% with a maximum of 201%, indicating a steep terrain, as seen in Figure 2d.
The climatic variables present mean precipitation and temperature of 749 ± 81.45 mm and 14.50 ± 2.80 °C, respectively. The highest precipitation (>800) is distributed in the forested areas located in the southern portion (Figure 2e), associated with higher NDVI values, higher altitude, and lower temperature (<12.1 °C) (Figure 2f).
As for the LS-factor in forest areas, values of up to 310 are found in the highest parts (Figure 2g). Rainfall erosivity (R-factor), on the other hand, ranges between 1856.80 and 3965 MJ mm ha−1 yr−1 with a mean of 3039.90 and is observed to be associated with higher altitudes (Figure 2h).

3.3. Statistical Analyses of Soil Organic Matter and Environmental Variables

The data in Table 4 shows the summary statistics of the SOM and environmental variables at the sampling sites. The SOM ranged from 1.10 to 10.28, with an arithmetic mean of 5.41 ± 2.8, meeting the normality assumption via the Shapiro-Wilk test (p > 0.05). The mean values of the topographic variables were 2670 masl, 31.93 % and 5.32 for altitude, slope, and LS-factor. The climatic variables of precipitation, temperature, and R-factor had a mean of 773.6 mm, 14.26 °C and 3088 MJ mm ha−1 yr−1, respectively. Finally, NDVI presented a mean of 0.23.

3.4. Relationship between Soil Organic Matter and Environmental Variables

The Pearson correlation analysis between the environmental variables and SOM is presented in Figure 3. Pearson’s correlation coefficient (Pearson’s r) evaluates the degree of the linear relationship between the variables. Pearson’s r shows that SOM presents a strong positive and significant relationship with precipitation > altitude > R-factor > NDVI, while a strong negative correlation is observed with temperature. Significant correlations also exist between environmental variables such as altitude with precipitation, R-factor, and temperature. In contrast, NDVI shows weak correlations with the rest of the variables. As expected, precipitation with the R-factor and slope with the LS-factor shows a highly strong correlation.

3.5. Multiple Linear Regression Model Analysis

The MLR model using the stepwise technique, obtained a coefficient of determination (R2) of 0.78. However, the cross-validation shows that the accurate coefficient of determination (R2cv) is 0.69 (Table 5). The model was constructed with the variables LS-factor, NDVI and precipitation. All the model variables showed coefficients significantly influencing SOM (p < 0.01). The above indicates that these environmental variables explain 69% of the variability of SOM in soil samples. The model presented an RMSEcv of 1.53% and an RPD of 1.83.
The absolute value of the coefficients shows the degree of influence of the environmental variables on the prediction of SOM; in this sense, it is observed that the NDVI coefficient (22.20) had the highest absolute value compared to the coefficients of the LS-factor (0.2839) and precipitation (0.0456). At the same time, the positive sign of the covariates indicates a positive correlation between SOM and the explanatory variables, as observed in the correlation analysis.
The model met the normality assumption of the residuals by the Shapiro-Wilk test (p > 0.05) and that of homoscedasticity by the Breusch–Pagan test (p > 0.05). Likewise, all the variables presented a variance inflation factor of less than 2, so multicollinearity was not a problem in the model. Furthermore, according to Cook’s distance test, outliers were absent (<0.03).
The scatter plot between the observed SOM and that predicted by the MLR model can be seen in Figure 4, which showed a good correlation (r = 0.88).

3.6. Mapping Soil Organic Matter in Forest Cover

The estimated SOM in the forest areas showed a mean of 5.2%, with a coefficient of variation of 80% and a range of values between the first and third quartiles from 2.42% to 7.76%. However, the presence of abnormal values in the prediction was identified.
Figure 5 shows the distribution of SOM estimation in forest areas. SOM was classified according to the Mexican Official Standard [46] as organic matter content: very low (<4%), low (4.1%–6%), medium (6.1%–10.9%), high (11%–16%) and very high (>16%). In this sense, about 38% of the forest areas presented a very low SOM content. In comparison, 30% showed a range considered medium, followed by low, high, and very high content, with an extension in the forest areas of 20, 11 and 1%, respectively.
Higher SOM values were related to higher altitudes, as well as to higher precipitation and lower temperature. High slopes but at lower altitude levels do not have high organic matter values. Likewise, forest density and health represented by NDVI were associated with higher SOM content.

4. Discussion

4.1. Forest Cover and Environmental Variables

This research shows that Landsat 8 L1TP satellite images with an FCC and a supervised classification are an efficient resource at the subnational level to know the distribution of forest cover at a high spatial resolution (30 m). It should be noted that the data used are accessible and free. Moreover, the easily applicable methodology has been widely used to identify land uses, forest areas and their boundaries [14,48,50,76,77]. Therefore, the above shows a reasonable basis for overcoming the challenges concerning identifying forest areas at local scales, where knowledge of forest resources is limited.
According to our study, forest cover represented nearly 41% of the municipality. Compared to other research, our results present an area between 3 and 4.5% higher in forest extension. For example, Arroyo et al. [78] reported an extent of forested surfaces of 36.5% in 2000 and 2015 (spatial resolution 250 m), while official Mexican data on land use and vegetation [15] define a forested surface in Ixtacamaxtitlan of 37.8% as of 2018 (scale 1:250,000). However, our results are more accurate, as we found an overall accuracy of 97.7%. In contrast, the land use and vegetation data at national-scale investigations place their accuracy between 73 and 89% [79,80] and in the definition of land cover in Ixtacamaxtitlan [78], their accuracy was 85%.
The difference in our results compared to the above data is not very wide. In contrast, the Forest and Soil Inventory of the state of Puebla, Mexico [44] indicates a forest extension of primary and secondary coniferous and broadleaf forests in only 17% of the territory as of 2013, an extension much smaller than that found in our research. However, making temporal comparisons regarding the increase and decrease of forest cover is a limitation due to methodological differences and spatial resolution. This situation highlights the need to produce accurate information on the distribution of forest resources and their spatial and temporal dynamics at major scales. In this regard, some authors [81,82] mention that official data from Mexico do not show forest degradation processes, but at major scale, it is possible to identify severely degraded or deforested areas; or, in their case, recovered or restored. Therefore, it is essential to generate updated information at the subnational level to detect degradation or recovery processes associated with changes in carbon fluxes in an effort to promote their conservation, reforestation and sustainable management.
The presence of forests was related to the climate with humid and sub-humid regimes and the ease of forest vegetation to develop on steep slopes and high altitudes. This situation implies a shorter duration in the runoff water concentration related to higher precipitation and lower temperature. Moreover, erosivity, considered moderate [83], is related to soil loss by entrainment of particles from the superficial horizons where there is a higher organic matter content, a process that is magnified with higher values of the LS-factor. However, the presence of forest vegetation reduces the impact of precipitation on the soil, and roots—which are more significant in temperate forests [84]—retain the soil, preventing it from runoff. This scenario enhances the importance of forest cover in an area characterized by diverse biophysical contrasts.
On the other hand, analyzing NDVI, low values are found in forest areas. This condition could be associated with a lower tree density represented by areas of secondary forest vegetation. In addition, fungal diseases have been reported in Pinus sp. in the Sierra Norte de Puebla [85], to which Ixtacamaxtitlan belongs, which might be promoting a decrease in forest health. Gutiérrez-Flores et al. [86] states that fungal diseases in forest areas are one of Mexico’s leading causes of decreased forest biomass. Likewise, an increase in droughts has been reported since 2019 in the central region of Mexico [87] that could promote water stress in forest species and, therefore, a decrease in primary productivity.
Several authors have reported that altitude, slope and precipitation present high standard deviations [19,26,88], as it was found in our research, implying greater spatial variability in the forest areas. On the other hand, the covariates mean at the sampling points are akin to those in the forest, except for the LS-factor, which is higher in the sampling mean. In this sense, although the sampling points covered a wide range of values of the environmental variables, some forest areas presented steep slopes (>101%), areas of higher (>2808 masl) or lower altitude (<2578 masl) and LS-factor values superior to those sampled (>14.25), showing values of environmental variables that were not covered by the soil samples. This variability can be verified by comparing the coefficients of variation, where DEM, R-factor and precipitation present a coefficient of variation with a difference greater than 40% in the forest zones.

4.2. Soil Organic Matter and Environmental Variables

The variables with the most significant influence on SOM content were altitude, NDVI, precipitation and temperature; this relationship has also been highlighted in previous research [18,20,23,28]; in addition to variables such as geology, soil texture and bulk density, which also have a significant effect on SOM/SOC. However, these variables were not measured in the present research. In the case of geology, this is because most of the territory has a volcanic geology with predominantly Andosol soils [43]. While in the case of soil texture and bulk density, there is no GIS information on these environmental variables in the municipality; consequently, it would have to be generated through field sampling and modeling.
On the other hand, despite the strong relationships between the covariates with SOM content; those selected by the stepwise technique were NDVI, LS-factor and precipitation; even though the LS-factor presented a weak positive association with SOM; the contribution of the covariate to the model was significant and independent. Several authors [19,23,24,25] mention that the explanatory variables selected with this technique are considered the most appropriate and with the highest weight. Likewise, it highlights that this regression technique set a topographic variable, a climatic variable and a vegetation variable which can help us better understand the dynamics of organic matter and SOC concerning these environmental variables [25].
The importance of NDVI as an indicator of organic matter input has been supported previously [89,90,91], being the covariable with the most significant weight in many investigations [23,24,31,37,41], as found in this research. Furthermore, the NDVI relevance as an environmental variable in predicting SOM in the temperate forests of Ixtacamaxtitlán highlights the importance of the amount of biomass and vegetation cover as drivers of organic matter accumulation [37,41].
Generally, topographic factors have more significance in prediction models [18,30]. In this research, LS-factor was the second most important variable in explaining changes in SOM with a positive effect on its content. Although similar results have been found in other investigations [23,92,93], our results differ from those reported by Kumar et al. [19], Liu et al. [94] and K. Wang et al. [26], who found a negative effect of slope on SOM in their prediction models. However, the positive effect between SOM content and LS-factor could be due to the difficulty of developing agricultural activities on steep slopes. This situation has prevented land use change, which promotes the permanence of forest vegetation cover, providing higher SOM contents. In addition, the analysis in this research was conducted only for forested areas so that the vegetation cover is assured. In contrast, in the aforementioned research, the SOM/SOC content is predicted in different land covers.
Climatic factors have also been essential in SOM predictions [35,95]. In this sense, precipitation was the third factor that explained the changes in SOM content; climatic factors are closely related to organic matter as they significantly influence the type of vegetation, its growth and microbial activity [95].

4.3. Soil Organic Matter Estimation Model

SOM content in the Ixtacamaxtitlán forests averaged 5.2%. This result is similar to that reported by Gallegos and Bautista [96], who generated a dataset listing the chemical and physical properties of the soils from a series of soil profiles. It was identified that the Andosols of the mountainous systems of the Transmexican Neovolcanic Axis—of which Ixcatamaxtitlan is a part—present a mean SOM content of 5.97% in the first 29 cm. On the other hand, our results differ from those reported by Galicia et al. [11] in Andosol soils of temperate pine forests in Durango, Mexico, where they detected SOM contents between 11%–17%; however, these estimates were made through soil profiles, so the depth of the superficial horizon is heterogeneous. Research in mountainous areas [24,26] generally detected a higher SOM content than in plain and lower altitude areas.
The model presented a good prediction (R2cv = 0.69). At the same time, the RMSEcv shows a value lower than the standard deviation of the SOM data at the sampling sites, with an RPD of 1.83, indicating that the model presents an acceptable fit [24,75,97]. Likewise, the mean SOM in the soil samples and the mean SOM predicted for the forest areas are similar. However, analyzing the observed values against the predicted ones, we can find that there is an underestimation of the low values and an overestimation of the highest values, this is reflected in the prediction at the level of forest areas, where the minimum and maximum values outside the range of SOM evaluation are observed. Regarding this, it is necessary to highlight that one of the disadvantages of the MLR is that it considers the global effect of environmental variables but ignores the impact of spatial heterogeneity [23]. Therefore, relatively large estimation errors are inevitable, especially in areas not represented through sampling.
Although different spatial modeling and prediction methods can show better performances (e.g., geographically weighted regression, random forests or artificial neural networks), MLR models and the stepwise technique are still used as the primary method to define the most important variables in the construction of prediction models [23,40,98,99]. Likewise, the MLR has shown a similar performance when comparisons are made between developed models. [23,25,93,94,100].
SOM exhibits greater spatial variability than SOC and a greater dependence on relationships with environmental variables [25], highlighting the need for a major number of samplings covering a more comprehensive range of the variability of covariates in forest areas to increase SOM predictions’ accuracy. Despite this, quantifying the variation in SOM explained by environmental factors in forested areas has important implications for sustainable land management and decision-making. Moreover, it leads to finding more accurate predictions of SOC accumulation. With the information on forest cover determined in this research, it will be possible to establish systematic sampling in these areas to allow the use of mixed models for prediction and obtain more accurate results. In addition, future research should focus on defining the composition of the forest, its density and the type of forest species present in Ixtacamaxtitlán.

5. Conclusions

The information reported in this research represents the distribution and extension, to the year 2020, of forest areas at the subnational level in the municipality of Ixcatacamaxtitlán, Puebla, Mexico, and the characterization of environmental variables related to topography, climate and vegetation generated through GIS. These data allowed the development of an empirical model to estimate the SOM content and define its spatial distribution in the study area.
The results show that forest areas in Ixtacamaxtitlan cover an extension of almost 41% of the municipality’s total area. Furthermore, these forest areas presented a wide range of topographic, climatic and vegetation environmental variables affecting SOM content. The stepwise technique allowed the development of an MLR model with significant covariates that provided information on their effect on the SOM, predicting its content and mapping its distribution.
The model presented a good fit; however, it is necessary to carry out a more significant number of samples to help us better understand the dynamics of the SOM in the study area. Likewise, it is recommended that the soil samples present a greater heterogeneity to represent the high variability of the biophysical factors in order to use prediction models associated with spatial variability.
Generating knowledge about forest areas at major scales is a critical step for developing conservation strategies and adequate management of forest resources. The results of this study provide an important basis for the search for a better understanding of the C cycle in forest areas at the subnational level.

Author Contributions

Conceptualization, I.A., V.T.-F. and R.C.; methodology, I.A., V.T.-F. and R.C.; formal analysis, I.A. and R.C.; investigation, I.A., V.T.-F. and R.C.; resources, V.T-F. and R.C.; data curation, I.A.; writing—original draft preparation, I.A. and R.C.; writing—review and editing, I.A. and R.C.; supervision, R.C.; funding acquisition, R.C. All authors have read and agreed to the published version of the manuscript.

Funding

This research project was funded by Consejo Nacional de Ciencia y Tecnología (CONACYT) (Grant Number 712144). The APC was funded by Benemérita Universidad Autónoma de Puebla.

Data Availability Statement

Not applicable.

Acknowledgments

We thank Instituto de Ciencias (ICUAP) for the administrative and technical support.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Landholm, D.M.; Pradhan, P.; Wegmann, P.; Sánchez, M.A.R.; Salazar, J.C.S.; Kropp, J.P. Reducing Deforestation and Improving Livestock Productivity: Greenhouse Gas Mitigation Potential of Silvopastoral Systems in Caquetá. Environ. Res. Lett. 2019, 14, 114007. [Google Scholar] [CrossRef]
  2. FAO. Carbono Orgánico Del Suelo: El Potencial Oculto; Lefévre, C., Rekik, F.V.A., Wiese, L., Eds.; Food and Agriculture Organization of the United Nations: Rome, Italy, 2017. [Google Scholar]
  3. Rodríguez Becerra, M.; Mance, H. Cambio Climático: Lo Que Está En Juego, Primera; Valderrama, J.A., Ed.; Foro Nacional Ambiental: Bogotá, Colombia, 2009; Volume s10-II. [Google Scholar]
  4. FAO. Global Forest Resource Assessment 2005: Progress towards Sustainable Forest Management; Food and Agriculture Organization of the United Nations: Rome, Italy, 2006; Volume 147. [Google Scholar]
  5. Pérez-Ramírez, S.; Ramírez, M.I.; Jaramillo-López, P.F.; Bautista, F. Contenido de Carbono Orgánico En El Suelo Bajo Diferentes Condiciones Forestales: Reserva de La Biosfera Mariposa Monarca, México. Rev. Chapingo Ser. Cienc. For. Del Ambient. 2013, 19, 157–173. [Google Scholar] [CrossRef] [Green Version]
  6. Zamora-Morales, B.P.; Mendoza-Cariño, M.; Sangerman-Jarquín, D.M.; Quevedo Nolasco, A.; Navarro Bravo, A. El Manejo Del Suelo En La Conservación de Carbono Orgánico. Rev. Mex. Cienc. Agrícolas 2018, 9, 1787–1799. [Google Scholar] [CrossRef] [Green Version]
  7. Hosciło, A.; Lewandowska, A. Mapping Forest Type and Tree Species on a Regional Scale Using Multi-Temporal Sentinel-2 Data. Remote Sens. 2019, 11, 929. [Google Scholar] [CrossRef] [Green Version]
  8. Zhang, D.; Wang, H.; Wang, X.; Lü, Z. Accuracy Assessment of the Global Forest Watch Tree Cover 2000 in China. Int. J. Appl. Earth Obs. Geoinf. 2020, 87, 102033. [Google Scholar] [CrossRef]
  9. Cano-Flores, O.; Vela-Correa, G.; Acevedo-Sandoval, O.A.; Valera-Pérez, M.Á. Organic Carbon Concentrations in the Woodland and Soils of the Protected Natural Area “El Faro” in Tlalmanalco, Estado de Mexico. Terra Latinoam. 2020, 38, 895–905. [Google Scholar] [CrossRef]
  10. Charro, E.; Hernández, S.; Martín, J.; Moyano, A.; Ruiz, N. Estimación Del Secuestro de Carbono En Suelos Bajo Masas Forestales de Pinus Halepensis En Castilla y León (España). Cuad. La SECF 2008, 130, 125–130. [Google Scholar]
  11. Galicia, L.; Gamboa, M.; Cram, S.; Vergara, B.; Ramírez, V.; Saynes, V.; Siebe, C. Almacén y Dinámica Del Carbono Orgánico Del Suelo En Bosques Templados de México. Terra Latinoam. 2016, 34, 1–29. [Google Scholar]
  12. Houghton, R.A.; Hall, F.; Goetz, S.J. Importance of Biomass in the Global Carbon Cycle. J. Geophys. Res. Biog. 2009, 114, 1–13. [Google Scholar] [CrossRef]
  13. Kumar, K.K.; Nagai, M.; Witayangkurn, A.; Kritiyutanant, K.; Nakamura, S. Above Ground Biomass Assessment from Combined Optical and SAR Remote Sensing Data in Surat Thani Province, Thailand. J. Geogr. Inf. Syst. 2016, 8, 506–516. [Google Scholar] [CrossRef] [Green Version]
  14. Song, D.-X.; Wang, Z.; He, T.; Wang, H.; Liang, S. Estimation and Validation of 30 m Fractional Vegetation Cover over China through Integrated Use of Landsat 8 and Gaofen 2 Data. Sci. Remote Sens. 2022, 6, 1–20. [Google Scholar] [CrossRef]
  15. INEGI. Conjunto de Datos Vectoriales de Uso de Suelo y Vegetación. Escala 1:250 000, Serie VII. Conjunto Nacional.’, Escala: 1:250 000. Edición: 1; Instituto Nacional de Estadística y Geografía: Aguascalientes, México, 2021.
  16. Guevara, M.; Vargas, R. Soil Organic Carbon Predictions across Mexico at 1 m of Soil Depth and 90 m of Spatial Resolution (1999–2009). Terra Latinoam. 2021, 39, 1–19. [Google Scholar] [CrossRef]
  17. Abdel-Fattah, M.K. Linear Regression Models to Estimate Exchangeable Sodium Percentage and Bulk Density of Salt Affected Soils in Sahl El-Hossinia, El-Sharkia Governorate, Egypt. Commun. Soil Sci. Plant Anal. 2019, 50, 2074–2087. [Google Scholar] [CrossRef]
  18. Van Huynh, C.; Pham, T.G.; Nguyen, L.H.K.; Nguyen, H.T.; Nguyen, P.T.; Le, Q.N.P.; Tran, P.T.; Nguyen, M.T.H.; Tran, T.T.A. Application GIS and Remote Sensing for Soil Organic Carbon Mapping in a Farm-Scale in the Hilly Area of Central Vietnam. Air Soil Water Res. 2022, 15, 4777. [Google Scholar] [CrossRef]
  19. Kumar, S.; Lal, R.; Liu, D. A Geographically Weighted Regression Kriging Approach for Mapping Soil Organic Carbon Stock. Geoderma 2012, 189–190, 627–634. [Google Scholar] [CrossRef]
  20. Mishra, U.; Lal, R.; Liu, D.; Van Meirvenne, M. Predicting the Spatial Variation of the Soil Organic Carbon Pool at a Regional Scale. Soil Sci. Soc. Am. J. 2010, 74, 906–914. [Google Scholar] [CrossRef]
  21. Pérez-Vázquez, Z.R.; Ángeles-Pérez, G.; Chávez-Vergara, B.; Valdez-Lazalde, J.R.; Ramírez-Guzmán, M.E. Spatial Approach for Modeling Litter Carbon in Forests under Management for Timber Production. Madera Y Bosques 2021, 27, 1–19. [Google Scholar] [CrossRef]
  22. Kuo, P.F.; Huang, T.E.; Putra, I.G.B. Comparing Kriging Estimators Using Weather Station Data and Local Greenhouse Sensors. Sensors 2021, 21, 1853. [Google Scholar] [CrossRef]
  23. Liu, T.; Zhang, H.; Shi, T. Modeling and Predictive Mapping of Soil Organic Carbon Density in a Small-Scale Area Using Geographically Weighted Regression Kriging Approach. Sustainbility 2020, 12, 9330. [Google Scholar] [CrossRef]
  24. Piccini, C.; Francaviglia, R.; Marchetti, A. Predicted Maps for Soil Organic Matter Evaluation: The Case of Abruzzo Region (Italy). Land 2020, 9, 349. [Google Scholar] [CrossRef]
  25. Costa, E.M.; Tassinari, W.d.S.; Pinheiro, H.S.K.; Beutler, S.J.; dos Anjos, L.H.C. Mapping Soil Organic Carbon and Organic Matter Fractions by Geographically Weighted Regression. J. Environ. Qual. 2018, 47, 718–725. [Google Scholar] [CrossRef] [PubMed]
  26. Wang, K.; Zhang, C.R.; Li, W.D.; Lin, J.; Zhang, D.X. Mapping Soil Organic Matter with Limited Sample Data Using Geographically Weighted Regression. J. Spat. Sci. 2014, 59, 91–106. [Google Scholar] [CrossRef]
  27. Minasny, B.; McBratney, A.B.; Malone, B.P.; Wheeler, I. Digital Mapping of Soil Carbon. Adv. Agron. 2013, 118, 1–47. [Google Scholar] [CrossRef]
  28. Hu, P.L.; Liu, S.J.; Ye, Y.Y.; Zhang, W.; Wang, K.L.; Su, Y.R. Effects of Environmental Factors on Soil Organic Carbon under Natural or Managed Vegetation Restoration. L. Degrad. Dev. 2018, 29, 387–397. [Google Scholar] [CrossRef]
  29. Orr, B.J.; Cowie, A.L.; Castillo Sánchez, V.M.; Chasek, P.; Crossman, N.D.; Erlewein, A.; Louwagie, G.; Maron, M.; Metternicht, G.I.; Minelli, S.; et al. Scientific Conceptual Framework for Land Degradation Neutrality. A Report of the Science-Policy Interface; United Nations Convention to Combat Desertification (UNCCD): Bonn, Germany, 2017; p. 140. [Google Scholar]
  30. Pouladi, N.; Møller, A.B.; Tabatabai, S.; Greve, M.H. Mapping Soil Organic Matter Contents at Field Level with Cubist, Random Forest and Kriging. Geoderma 2019, 342, 85–92. [Google Scholar] [CrossRef]
  31. Mirzaee, S.; Ghorbani-Dashtaki, S.; Mohammadi, J.; Asadi, H.; Asadzadeh, F. Spatial Variability of Soil Organic Matter Using Remote Sensing Data. Catena 2016, 145, 118–127. [Google Scholar] [CrossRef]
  32. Zeng, C.; Yang, L.; Zhu, A.X.; Rossiter, D.G.; Liu, J.; Liu, J.; Qin, C.; Wang, D. Mapping Soil Organic Matter Concentration at Different Scales Using a Mixed Geographically Weighted Regression Method. Geoderma 2016, 281, 69–82. [Google Scholar] [CrossRef]
  33. Nadporozhskaya, M.A.; Mohren, G.M.J.; Chertov, O.G.; Komarov, A.S.; Mikhailov, A.V. Dynamics of Soil Organic Matter in Primary and Secondary Forest Succession on Sandy Soils in the Netherlands: An Application of the ROMUL Model. Ecol. Modell. 2006, 190, 399–418. [Google Scholar] [CrossRef]
  34. López-Castañeda, A.; Zavala-Cruz, J.; Palma-López, D.J.; Rincón-Ramírez, J.A.; Bautista, F. Digital Mapping of Soil Profile Properties for Precision Agriculture in Developing Countries. Agronomy 2022, 12, 353. [Google Scholar] [CrossRef]
  35. Pérez-Rodríguez, G.; López-Santos, A.; Velásquez-Valle, M.A.; Villanueva-Díaz, J.; García-Rodríguez, J.L. Spatial Distribution of Soil Organic Carbon by Digital Mapping: The Case of the Medio Aguanaval River Sub-Basin. Ing. Agrícola Biosist. 2021, 13, 227–245. [Google Scholar] [CrossRef]
  36. Guevara, M.; Federico Olmedo, G.; Stell, E.; Yigini, Y.; Aguilar Duarte, Y.; Arellano Hernández, C.; Arévalo, G.E.; Eduardo Arroyo-Cruz, C.; Bolivar, A.; Bunning, S.; et al. No Silver Bullet for Digital Soil Mapping: Country-Specific Soil Organic Carbon Estimates across Latin America. Soil 2018, 4, 173–193. [Google Scholar] [CrossRef] [Green Version]
  37. Gomes, L.C.; Faria, R.M.; de Souza, E.; Veloso, G.V.; Schaefer, C.E.G.R.; Filho, E.I.F. Modelling and Mapping Soil Organic Carbon Stocks in Brazil. Geoderma 2019, 340, 337–350. [Google Scholar] [CrossRef]
  38. Peri, P.; Maradei, D.; Lupi, A.; Vazquez, C.; Gyenge, J.; Gatica, M.; Sandoval, M.; Gaute, M. Estimación de Las Reservas de Carbono Orgánico Del Suelo Con Plantaciones Forestales y Otros Usos de La Tierra En Distintas Regiones de Argentina; INTA-Dirección Nacional de Desarrollo Foresto-Industrial del Ministerio de Agricultura Ganadería y Pesca: Buenos Aires, Argentina, 2022.
  39. Padilha, M.C.d.C.; Vicente, L.E.; Demattê, J.A.M.; dos Santos Wendriner Loebmann, D.G.; Vicente, A.K.; Salazar, D.F.U.; Guimarães, C.C.B. Using Landsat and Soil Clay Content to Map Soil Organic Carbon of Oxisols and Ultisols near São Paulo, Brazil. Geoderma Reg. 2020, 21, 253. [Google Scholar] [CrossRef]
  40. Zhang, H.; Wu, P.; Yin, A.; Yang, X.; Zhang, M.; Gao, C. Prediction of Soil Organic Carbon in an Intensively Managed Reclamation Zone of Eastern China: A Comparison of Multiple Linear Regressions and the Random Forest Model. Sci. Total Environ. 2017, 592, 704–713. [Google Scholar] [CrossRef] [PubMed]
  41. Duarte, E.; Zagal, E.; Barrera, J.A.; Dube, F.; Casco, F.; Hernández, A.J. Digital Mapping of Soil Organic Carbon Stocks in the Forest Lands of Dominican Republic. Eur. J. Remote Sens. 2022, 55, 213–231. [Google Scholar] [CrossRef]
  42. Martínez Pastur, G.; Aravena Acuña, M.C.; Silveira, E.M.O.; Von Müller, A.; La Manna, L.; González-Polo, M.; Chaves, J.E.; Cellini, J.M.; Lencinas, M.V.; Radeloff, V.C.; et al. Mapping Soil Organic Carbon Content in Patagonian Forests Based on Climate, Topography and Vegetation Metrics from Satellite Imagery. Remote Sens. 2022, 14, 5702. [Google Scholar] [CrossRef]
  43. FAO. Digital Soil Map of the World Version 3.6. Land and Water Development Division; FAO: Rome, Italy, 2003. [Google Scholar]
  44. CONAFOR. Inventario Estatal Forestal y de Suelos—Puebla 2013; Secretaria de Medio Ambiente y Recursos Naturales (SEMARNAT); Comisión Nacional Forestal (CONAFOR): Zapopan, Jalisco, México, 2014; p. 135. [Google Scholar]
  45. Walkley, A.; Black, I.A. An Examination of the Degtjareff Method for Determining Soil Organic Matter, and a Proposed Modification of the Chromic Acid Titration Method. Soil Sci. 1934, 37, 29–38. [Google Scholar] [CrossRef]
  46. SEMARNAT. NOM-021-RECNAT-2000; SEMARNAT: Ciudad de México, México, 2002. [Google Scholar]
  47. Young, N.E.; Anderson, R.S.; Chignell, S.M.; Vorster, A.G.; Lawrence, R.; Evangelista, P.H. A Survival Guide to Landsat Preprocessing. Ecology 2017, 98, 920–932. [Google Scholar] [CrossRef] [Green Version]
  48. Yadav, A.; Dodamani, B.M.; Dwarakish, G.S. Shoreline Analysis Using Landsat-8 Satellite Image. ISH J. Hydraul. Eng. 2021, 27, 347–355. [Google Scholar] [CrossRef]
  49. Comisión Nacional del Agua—Servicio Meteorológico Nacional. Red de Estaciones Climatológicas. Available online: https://smn.conagua.gob.mx/es/ (accessed on 15 December 2022).
  50. Khan, I.A.; Khan, M.R.; Baig, M.H.A.; Hussain, Z.; Hameed, N.; Khan, J.A. Assessment of Forest Cover and Carbon Stock Changes in Sub-Tropical Pine Forest of Azad Jammu & Kashmir (AJK), Pakistan Using Multitemporal Landsat Satellite Data and Field Inventory. PLoS ONE 2020, 15, e0226341. [Google Scholar] [CrossRef] [Green Version]
  51. ESRI. ArcGis Desktop: Versión 10.4; Environmental Systems Research Institute (ESRI): Redlands, CA, USA, 2016. [Google Scholar]
  52. Kogo, B.K.; Kumar, L.; Koech, R. Forest Cover Dynamics and Underlying Driving Forces Affecting Ecosystem Services in Western Kenya. Remote Sens. Appl. Soc. Environ. 2019, 14, 75–83. [Google Scholar] [CrossRef]
  53. Google Earth PRO; Google: Mountain View, CA, USA, 2022.
  54. Olofsson, P.; Foody, G.M.; Herold, M.; Stehman, S.V.; Woodcock, C.E.; Wulder, M.A. Good Practices for Estimating Area and Assessing Accuracy of Land Change. Remote Sens. Environ. 2014, 148, 42–57. [Google Scholar] [CrossRef]
  55. Branders, S.; Pereira, A.; Bernard, G.; Ernst, M.; Albert, A. Leveraging Historical Data for High-Dimensional Regression Adjustment, a Composite Covariate Approach. ArXiv 2021, arXiv:2103.14421. [Google Scholar]
  56. Wischmeier, W.H.; Smith, D.D. Predicting Rainfall Erosion Losses. A Guide to Conservation Planning; U.S Department of Agriculture. U.S. Government Printing Office: Washington, DC, USA, 1978; Volume 1. [Google Scholar]
  57. Rouse, J.W.; Haas, R.H.; Schell, J.A.; Deering, D.W.; Harlan, J.C. Monitoring the Vernal Advancements and Retrogradation of Natural Vegetation. In NASA/GSFC, Final Report, Greenbelt, MD, USA; Remote Sensing Center Texas A&M University: College Station, TA, USA, 1974; pp. 1–137. [Google Scholar]
  58. Baskan, O.; Dengiz, O.; Demirag, İ.T. The Land Productivity Dynamics Trend as a Tool for Land Degradation Assessment in a Dryland Ecosystem. Environ. Monit. Assess. 2017, 189, 1–21. [Google Scholar] [CrossRef]
  59. Cha, S.; Kim, C.; Kim, J.; Lee, A.L.; Park, K.; Koo, N.; Kim, Y.S.; Cha, S.; Kim, C.; Kim, J.; et al. Land-Use Changes and Practical Application of the Land Degradation Neutrality ( LDN ) Indicators : A Case Study in the Subalpine Forest Ecosystems. Republic of Korea. For. Sci. Technol. 2020, 16, 8–17. [Google Scholar] [CrossRef]
  60. Corradine, M.F.; Moreno, T. Plan de Ordenación y Manejo de La Cuenca Hidrográfica. In Actualización POMCA. Rio Garagoa; Consorcio Río Garagoa: Bogotá, Colombia, 2017; p. 131. [Google Scholar]
  61. ESRI. Cómo Funciona Pendiente. Available online: https://desktop.arcgis.com/es/arcmap/10.5/tools/3d-analyst-toolbox/how-slope-works.htm (accessed on 15 November 2021).
  62. Foster, G.R.; Meyer, L.D.; Onstad, C.A. Runoff Erosivity Factor and Variable Slope Length Exponents for Soil Loss Estimates. Trans. Am. Soc. Agric. Eng. 1977, 20, 683–687. [Google Scholar] [CrossRef]
  63. Desmet, P.J.J.; Govers, G. A GIS Procedure for Automatically Calculating the USLE LS Factor on Topographically Complex Landscape Units. J. Soil Water Conserv. 1996, 51, 427–433. [Google Scholar]
  64. McCool, D.K.; Brown, L.C.; Foster, G.R.; Mutchler, C.K.; Meyer, L.D. Revised Slope Steepness Factor for the Universal Soil Loss Equation. Trans. Am. Soc. Agric. Eng. 1987, 30, 1387–1396. [Google Scholar] [CrossRef]
  65. Pena, S.B.; Abreu, M.M.; Magalhães, M.R.; Cortez, N. Water Erosion Aspects of Land Degradation Neutrality to Landscape Planning Tools at National Scale. Geoderma 2020, 363, 114093. [Google Scholar] [CrossRef]
  66. Garcia Calabrese, M.V.; Paniagua, M.T.; Chicaiza, E. Análisis de La Distribución Espacial de La Precipitación Anual (1979–2014) Aplicando Métodos de Interpolación En La Región Occidental Del Paraguay. Rev. Geog. 2022, 164, 63–81. [Google Scholar] [CrossRef]
  67. Arnoldus, H.M.J.; Boodt, M.D.; Gabriels, D. An Approximation of the Rainfall Factor in the Universal Soil Loss Equation. Assess. Eros. 1980, 1, 127–132. [Google Scholar]
  68. Castelan Vega, R.; Tamariz Flores, V.; Linares Fleites, G.; Cruz Montalvo, A. Agresividad de Las Precipitaciones En La Subcuenca Del Río San Marcos, Puebla, México. Investig. Geog. 2015, 83, 28–40. [Google Scholar] [CrossRef] [Green Version]
  69. Royston, J. Algorithm AS 181: The W Test for Normality. J. Appl. Stat. 1982, 1, 176–180. [Google Scholar] [CrossRef]
  70. Cavanaugh, J.E.; Neath, A.A. The Akaike Information Criterion: Background, Derivation, Properties, Application, Interpretation, and Refinements. Wiley Interdiscip. Rev. Comput. Stat. 2019, 11, 1460. [Google Scholar] [CrossRef]
  71. Tziachris, P.; Aschonitis, V.; Chatzistathis, T.; Papadopoulou, M. Assessment of Spatial Hybrid Methods for Predicting Soil Organic Matter Using DEM Derivatives and Soil Parameters. Catena 2019, 174, 206–216. [Google Scholar] [CrossRef]
  72. Feng, C.; Feng, C.; Li, L.; Sadeghpour, A. A Comparison of Residual Diagnosis Tools for Diagnosing Regression Models for Count Data. BMC Med. Res. Methodol. 2020, 20, 1–21. [Google Scholar] [CrossRef] [PubMed]
  73. Kpienbaareh, D.; Mohammed, K.; Luginaah, I.; Wang, J.; Bezner Kerr, R.; Lupafya, E.; Dakishoni, L. Estimating Groundnut Yield in Smallholder Agriculture Systems Using PlanetScope Data. Land 2022, 11, 1752. [Google Scholar] [CrossRef]
  74. Kweon, G.; Lund, E.; Maxton, C. Soil Organic Matter and Cation-Exchange Capacity Sensing with on-the-Go Electrical Conductivity and Optical Sensors. Geoderma 2013, 199, 80–89. [Google Scholar] [CrossRef]
  75. Peón, J.; Recondo, C.; Fernández, S.; Calleja, J.F.; De Miguel, E.; Carretero, L. Prediction of Topsoil Organic Carbon Using Airborne and Satellite Hyperspectral Imagery. Remote Sens. 2017, 9, 1211. [Google Scholar] [CrossRef] [Green Version]
  76. Kiani-harchegani, M.; Hamidreza, S. Practicing Land Degradation Neutrality (LDN) Approach in the Shazand Watershed, Iran. Sci. Total Environ. 2020, 698, 134319. [Google Scholar] [CrossRef]
  77. Zhang, D.; Liu, G.; Hu, W. Mapping Tidal Flats with Landsat 8 Images and Google Earth Engine: A Case Study of the China’s Eastern Coastal Zone circa 2015. Remote Sens. 2019, 11, 924. [Google Scholar] [CrossRef] [Green Version]
  78. Arroyo, I.; Cervantes, V.; Tamaríz-Flores, V.; Castelán, R. Land Degradation Neutrality: State and Trend of Degradation at the Subnational Level in Mexico. Land 2022, 11, 562. [Google Scholar] [CrossRef]
  79. Flores Cesareo, J.C.; Bustamante González, A.; Vargas López, S.; Cajuste, L.; Escobedo, F.J.; Ramírez Valadez, M. Cartografía Del Uso Del Suelo En La Subcuenca Huaquechula, Puebla, México, Con Un Índice Combinado de Imágenes de Satélite. Investig. Geog. 2020, 101, 1–14. [Google Scholar] [CrossRef]
  80. Mas, J.-F.; Velázquez, A.; Couturier, S. La Evaluación de Los Cambios de Cobertura/Uso Del Suelo En La República Mexicana. Investig. Ambient. 2009, 1, 23–39. [Google Scholar]
  81. Leyva-Ovalle, Á.; Valdez-Lazalde, J.R.; de los Santos-Posadas, H.M.; Martínez-Trinidad, T.; Herrera-Corredor, J.A.; Lugo-Espinosa, O.; García-Nava, J.R. Monitoreo de La Degradación Forestal En México Con Base En El Inventario Nacional Forestal y de Suelos (Infys). Madera Bosques 2017, 23, 69–83. [Google Scholar] [CrossRef] [Green Version]
  82. Rosete, F.; Pérez, J.; Villalobos, M.; Navarro, E.; Salinas, E.; Remond, R. El Avance de La Deforestación En México 1976–2007. Madera Y Bosques 2014, 20, 21–35. [Google Scholar] [CrossRef] [Green Version]
  83. Vega, M.B.; Febles, J.M. La Agresividad de La Lluvia En Áreas Rurales de La Provincia La Habana Como Factor de Presión En La Sostenibilidad Agroambiental. In II Seminario Internacional de Cooperación y Desarrollo en Espacios Rurales Iberoamericanos; Ministerio de Educación y Ciencia, España: Almeria, Spain, 2008. [Google Scholar]
  84. An, J.Y.; Park, B.B.; Chun, J.H.; Osawa, A. Litterfall Production and Fine Root Dynamics in Cool-Temperate Forests. PLoS ONE 2017, 12, e0180126. [Google Scholar] [CrossRef] [Green Version]
  85. Gutiérrez-Flores, L.M.; Mauricio-Gutiérrez, A.; Carcaño-Montiel, M.G.; Portillo-Manzano, E.; Gómez-Velázquez, L.; Sánchez-Alonso, P.; López-Reyes, L. Fungi Associated with Sick Trees of Pinus Patula in Tetela de Ocampo, Puebla, Mexico. Arch. Phytopathol. Plant Prot. 2020, 53, 591–611. [Google Scholar] [CrossRef]
  86. Gutiérrez-Flores, L.M.; López-Reyes, L.; Hipólito-Romero, E.; Torres-Ramírez, E.; Castañeda-Roldán, E.I.; Mauricio-Gutiérrez, A. Biological Control Perspectives in the Pine Forest (Pinus Spp.), an Environmentally Friendly Alternative to the Use of Pesticides. Rev. Mex. Fitopatol. Mex. J. Phytopathol. 2022, 40, 1–24. [Google Scholar] [CrossRef]
  87. Ramirez-Huerta, M.; Juárez-Sánchez, J.P.; Ramírez-Valverde, B.; Martínez-Carrera, D.C.; Morales-Acoltzi, T. Adaptación Alimentaria de Campesinos Productores de Maíz Ante La Variabilidad Climática En El Centro Oriente Del Estado de Puebla, México. Estud. Soc. Rev. Aliment. Contemp. Desarro. Reg. 2021, 222, 1–30. [Google Scholar] [CrossRef]
  88. 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. Environ. Monit. Assess. 2006, 121, 595–611. [Google Scholar] [CrossRef] [PubMed]
  89. Mora-Vallejo, A.; Claessens, L.; Stoorvogel, J.; Heuvelink, G.B.M. Small Scale Digital Soil Mapping in Southeastern Kenya. Catena 2008, 76, 44–53. [Google Scholar] [CrossRef]
  90. Lamichhane, S.; Kumar, L.; Wilson, B. Digital Soil Mapping Algorithms and Covariates for Soil Organic Carbon Mapping and Their Implications: A Review. Geoderma 2019, 352, 395–413. [Google Scholar] [CrossRef]
  91. Mahmoudabadi, E.; Karimi, A.; Haghnia, G.H.; Sepehr, A. Digital Soil Mapping Using Remote Sensing Indices, Terrain Attributes, and Vegetation Features in the Rangelands of Northeastern Iran. Environ. Monit. Assess. 2017, 189, 500. [Google Scholar] [CrossRef] [PubMed]
  92. Wang, K.; Zhang, C.; Li, W. Comparison of Geographically Weighted Regression and Regression Kriging for Estimating the Spatial Distribution of Soil Organic Matter. GIScience Remote Sens. 2012, 49, 915–932. [Google Scholar] [CrossRef]
  93. Guo, P.T.; Wu, W.; Sheng, Q.K.; Li, M.F.; Liu, H.B.; Wang, Z.Y. Prediction of Soil Organic Matter Using Artificial Neural Network and Topographic Indicators in Hilly Areas. Nutr. Cycl. Agro. 2013, 95, 333–344. [Google Scholar] [CrossRef]
  94. Liu, Y.; Guo, L.; Jiang, Q.; Zhang, H.; Chen, Y. Comparing Geospatial Techniques to Predict SOC Stocks. Soil Tillage Res. 2015, 148, 46–58. [Google Scholar] [CrossRef]
  95. Wang, B.; Gray, J.M.; Waters, C.M.; Rajin Anwar, M.; Orgill, S.E.; Cowie, A.L.; Feng, P.; Li Liu, D. Modelling and Mapping Soil Organic Carbon Stocks under Future Climate Change in South-Eastern Australia. Geoderma 2022, 405, 115442. [Google Scholar] [CrossRef]
  96. Gallegos, Á.; Bautista, F. Soil Profile Photograph Dataset from Central Mexico to Delineate Horizons and Quantify Coarse Fragments. Data Br. 2022, 40, 7749. [Google Scholar] [CrossRef]
  97. Veum, K.S.; Goyne, K.W.; Kremer, R.J.; Miles, R.J.; Sudduth, K.A. Biological Indicators of Soil Quality and Soil Organic Matter Characteristics in an Agricultural Management Continuum. Biogeochemistry 2014, 117, 81–99. [Google Scholar] [CrossRef]
  98. Handal-Silva, A.; Pérez-Castresana, G.; Morán-Perales, J.; García-Suastegui, W. Historia de La Contaminación Hídrica Del Alto Balsas. Rev. Del Desarro. Urbano Sustentable 2017, 3, 10–23. [Google Scholar]
  99. Odhiambo, B.O.; Kenduiywo, B.K.; Were, K. Spatial Prediction and Mapping of Soil PH across a Tropical Afro-Montane Landscape. Appl. Geogr. 2020, 114, 102129. [Google Scholar] [CrossRef]
  100. Bouasria, A.; Ibno Namr, K.; Rahimi, A.; Ettachfini, E.M.; Rerhou, B. Evaluation of Landsat 8 Image Pansharpening in Estimating Soil Organic Matter Using Multiple Linear Regression and Artificial Neural Networks. Geo-Spatial Inf. Sci. 2022, 25, 353–364. [Google Scholar] [CrossRef]
Forests 14 00539 g001 550
Figure 1. Location of the study site, sampling areas and sampling points in Mexico. The color composite (Bands 5,4,3) of Landsat 8 from 7 October 2020.
Figure 1. Location of the study site, sampling areas and sampling points in Mexico. The color composite (Bands 5,4,3) of Landsat 8 from 7 October 2020.
Forests 14 00539 g001
Forests 14 00539 g002 550
Figure 2. (a) Forest/Non-forest cover spatial distribution in Ixtacamaxtitlan, Puebla, México; (b) NDVI in forest cover; (c) DEM; (d) slope; (e) precipitation; (f) temperature; (g) LS-factor; (h) R-factor.
Figure 2. (a) Forest/Non-forest cover spatial distribution in Ixtacamaxtitlan, Puebla, México; (b) NDVI in forest cover; (c) DEM; (d) slope; (e) precipitation; (f) temperature; (g) LS-factor; (h) R-factor.
Forests 14 00539 g002
Forests 14 00539 g003 550
Figure 3. Correlation analysis between SOM and environmental variables. * Correlation is significant at the 0.05 level.
Figure 3. Correlation analysis between SOM and environmental variables. * Correlation is significant at the 0.05 level.
Forests 14 00539 g003
Forests 14 00539 g004 550
Figure 4. Scatter plot of observed SOM vs. predicted SOM by the MLR model.
Figure 4. Scatter plot of observed SOM vs. predicted SOM by the MLR model.
Forests 14 00539 g004
Forests 14 00539 g005 550
Figure 5. Map of SOM content in forest areas of Ixtacamaxtitlan, Puebla, Mexico.
Figure 5. Map of SOM content in forest areas of Ixtacamaxtitlan, Puebla, Mexico.
Forests 14 00539 g005
Table
Table 1. Land cover types and reclassification to forest/non-forest cover.
Table 1. Land cover types and reclassification to forest/non-forest cover.
Land CoverDescriptionReclassification
ForestForest vegetation with cover > 40%.Forest cover
Secondary forestForest vegetation with cover < 40%.Forest cover
CroplandCrop areas and bare landNon-forest cover
GrasslandGrassed and shrub areasNon-forest cover
SettlementUrban areasNon-forest cover
Table
Table 2. Accuracy assessment for the forest and non-forest cover classification.
Table 2. Accuracy assessment for the forest and non-forest cover classification.
Forest CoverNon-Forest Cover
User’s accuracy (%)98.7597.00
Producer’s accuracy (%)95.6499.15
Rate(%)
Overall accuracy0.97797.7
Kappa0.95295.2
True positive rate0.98898.8
False positive rate0.0303.0
Table
Table 3. Summary of the environmental variables in forested areas.
Table 3. Summary of the environmental variables in forested areas.
VariableMinimumMaximumMeanSDCV (%)
NDVI0.010.450.280.0517.86
DEM (masl)2056.003506.002755.00297.6610.80
Slope (%)0.00201.0039.8519.4448.78
Precipitation (mm)589.00993.00749.0081.4510.87
Temperature (°C)9.8019.2014.502.8019.31
LS-factor0.02309.927.023.9456.12
R-factor (MJ mm ha−1 yr−1)1856.803965.803039.90576.1018.95
SD = Standard Deviation; CV = Coefficient of variation.
Table
Table 4. Descriptive statistics of SOM and environmental variable at the sampling points.
Table 4. Descriptive statistics of SOM and environmental variable at the sampling points.
VariableMinimumMaximumMeanSDCV (%)
SOM (%)1.1010.285.032.8055.66
NDVI0.140.300.230.0521.73
DEM (masl)2578.002808.00267088.993.33
Slope (%)11.44101.1931.9320.4363.98
Precipitation (mm)720.30839.40773.639.705.13
Temperature (°C)11.5816.2214.261.9713.81
LS-factor0.5514.255.323.3062.03
R-factor (MJ mm ha−1 yr−1)2801.003587.003088.00322.4310.44
SD = Standard Deviation; CV = Coefficient of variation.
Table
Table 5. Performance of the MLR model and LOOCV for SOM.
Table 5. Performance of the MLR model and LOOCV for SOM.
FunctionMLRLOOCV
R2RMSEAICR2cvRMSEcvRPD
SOM = −36.96 + NDVI [22.20] + LS-factor [0.2839] + Precipitación [0.0456]0.781.2783.1170.691.531.83
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Arroyo, I.; Tamaríz-Flores, V.; Castelán, R. Mapping Forest Cover and Estimating Soil Organic Matter by GIS-Data and an Empirical Model at the Subnational Level in Mexico. Forests 2023, 14, 539. https://doi.org/10.3390/f14030539

AMA Style

Arroyo I, Tamaríz-Flores V, Castelán R. Mapping Forest Cover and Estimating Soil Organic Matter by GIS-Data and an Empirical Model at the Subnational Level in Mexico. Forests. 2023; 14(3):539. https://doi.org/10.3390/f14030539

Chicago/Turabian Style

Arroyo, Itzel, Víctor Tamaríz-Flores, and Rosalía Castelán. 2023. "Mapping Forest Cover and Estimating Soil Organic Matter by GIS-Data and an Empirical Model at the Subnational Level in Mexico" Forests 14, no. 3: 539. https://doi.org/10.3390/f14030539

APA Style

Arroyo, I., Tamaríz-Flores, V., & Castelán, R. (2023). Mapping Forest Cover and Estimating Soil Organic Matter by GIS-Data and an Empirical Model at the Subnational Level in Mexico. Forests, 14(3), 539. https://doi.org/10.3390/f14030539

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop