Predicting soil chemical properties such as soil organic carbon (SOC) and available phosphorus (Ava-P) content is critical in areas where different land uses exist. The distribution of SOC and Ava-P is influenced by both natural and anthropogenic factors. This study aimed at (1) predicting SOC and Ava-P in a piedmont plain of Northeast Iran using the Random Forests (RF) and Cubist mathematical models and hybrid models (Regression Kriging), (2) comparing the models’ results, and (3) identifying the key variables that influence the spatial dynamics of soil SOC and Ava-P under different agricultural practices. The machine learning models were trained with 201 composite surface soil samples and 24 ancillary data, including climate (C), organism (O), topography- relief (R), parent material (P) and key soil features (S) according to the SCORPAN digital soil mapping framework, which can predictively represent soil formation factors spatially. Clay, one of the most critical soil properties with a well-known relationship to SOC, was the most important predictor of SOC, followed by open-access multispectral satellite images-based vegetation and soil indices. Ava-P had a similar set of effective variables. Hybrid approaches did not improve model accuracy significantly, but they did reduce map uncertainty. In the validation set, Ava-P was calculated using the RF algorithm with a normalized root mean square (NRMSE) of 96.8, while SOC was calculated using the Cubist algorithm with an NRMSE of 94.2. These values did not change when using the hybrid technique for Ava-P; however, they changed just by 1% for SOC. The management of SOC content and the supply of Ava-P in agricultural activities can be guided by SOC and Ava-P digital distribution maps. Produced digital maps in which the soil scientist plays an active role can be used to identify areas where concentrations are high and need to be protected, where uncertainty is high and sampling is required for further monitoring.

Soil mapping is essential for sustainable land management, precision agriculture, and environmental monitoring. This synoptic review evaluates the evolution of classical and modern geostatistical methods, spanning 2000 to 2024, and their integration with machine learning (ML) and remote sensing (RS) technologies. Key techniques such as spatial autocorrelation, variogram analysis, kriging and its variants, Bayesian models, Markov Chains, geostatistical simulations, and ensemble modeling were assessed in the context of digital soil mapping (DSM). Emphasis was placed on hybrid approaches that fuse geostatistics with ML algorithms including Random Forest (RF), Extreme Gradient Boosting (XGBoost), Support Vector Machines (SVM), and Artificial Neural Networks (ANN), along with the enrichment of spatial models using RS data. Hybrid approaches combining geostatistics with ML algorithms (e.g., RF, Boost, SVM, ANN) demonstrate promise in addressing spatial uncertainty, while RS data enhances covariate enrichment and near-real-time applications. Although advancements in variogram estimation and kriging techniques have optimized sampling strategies, and improved prediction accuracy, challenges persist in computational efficiency and uncertainty quantification. However, emerging trends such as, Bayesian methods, multiscale modeling, and sensor networks show the potential to mitigate these limitations. Future research should focus on improving model accuracy, interpretability, transparency, and explainability. Additionally, addressing data quality issues, scalability, and increasing sampling density to support robust spatial analysis and uncertainty quantification are major research gaps and future directions for soil mapping and information systems. This review is a valuable resource for researchers, practitioners, and policymakers engaged in data-driven land and soil management.