Machine Learning for Shield-Scale Soil Health Mapping: A Systematic Review

2.1 Research design

This study employed a systematic literature review (SLR) design following the Preferred Reporting Items for Systematic Reviews and Meta-Analyses (PRISMA) 2020 guidelines (Page et al., 2021). The review protocol was developed a priori and documented internally. The review was not registered in PROSPERO as the protocol pertains to a methodological review without health-related outcomes. The PRISMA 2020 checklist was completed to ensure reporting transparency.

2.2 Eligibility criteria

Eligibility criteria were defined based on the PICOC framework and operationalized for consistent application. Two reviewers independently applied the criteria; disagreements were resolved through discussion and consensus. Table 1 presents the complete inclusion and exclusion criteria.

2.3 Search strategy and information sources

We systematically searched six electronic databases: Scopus, Springer, ScienceDirect, Nature, IEEE Xplore, and Wiley. The search was conducted between January and March 2026, covering publications from January 2016 through March 2026. Hand-searching of reference lists supplemented the database searches.

The search strategy combined four thematic blocks using Boolean operators. The primary search string for Scopus (TITLE-ABS-KEY) was:

Final Boolean String: (Block 1) AND (Block 2) AND (Block 3) AND (Block 4).

2.4 Study selection process

The study selection process followed the PRISMA 2020 flow diagram, comprising four stages: identification, screening, eligibility, and inclusion. All records were exported to EndNote X9 for duplicate removal. Two reviewers independently screened titles and abstracts using Rayyan (Ouzzani et al., 2016). Inter-reviewer agreement was κ = 0.83. Full-text articles were retrieved and assessed independently. Reasons for exclusion were documented. Studies that satisfied all inclusion criteria were included in the qualitative and quantitative synthesis.

2.5 Data extraction

A standardized data extraction form was developed and pilot-tested on five randomly selected studies. Two reviewers independently extracted data; discrepancies were resolved through consensus. For multi-study papers, all relevant model runs were extracted separately. Data extraction was performed using Microsoft Excel with validation rules. Authors were contacted for missing data where necessary (n = 4 studies responded).

2.6 Data items

The following data items were extracted:

2.7 Quality assessment

Quality assessment used a modified JBI Critical Appraisal Checklist with seven criteria (Munn et al., 2015):

Scoring: ≥6 Yes = High quality; 4–5 Yes = Medium quality; ≤3 Yes = Low quality. Two reviewers independently assessed each study. Results are presented in Table 2.

2.8 Data synthesis

This review employed narrative synthesis as the primary approach for data integration and interpretation, structured around the four research questions to ensure alignment between extracted evidence and review objectives (Popay et al., 2006). The synthesis process involved three interrelated stages. First, in the preliminary synthesis stage, we organized extracted data from all 46 included studies into descriptive summaries grouped by ML algorithms, data sources, target variables, and geographic locations to identify patterns and frequencies across the literature. Second, during the exploration of relationships stage, we examined how contextual factors such as sample size, data source, validation method, and geographic location influenced study outcomes and performance variability. Third, the assessment of robustness stage involved evaluating methodological quality, consistency of findings, and potential for publication bias to inform confidence in conclusions drawn.

For quantitative synthesis, we extracted reported R² values and calculated descriptive statistics including means, medians, and ranges grouped by algorithm type and target property. Boxplots were generated to visualize R² distributions across algorithms, with Figure 3 presented and discussed in Section 3.3. Where studies compared ML with conventional methods, we calculated percentage improvement using the formula: Improvement (%) = [(ML_metric - Conventional_metric) /Conventional_metric] × 100. Subgroup analyses were conducted by target property, data source, validation method, and sample size. Geographic analysis was performed by geocoding study locations and presenting the distribution in Figure 4, which is presented and discussed in Section 3.2. The narrative synthesis is presented through text, tables, and figures to facilitate understanding and transparency.

3.1 Selection results

The systematic search identified 921 records from six databases (Scopus: 267, Springer: 184, ScienceDirect: 156, Nature: 98, IEEE Xplore: 124, Wiley: 92). After removing duplicates, 821 records underwent title and abstract screening. Of these, 316 records were excluded based on title and abstract screening. The primary reasons for exclusion were: publication year outside 2015–2026 (89 records), non-English language (42 records), book chapters, retracted articles, short surveys, notes, editorials, data papers, or other non-peer-reviewed documents (97 records), and articles not relevant to ML or soil mapping (88 records). The remaining 505 full-text articles were assessed for eligibility.

At full-text screening, 459 articles were excluded. The reasons were: studies conducted outside agricultural field-scale context or without clear field-scale delineation (<100 Ha) (178 studies); studies that did not substantially address ML algorithms or digital soil mapping (96 studies); unclear methodology or inaccessible full text (67 studies); and studies not relevant to the research questions, such as those only focusing on soil property measurement without spatial prediction or mapping, or using conventional statistics without ML comparison (118 studies). The remaining 46 studies met all inclusion criteria and were included in the qualitative and quantitative synthesis, as shown in Figure 1.

Figure 1. PRISma flow diagram.

Source: Researcher’s Visualisation (2026).

3.2 Geographic distribution

China dominated the geographic representation, contributing 24 studies (52.2%), with most studies conducted in temperate and subtropical agricultural regions such as the North China Plain, Sichuan Basin, and Yangtze River Delta. Studies from China include L. Wang et al. (2026), Guo et al. (2026), Bao et al. (2025), L. Wang & Zhou (2023), Z. Wang et al. (2025b), F. Zhou et al. (2025), Ning et al. (2026), He et al. (2024), Deng et al. (2025), Luo et al. (2026), Xiang et al. (2025), Zheng et al. (2023), J. Zhou et al. (2023), R. Wang et al. (2025a), N. Chen et al. (2026), (Li et al., 2026), Qiu et al. (2025), X. Zhang et al. (2025), Wei et al. (2025), Wu et al. (2022), Huang et al. (2024), Dong et al. (2025), Han et al. (2025), K. Wang et al. (2021) and W. chun Zhang et al. (2022).

The remaining studies were distributed across multiple countries: the United States contributed 3 studies (J. Chen et al., 2025; Xu et al., 2025); Turkey contributed 3 studies (Abiye & Dengiz, 2025; Acir, 2025; Kaya et al., 2023), the latter as Turkey/India collaboration; Peru contributed 3 studies (Carbajal-Llosa et al., 2025; Enriquez et al., 2025; Pizarro et al., 2023); Spain contributed 2 studies (El-Jamaoui, Delgado-Iniesta, et al., 2025a; El-Jamaoui, José Martínez Sánchez, et al., 2025b); Canada contributed 2 studies (Naeimi et al., 2025, 2026); the Netherlands contributed 2 studies (Schmidinger et al., 2026; Vavlas et al., 2024) as Germany/Netherlands collaboration; India contributed 2 studies (C et al., 2026; Kumaraperumal et al., 2022); while Germany (Schmidinger et al., 2026), Russia (Suleymanov et al., 2026), Belgium (De Bièvre et al., 2026), and Iran (Hosseinpour-Zarnaq et al., 2025) each contributed 1 study. Notably, only 8 studies (17.4%) were conducted in tropical regions, including Peru (Carbajal-Llosa et al., 2025; Enriquez et al., 2025; Pizarro et al., 2023), India (C et al., 2026; Kumaraperumal et al., 2022), Turkey (Abiye & Dengiz, 2025; Kaya et al., 2023), and Iran (Hosseinpour-Zarnaq et al., 2025). A striking observation was the complete absence of studies from Southeast Asian countries, including Indonesia, Thailand, Vietnam, Malaysia, and the Philippines, despite their significant agricultural importance and vulnerability to soil degradation. Figure 4 illustrates this geographic concentration and the substantial gaps in tropical region representation, particularly in Southeast Asia.

3.3 Study scale and sample sizes

The study areas ranged from 40 Ha in J. Zhou et al. (2023) to 16,258 Ha in De Bièvre et al. (2026), with a median area of 198 Ha and interquartile range of 103 to 397 Ha. Studies with smaller areas included El-Jamaoui, Delgado-Iniesta, et al. (2025a) with 76 Ha, Abiye & Dengiz (2025) with 50 Ha, and Pizarro et al. (2023) with 46 Ha, all of which utilized UAV data to achieve high-resolution mapping. The sample sizes, representing the number of soil sampling points, exhibited considerable variation from 40 samples in J. Zhou et al. (2023) to 16,258 samples in De Bièvre et al. (2026), with a median of 214 samples and interquartile range of 144 to 500. Studies with large sample sizes included J. Chen et al. (2025) with 98,211 samples, De Bièvre et al. (2026) with 16,258 samples, and Kumaraperumal et al. (2022) with 4,479 samples. Notably, 12 studies (26.1%) had sample sizes below 100, including J. Zhou et al. (2023) with 40 samples, Lyu et al. (2025) with 45 samples, Pizarro et al. (2023) with 46 samples, Abiye & Dengiz (2025) with 50 samples, Han et al. (2025) with 60 samples, He et al. (2024) with 62 samples, Acir (2025) with 64 samples, El-Jamaoui, José Martínez Sánchez, et al. (2025b) with 70 samples, El-Jamaoui, Delgado-Iniesta, et al. (2025a) with 76 samples, Luo et al. (2026) with 82 samples, Zheng et al. (2023) with 100 samples, and Enriquez et al. (2025) with 103 samples. This raises concerns about model overfitting and limited generalizability for these small-sample studies. Table 3 provides a comprehensive summary of all 46 included studies.

3.4 Quality assessment results

The quality assessment revealed that 28 studies (60.9%) were rated as high quality, 15 studies (32.6%) were rated as medium quality, and 3 studies (6.5%) were rated as low quality. The most common weakness identified across all studies was Q5 (spatial validation), with only 14 studies (30.4%) implementing spatial cross-validation, geostatistical-based validation, or geographically independent test sets. Studies that employed spatial validation included Guo et al. (2026), J. Chen et al. (2025), F. Zhou et al. (2025), Carbajal-Llosa et al. (2025), Schmidinger et al. (2026), De Bièvre et al. (2026), Xu et al. (2025), and W. chun Zhang et al. (2022), among others. All 46 studies clearly stated their objectives (Q1), used appropriate study designs (Q2), and reported adequate performance metrics including R² and RMSE (Q6). However, the low-quality studies (n = 3) failed to report reproducible methodology or adequate validation metrics, which compromised their overall quality scores. The predominance of high-quality studies (60.9%) indicates that the included literature generally meets acceptable methodological standards, though the spatial validation gap remains a significant concern.

3.5 Dominant algorithms and data sources

The frequency of ML algorithms across the 46 studies is shown in Figure 2. Random Forest was the most common (34 studies, 73.9%), used in (El-Jamaoui, Delgado-Iniesta, et al., 2025a; El-Jamaoui, José Martínez Sánchez, et al., 2025b), L. Wang & Zhou (2023), Kumaraperumal et al. (2022), Enriquez et al. (2025), F. Zhou et al. (2025), He et al. (2024), Deng et al. (2025), C et al. (2026), Qiu et al. (2025), X. Zhang et al. (2025), Wu et al. (2022), Suleymanov et al. (2026), Han et al. (2025), K. Wang et al. (2021), Hosseinpour-Zarnaq et al. (2025), Kaya et al. (2023), and many others. SVM/SVR appeared in 23 studies (50.0%), XGBoost in 14 (30.4%), ANN in 9 (19.6%), CNN in 5 (10.9%), Cubist in 8 (17.4%), PLSR in 6 (13.0%), and LSTM in 2 (4.3%). RF was included as a benchmark in 80% of comparative studies. Common combinations were RF + SVM (n = 15), RF + XGBoost (n = 10), RF + ANN (n = 7), and RF + Cubist (n = 6). Ensemble methods (RF, XGBoost, hybrids) appeared in 37 studies (80%). Temporal trends show increasing XGBoost adoption from 2023 onward, and all CNN studies were published in 2025–2026, indicating a shift toward more complex models.

Figure 2. Distribution of machine learning algorithms.

A comparison of remote sensing platforms is given in Table 5. Sentinel-2 was the most frequently used platform (25 studies, 54.3%), including El-Jamaoui, José Martínez Sánchez, et al. (2025b), Bao et al. (2025), L. Wang & Zhou (2023), F. Zhou et al. (2025), He et al. (2024), N. Chen et al. (2026), De Bièvre et al. (2026), Suleymanov et al. (2026), K. Wang et al. (2021), and others. Landsat 8 was used in 18 studies (39.1%), e.g., Kaya et al. (2023), Deng et al. (2025), Carbajal-Llosa et al. (2025), Qiu et al. (2025), Wei et al. (2025), Huang et al. (2024). UAV/drone platforms were employed in 10 studies (21.7%), including El-Jamaoui, José Martínez Sánchez, et al. (2025b) & El-Jamaoui, Delgado-Iniesta, et al. (2025a), Bao et al. (2025), Vavlas et al. (2024), Enriquez et al. (2025), J. Zhou et al. (2023), Pizarro et al. (2023), and Lyu et al. (2025). MODIS appeared in 5 studies (10.9%) – J. Chen et al. (2025), Ning et al. (2026), W. chun Zhang et al. (2022) – and hyperspectral sensors in 6 (13.0%) – Bao et al. (2025), Hosseinpour-Zarnaq et al. (2025), Song et al. (2018), Han et al. (2025), Luo et al. (2026). PlanetScope was used in two studies (Lyu et al., 2025; Z. Wang et al., 2025b). UAV-based studies yielded the highest R² (0.80–0.92), followed by hyperspectral (0.70–0.95), Sentinel-2 (0.56–0.85), and Landsat 8 (0.41–0.86). Environmental covariates (DEM, climate, legacy soil) were integrated in 14, 6, and 5 studies respectively. Data fusion (multi-source) was practiced in 28 studies (60.9%), with common combinations being Sentinel-2 + Landsat, Sentinel-2 + UAV, and satellite+DEM + climate.

Performance by soil property is summarized in Table 4. SOC/SOM were the most frequent targets (33 studies, 71.7%), with R² ranging from 0.41 to 0.95 and RMSE from 0.20 to 4.25 g/kg. The highest SOC accuracies came from UAV-based RF models (El-Jamaoui, Delgado-Iniesta, et al., 2025a) R² = 0.92; (J. Zhou et al., 2023) R² = 0.91. pH was studied in 10 studies (R² 0.45–0.79), with best results from weighted ensemble (Carbajal-Llosa et al., 2025) >0.93 and XGBoost (Z. Wang et al., 2025b) 0.72. Total nitrogen (TN/STN) was predicted in 8 studies (R² 0.56–0.95), with (He et al., 2024) reaching 0.95 and (Luo et al., 2026) 0.77. Phosphorus and potassium (each in 4–5 studies) showed moderate R² (0.42–0.89 and 0.51–0.87). Salinity/EC (4 studies) achieved R² 0.71–0.88 using MB-CNN (Dong et al., 2025) and FOPW-EL (L. Wang et al., 2026). CEC was the most challenging R² 0.41–0.56(Kaya et al., 2023; Kumaraperumal et al., 2022). Texture (3 studies) ranged from 0.42 to 0.82 (Hosseinpour-Zarnaq et al., 2025; Qiu et al., 2025). SOM (5 studies) and soil moisture (2 studies) also showed good performance with RF and ADF respectively.

3.6 Ensemble vs. Single classifier performance

The boxplot comparison of R² distributions between ensemble and single-classifier models is presented in Figure 3. Ensemble methods (RF, XGBoost, Cubist, hybrids) had a median R² of 0.75 (IQR 0.62–0.85), while single classifiers (SVM, ANN, PLSR) had a median of 0.70 (IQR 0.58–0.82). The difference was statistically significant (Wilcoxon, p < 0.001), with ensembles outperforming singles by 0.05–0.10 R². The boxplot also shows narrower interquartile ranges for ensembles, indicating greater consistency. An outlier (RF, R² = 0.95) from He et al. (2024) illustrates the potential of ensemble methods under favorable conditions.

Figure 3. Performance comparison of ML algorithms (R² Boxplot).

Source: Researcher’s Process (2026).

Figure 4. Geographic distribution of included studies.

Source: Researcher’s Process (2026).

Across 83 model runs, RF achieved mean R² = 0.74 (SD 0.15, median 0.75, range 0.41–0.95), XGBoost 0.72 (0.14, 0.73, 0.57–0.87), SVM 0.68 (0.16, 0.69, 0.42–0.92), ANN 0.70 (0.18, 0.73, 0.40–0.95), CNN 0.68 (0.12, 0.67, 0.55–0.80), Cubist 0.70 (0.14, 0.72, 0.47–0.85), and ensemble hybrids 0.76 (0.13, 0.78, 0.55–0.88). RF showed the lowest variability, confirming its robustness. ANN had the widest range, indicating sensitivity to data and tuning. Ensemble hybrids achieved the highest mean R², suggesting added value from combining diverse learners.

For SOC/SOM, ensemble methods (RF/XGBoost) averaged R² = 0.76 (n = 28) versus 0.70 for singles (SVM/ANN, n = 20), a + 0.06 advantage. For pH, ensembles averaged 0.87 (n = 6) vs 0.82 (n = 4), +0.05. For total N, 0.82 (n = 5) vs 0.78 (n = 3), +0.04. For salinity, 0.79 (n = 4) vs 0.74 (n = 3), +0.05. The ensemble advantage was largest for SOC, the most complex property, and smaller for pH and nitrogen, suggesting that non-linearity and interaction complexity moderate the benefit of ensemble methods.

RF was the best algorithm in 22 studies (47.8%), including El-Jamaoui, Delgado-Iniesta, et al. (2025a) and El-Jamaoui, José Martínez Sánchez, et al. (2025b), L. Wang & Zhou (2023), Kumaraperumal et al. (2022), Enriquez et al. (2025), F. Zhou et al. (2025), He et al. (2024), Deng et al. (2025), C et al. (2026), Qiu et al. (2025), X. Zhang et al. (2025), Wu et al. (2022), Han et al. (2025), K. Wang et al. (2021), and others. XGBoost was best in 8 studies (17.4%) – Z. Wang et al. (2025b), Xiang et al. (2025), Wei et al. (2025), Huang et al. (2024). Ensemble hybrids were best in 6 studies (13.0%) – Carbajal-Llosa et al. (2025), Lyu et al. (2025), Bao et al. (2025). SVM best in 5 (10.9%) – Abiye & Dengiz (2025), Hosseinpour-Zarnaq et al. (2025). ANN best in 3 (6.5%) – Acir (2025), Kaya et al. (2023). CNN best in 2 (4.3%) – Ning et al. (2026), Dong et al. (2025). In 14 studies (30.4%), RF was among the top two even when not the absolute best, confirming its reliability as a benchmark.

3.7 ML vs. Conventional geostatistics

Eight studies (17.4%) directly compared ML with conventional methods (OK, IDW, SLR, PLSR). These include J. Zhou et al. (2023), K. Wang et al. (2021), Song et al. (2018), Wu et al. (2022), Qiu et al. (2025), W. chun Zhang et al. (2022), L. Wang & Zhou (2023), and El-Jamaoui, José Martínez Sánchez, et al. (2025b). Conventional methods included OK (5 studies), IDW (2), SLR (3), and PLSR (2).

ML improved R² by 11% to 57% (mean 35.3%; weighted by sample size 32.6%). The largest gains were seen in W. chun Zhang et al. (2022) with RFRK over OK (+57% NSE), Song et al. (2018) with BPNNOK over OK (+56%), and Wu et al. (2022) with Cubist over SLR (+52%). K. Wang et al. (2021) reported +34%, L. Wang & Zhou (2023) +29%, Qiu et al. (2025) +26%, El-Jamaoui, José Martínez Sánchez, et al. (2025b) +18%, and J. Zhou et al. (2023) +11%. RMSE reductions ranged from 12% to 43% (mean 24.7%), with SOC showing the largest reductions (~31%) and texture the smallest (~18%). These improvements are practically meaningful; for example, K. Wang et al. (2021) reduced SOC RMSE from ~6.2 g/kg (PLSR) to ~3.9 g/kg (RF).

Sample size strongly modulated improvement: studies with n > 300 showed >30% R² gain, while n < 100 showed ≤15% gain. Spatial autocorrelation strength also mattered; in low-relief plains (Wu et al., 2022), kriging already captured spatial structure, reducing ML’s edge. Data richness (multi-source integration) led to larger improvements (>40%) compared to spectral-only data (≤25%). Target variable also influenced results: SOC and pH improved most, while texture and CEC improved least.

Conventional methods remained competitive in three scenarios: (i) very small samples (n < 50), as in J. Zhou et al. (2023) with only 11% improvement; (ii) low-relief, homogeneous landscapes with strong spatial autocorrelation; and (iii) when interpretability and simplicity are paramount, as shown by Acir (2025) using Lasso with R² = 0.75.

3.8 Methodological challenges and research gaps

As noted in Table 2, only 14 studies (30.4%) used spatial cross-validation. The remaining 32 used random CV, which likely inflates R² by 0.04–0.05 (observed mean R²: spatial CV 0.71, random CV 0.75). This gap is critical because random CV underestimates prediction error in spatially autocorrelated data, leading to over-optimistic performance claims. Twelve studies (26.1%) had n < 100, and these reported higher mean R² (0.76) than studies with n > 300 (0.71), indicating possible overfitting or publication bias. Only 28% of studies achieved sampling density ≥ 1 sample per 5 ha, considered adequate for field-scale mapping. The median density was 1 sample per 8.5 ha, which may be insufficient for heterogeneous tropical soils.

Only 6 studies (13.0%) used XAI techniques: SHAP in Guo et al. (2026), Xiang et al. (2025), and Acir (2025); LIME in one study; and feature importance plots in two others. The remaining 40 studies treated ML as black boxes, limiting trust and practical adoption, especially in farmer-centric applications. China’s dominance (53.2%) and the severe underrepresentation of tropical regions (only 17.4%) – with zero studies from Southeast Asia – constitute a major geographic bias. Tropical soils have different mineralogy, spectral properties, and cloud cover regimes, so models developed in temperate China may not transfer. Cloud-resistant methods (SAR/radar) were absent in tropical studies, and bare-soil composites were used in only two tropical studies.

Only 3 studies (6.5%) used time-series ML (J. Chen et al., 2025; L. Wang & Zhou, 2023; F. Zhou et al., 2025). No study tested model transferability to geographically independent fields, which is essential for operational deployment. Hyperparameter tuning was inadequately reported in 43% of studies. Data preprocessing was often superficially described. Uncertainty quantification (prediction intervals) was missing in 89% of studies. A synthesis of all research gaps and proposed future directions is provided in Table 6.

4.1 Dominance of random forest and data fusion

The overwhelming dominance of Random Forest (74% of studies) corroborates earlier reviews (Hengl et al., 2017; Padarian et al., 2020). RF’s practical advantages – handling high-dimensional spectral data, robustness to noise and missing values, simple hyperparameter tuning, built-in feature importance, scalability, and mature implementations – make it the go-to algorithm for field-scale soil mapping. Studies like El-Jamaoui, Delgado-Iniesta, et al. (2025a) achieved R² = 0.92 using RF with UAV hyperspectral data, and K. Wang et al. (2021) reached R² = 0.86 with Sentinel-2A. The growing use of XGBoost (30% of studies, all post-2023) reflects its strength in handling large datasets and its regularization capabilities, though it requires more careful tuning.

The dominance of Sentinel-2 (54%) over Landsat (39%) is due to its finer spatial resolution and 5-day revisit, but Landsat’s long archive remains valuable. Data fusion (61% of studies) is widely recognized as essential for capturing the full set of soil-forming factors (McBratney et al., 2003). Multi-scale integrations, such as Bao et al. (2025) combining UAV, ZY1-02D, and Sentinel-2 (R² = 0.75), and Qiu et al. (2025) integrating Landsat 8 and DEM (R² = 0.82 for silt), exemplify the power of fusing complementary data. However, simple concatenation remains the dominant fusion method, and more sophisticated approaches (wavelet, deep fusion) are rarely explored.

4.2 Ensemble advantage and model selection

The consistent 0.05–0.10 R² advantage of ensemble methods over single classifiers can be attributed to variance reduction via bagging (RF), sequential error correction via boosting (XGBoost), and automatic capture of non-linear interactions. Carbajal-Llosa et al. (2025) showed that weighted ensembles of RF, SVM, ANN, and XGBoost can exceed R² = 0.93 for pH and EC. The high variability of ANN performance (0.40–0.95) indicates sensitivity to sample size and tuning; Li et al. (2026) with n = 995 achieved 0.76 for SOM, but small-sample ANN studies often overfit. CNN applications (Dong et al., 2025; Ning et al., 2026) are promising but remain data-hungry. For typical field-scale studies with n = 100–500, RF offers the best balance of accuracy, robustness, interpretability, and ease of use. XGBoost is recommended when n > 500, and deep learning only when n > 1000 or when raw image input is used.

4.3 ML vs. Conventional geostatistics

The 11–57% R² improvement of ML over kriging/regression is consistent with the literature (Padarian et al., 2020; Wadoux et al., 2020). ML outperforms because it integrates auxiliary covariates (spectral, topographic, climatic) that capture soil-forming factors beyond spatial proximity, and it handles non-linearities and interactions. Hybrid approaches (RF-RK, BPNNOK, RFRK) as in Song et al. (2018), Qiu et al. (2025), and W. chun Zhang et al. (2022) are particularly effective because they combine ML’s predictive power with kriging’s ability to model residual spatial autocorrelation. However, ML’s advantage is context-dependent: small samples (n < 50), strong spatial autocorrelation, and simple linear relationships reduce its edge. In such cases, conventional methods may be nearly as accurate and more interpretable a point underscored by Acir (2025) using Lasso for soil fertility index with R² = 0.75.

4.4 Methodological gaps and their implications

The spatial validation deficit (70% of studies lacking spatial CV) is the most critical flaw. Random CV overestimates performance because nearby points are not independent (Wadoux et al., 2020). The observed 0.04–0.05 R² gap between spatial and random CV in our dataset likely understates the true inflation. This over-optimism undermines trust and operational reliability models that report R² = 0.80 may fail in new fields. The geographic bias toward China (53%) and the absence of Southeast Asian studies mean that most models are developed in temperate systems and may not transfer to tropical conditions (Minasny & McBratney, 2016). Tropical soils differ in mineralogy, clay types, organic matter dynamics, and spectral behavior, and cloud cover limits optical data. The few tropical studies (Carbajal-Llosa et al., 2025; Enriquez et al., 2025) show that UAV and ensemble methods work well, but they are too scarce. The interpretability gap (only 13% XAI) is another major barrier; without SHAP or LIME, farmers and agronomists cannot understand predictions, limiting adoption. Finally, the lack of transferability testing and the scarcity of uncertainty reporting further constrain practical use.

4.5 Comparison with existing literature

Our findings align with Padarian et al. (2020) in terms of RF dominance and increasing ensemble use, but we extend the timeline to 2026 and confirm that XGBoost has emerged as a strong competitor. Compared to global mapping efforts (Hengl et al., 2017), field-scale studies achieve higher R² (0.70–0.85 vs. 0.40–0.60) due to more homogeneous conditions and better local data. Unlike Jenčo et al. (2020), who found PLSR common, our review shows a decisive shift toward ML, with PLSR now used mainly as a baseline. The spatial validation gap, noted earlier by Wadoux et al. (2020), persists at 70%, indicating that awareness has not translated into widespread practice.

4.6 Implications for practice and policy

For researchers, we recommend mandatory spatial cross-validation (spatial block CV) and reporting both random and spatial CV results to quantify over-optimism. Minimum reporting standards should include sampling design, preprocessing steps, hyperparameters, feature importance, and uncertainty metrics. For practitioners, RF remains the default choice for n = 100–500; XGBoost for n > 500; deep learning only for large datasets. Sentinel-2 offers the best cost-performance trade-off for most applications. For policymakers, investment in tropical soil spectral libraries, open data infrastructure, capacity building, and cloud-resistant methods (SAR/radar) is urgent. Southeast Asia should be a priority region for establishing soil health monitoring networks.

4.7 Research agenda for tropical field-scale soil mapping

We prioritize: (1) creation of open tropical soil spectral libraries and benchmark datasets; (2) development of standardized spatial validation protocols tailored to smallholder systems; (3) co-design of XAI tools with farmers and extension agents; (4) optical-SAR fusion for cloud-resistant mapping; (5) transfer learning and few-shot learning to reduce ground-truth needs; (6) smartphone-based proximal sensing (Naeimi et al., 2025, 2026) coupled with edge ML; and (7) multi-property, integrated soil health frameworks that link to crop yield and ecosystem services.

4.8 Limitations of this review

This review has potential publication bias (high-R² studies overrepresented), geographic bias (Chinese dominance), language bias (English only), and database coverage limitations. Heterogeneity in study designs, target variables, and metrics complicates direct meta-analysis. Our modified JBI checklist has not been validated for ML studies. Despite these, the review’s strengths include its field-scale focus, systematic rigor, comprehensive data extraction, quantitative synthesis, and identification of actionable gaps.