Optimized machine learning-based enhanced modeling of pile bearing capacity in layered soils using random and grid search techniques

Abstract

The bearing capacity of a pile is a critical factor in geotechnical design, necessitating extensive testing procedures that often increase both the time and cost of earthwork. Consequently, there is a growing demand for efficient and reliable methods to determine pile bearing capacity. This study aims to propose optimized machine learning based models through the application of Random Search (RS) and Grid Search (GS) optimization techniques for the prediction of pile-bearing capacity in layered soils. For this purpose, an extensive dataset is sourced from literature, and various machine learning algorithms including Random Forest (RF), Support Vector Machine (SVM), and XGBoost are investigated. Through a systematic modeling approach, multiple models are generated, and the performance of machine learning algorithms is refined using RS and GS cross validation (CV) using a customized code in Python. Optimized models are further assessed based on comprehensive evaluation criteria using key statistical performance indices. The results demonstrate that both RS and GS-tuned machine learning models achieve high accuracy, with R² values exceeding 0.9 and a low error index score across testing and training datasets. Notably, GS exhibits slightly superior statistical performance compared to RS. Furthermore, the tuned models with RS and GS showcase high performance on the validation dataset, with GS consistently outperforming RS. XGBoost emerges as the top performer among the machine learning models, followed by RF and SVM, highlighting the efficacy of tree-based algorithms in capturing the geotechnical variability inherent within pile bearing data. The proposed models offer valuable insights for predicting the preliminary evaluation of pile bearing capacity, facilitating swift and cost-effective geotechnical characterization within an acceptable error margin. This study introduces advancements in predictive modeling for geotechnical engineering, highlighting the transformative potential of optimization methodologies to enhance the machine learning models used for decision-making processes in civil engineering applications.

Introduction

Pile foundations are frequently utilized to support structures by penetrating deep into the ground and transferring loads onto denser and more stable soils or rocks. It is crucial to accurately ascertain their load-bearing capacity (Momeni et al. 2015). The maximum load that can be applied to a pile without causing failure in the soil beneath the pile or pile itself is known as the bearing capacity of the pile. Owing to this, pile bearing capacity is the most important parameter while designing piles (Drusa et al. 2016). For safe and economical construction of the project, an accurate estimation of bearing capacity is needed (Pham et al. 2020). Throughout years of research, various approaches have been adopted to determine the pile-bearing capacity such as empirical formulations and field testing. Empirical formulations offer insights based on historical data and observations, while static and dynamic field tests provide direct and real-time evaluations of pile performance under specific conditions. Empirical formulas often yield inaccurate estimates due to the complex interplay of external factors like soil mechanics and geotechnical conditions, as well as internal considerations such as pile dimensions (Liu et al. 2019). For the meticulous analysis, the presence of uncertainties in the impact of the hammer and variations in soil strength between the time of pile driving and the time of loading lead to challenges in accurately determining bearing capacity (Nguyen et al. 2024). The pile load test, which measures full-scale pile settlement under static load (Nguyen et al. 2023), is considered to offer highly accurate results (Shooshpasha et al. 2013). The reason for the accuracy is due to the method's similarity to the actual installation of driven piles. Despite its precision, the method's expense and time intensity make it less suited for small-scale projects, relegating its application primarily to larger ventures (Hoang et al. 2022). The measurement of soil properties using in-situ test methods has advanced; along with this progress, there has been a growing trend of using this data to predict the bearing capacity of piles (Nguyen et al. 2023). Out of in-situ tests, the Standard Penetration Test (SPT) is not only the most used test for determining the bearing capacity (Heidarie Golafzani et al. 2020), but also is a feasible alternative to pile load testing (Kozłowski and Niemczynski 2016). A number of past studies (Bazaraa and Kurkur 1986; Biarez and Foray 1977; Shioi and Fukui 2021) have used SPT results and proposed several empirical formulas to determine bearing capacity in sandy soil. Also, some studies have proposed several formulas to determine the bearing capacity for various types of soils including sand and clay. However, previous predictive models relied on basic regression techniques, resulting in the unreliability of predictions. Given advancements in data science and artificial intelligence, there is a pressing requirement to construct predictive models for pile bearing capacity utilizing machine learning algorithms.

Machine learning, a subset of artificial intelligence, draws insights from historical data by adaptively learning patterns, much like human brain functions (Mirrashid and Naderpour 2021; ur Rehman et al. 2024; Rehman et al. 2022). Recent research suggests that machine learning could serve as an alternative to traditional empirical methods and field tests for estimating the maximum capacity of structures (Alwanas et al. 2019; Mangalathu et al. 2022; Mohammed et al. 2020; Sarothi et al. 2022; Todorov and Billah 2022). Artificial neural networks (ANN) have been employed to estimate the friction capacity of driven piles in clay, supported by literature data for training the model (Nguyen et al. 2020). Armaghani et al., (Armaghani et al. 2020) formulated a hybrid intelligence model, combining ANN and particle swarm optimization (PSO), to predict pile settlement. Using data from various piles embedded in rock mass, the model was designed. Sensitivity analyses were then conducted to pinpoint key PSO parameters influencing prediction accuracy. Ghorbani et al. (2018) combined ANN and fuzzy inference system to predict the ultimate bearing capacity of axially loaded piles, using 108 load test data from the literature, including the Cone Penetration Test (CPT) measurements (Ghorbani et al. 2018). Momeni et al. (2015) introduced an ANN method to predict the axial bearing capacity of concrete piles. This model was trained using data from 36 field tests conducted with a pile driving analyzer, considering pile geometry and soil investigation results (Momeni et al. 2015). There have been few studies on utilizing pile load test data i.e., pile diameter, pile length, and along with SPT test results for establishing quick predictive models. Some studies have used single machine learning algorithms like Random Forest (RF), ANN and RF, Support Vector Regression (SVR), XGBoost, AdaBoost, and Decision Tree (DT) and compared their predicting performances. The others have focused on hybridized algorithms for example, the utilization of Genetic Algorithms (GA) to tune Deep Learning Neural Networks (Pham et al. 2020), the hybridization of Random Forest through Particle Swarm Optimization (Pham and Tran 2022), enhancement of Least Squares Support Vector Regression through ordinary differential flower pollination (ODFP) metaheuristic (Hoang et al. 2022), and the exploration of XGBoost parameter optimization through the Whale Optimization Algorithm, These studies not only laid the foundation for using pile load test data for bearing capacity determination but also demonstrate the effectiveness of hybridization and fine-tuning with algorithms like GA, PSO, and Whale Optimization Algorithms (WOA). These adaptations have been shown to significantly surpass the performance of standard models. However, relying solely on individual machine learning tools or single metaheuristic optimizers, such as GA, PSO, or WOA, risks missing the mark for global optimization. This risk may stem from falling into local minima traps or from the fact that these optimizers themselves require hyperparameter tuning, sometimes complicating simpler tasks. Moreover, it is not guaranteed that these models perform well in all scenarios; they might exhibit promising accuracy in one setting but falter in another, particularly as different terrains present unique challenges. Despite academic research limitations, the implementation of these models in the field remains unexplored, thus reducing the urgency for quick methods. Considering this discussion, there is a need for continued exploration beyond a single optimization technique. Proposing optimized algorithms that can be applied to predict the load-bearing capacity of driven piles is essential. Based on the findings of studies (Nguyen et al. 2023; Ge et al. 2023; Li and Li 2024; Yu 2024; Chala and Ray 2023), RF, SVR, and XGBoost are selected for the current study due to their proven ability to address complex and nonlinear problems, as well as their successful application in similar pile foundation tasks. These algorithms are particularly effective for high-dimensional and noisy datasets commonly encountered in pile load-bearing capacity prediction. RF and XGBoost, as tree-based ensemble methods, model relationships by iteratively partitioning the input space through decision trees, while SVR approximates the input–output relationship using a kernel-based optimization framework. However, these algorithms require further optimization to ensure reliable model performance, minimize bias in training on noisy data, and improve their generalizability and applicability. However, to the authors’ best knowledge, no previous study has applied Random Search (RS) and Grid Search (GS) optimization techniques to enhance the performance of these widely adopted algorithms for enhanced modeling of pile bearing capacity for quick and reliable determination without conducting extensive pile-load tests. This study aims to fill this gap by employing these advanced machine learning algorithms to model the pile bearing capacity and fine tuning these algorithms using RS and GS techniques, incorporating cross-validation (CV) to enhance prediction accuracy. Validation of these models is rigorously conducted through learning curves, cross-validation scores, and tests on independent validation data. To facilitate the preliminary evaluation of pile bearing capacity before conducting actual pile load tests, or for small projects where pile load test may be cost-prohibitive, a web application can serve as a practical and efficient tool. This approach refines the accuracy of the models and ensures their practical reliability, underscoring a commitment to clarity and reproducibility.

Methodology

Modeling data

Data description

The schematic illustration of the methodology for the current study is shown in Fig. 1. The dataset employed for this research was sourced from the work of Pham and Tran (Pham et al. 2020). This dataset comprised the results of field tests conducted on 472 pre-cast reinforced concrete piles in Ha Nam Province, Vietnam. These piles, characterized by their square cross-sections and closed tips, were installed using a hydraulic jack-in apparatus operating at a uniform velocity. The structure of the pile along with its geometric variables and the layers of soil surrounding it are shown in Fig. 2. It is noteworthy that the piles underwent testing no sooner than 7 days after installation. The testing procedure entailed progressively increasing vertical loads to 100%, 150%, and 200% of the design load, spanning periods of approximately 6 h, 12 h, and 24 h, respectively.

Fig. 1

The schematic illustration of methodology for current study

Fig. 2

The demonstration of pile load test in multi layered soil

The pile bearing capacity is determined based on one or both of the following criteria: (1) The pile's settlement at the current load matches or exceeds five times the settlement observed at the previous load level, or (2) The load-settlement curve is nearly linear. The pile bearing capacity corresponds to the load causing a settlement exceeding 10% of the pile's diameter. The pile bearing capacity corresponds to the load causing a settlement exceeding 10% of the pile's diameter. There are 10 input parameters i.e., pile diameter (D), depth of first soil layer (Z₁), depth of second soil layer (Z₂), depth of third soil layer (Z₃), elevation of pile top (Z_p), elevation of ground level (Z_g), elevation of pile extra segment of pile top (Z_t), elevation of pile tip (Z_m), average SPT count along pile shaft (N_sh), average SPT count at pile tip (N_t) and one output parameter the bearing capacity (P_u). Meanwhile, the average SPT value (N_t) is calculated according to Pham et al. (2020): considering values from 8 times the pile diameter (8D) above and 3 times the pile diameter (3D) below the pile tip, where D represents the pile's diameter. Additionally, N_sh is calculated by measuring the SPT blow count (N) at intervals of 300 mm (about 11.81 inches) along the embedded length of the pile. The histograms of inputs and outputs are shown in Fig. 3, and Table S1 provides insights into the data.

Fig. 3

The histogram plots of inputs and output (a) D; (b) Z; (c) Z_p and Zg; (d) Z_t and Z_m; (e) N_sh and N_t; (f) P_u

Input analysis

The specific selection of this data is based on its parameters, which directly influence pile bearing capacity. These parameters encompass pile geometry, thickness of soil layers, and soil characteristics obtained through SPT results. Furthermore, SPT blow counts are measured along the pile length to assess soil compressibility, a crucial factor in our analysis. This dataset is meticulously structured to encompass critical aspects of both pile and soil characteristics, aligning perfectly with the core objectives of our research and ensuring the accuracy of bearing capacity predictions.

The dataset comprises approximately 11 parameters, each with its statistical significance as depicted in Table 1. For parameters D and Z_g, low variability and a narrow range indicate that these features may not offer significant unique information for a model unless the target variable is closely associated with subtle changes in these features. Parameters Z₁, Z_p, Z_t, and N_t exhibit moderate variability, and their range suggests that these features may provide useful information to a model, particularly if they are correlated with the target variable. Similarly, parameters Z₂, Z_m, and N_sh display high variability, and a wide range suggests that these features may serve as important predictors, provided that their variability represents meaningful information rather than mere noise. Parameter Z₃ contains a large number of zeros because it indicates the depth of the third layer of soil.

Table 1 Statistical evaluation of data

In this study, feature selection was a crucial preliminary step before data analysis, recognizing that irrelevant or redundant features can adversely impact machine learning models. Such factors may lead to enhanced computational demands and reduced predictive accuracy. Therefore, for feature selection, recommendations from the study by Minh-Tu et al. (Cao et al. 2022) were taken into account. Irrelevant features, that do not contribute to understanding pile bearing capacity, exhibit statistical insignificance, or demonstrate dependent relationships, were eliminated. From the adopted dataset, features such as the elevation of the Z_t and Z_g were removed due to their weak correlation, low to moderate variability, and lack of domain relevance. Subsequently, the feature representing the elevation of the pile's lower head (Z_m) was eliminated despite its strong correlation and high variability because of its dependent nature and redundancy, which could adversely affect the model. Specifically, Z_m can be derived by summing the length of the pile embedded in each soil stratum (Z₁, Z₂, Z₃) with the Z_p. Additionally, further correlation coefficient analysis was conducted to reinforce the feature selection process, as illustrated in Fig. 4 and 7 input features were selected for modeling.

Fig. 4

The Correlation coefficient analysis of inputs

Preprocessing methods

Real-world data contains lots of noise, inconsistencies, and redundant information. The current dataset used in this study is extensive but exhibited significant skewness and outliers that could potentially bias the model. To address this, a detailed preprocessing procedure, including transformations and outlier handling, was implemented. These steps improved the symmetry of the data and reduced the presence of outliers, making the dataset more suitable for training machine learning models.

For the present data, a statistical approach named Gaussian approximation was implemented through a Python library Feature-engine (Galli 2021). It finds data points that lie further than 3 times standard deviations from the mean as outliers and removes outliers on both ends (tails) of observations.

$$Right\;tail = mean + 3\times std.$$

(1)

$$Left\;tail = mean - 3\times std.$$

(2)

Moreover, initially, the raw dataset underwent a preprocessing phase, where it is normalized using the MinMaxScaler function from scikit-learn. This preprocessing function re-scales each feature in the dataset to the range of 0 to 1, thereby ensuring that the features have the same scale. This is a crucial step, as it mitigates the potential for feature dominance and facilitates a more balanced learning process for the algorithms used.

Modeling and evaluation approaches

Three distinct machine learning models were used in this study i.e., RF, SVR, and Boost. RF aggregates the decision trees and forms a forest of them (Xue et al. 2021). Random forests generally outperform decision trees in terms of speed, accuracy, and computational cost (Gong et al. 2018). In this study, Python’s scikit-learn library is used for RF implementation, with hyperparameters tuned for optimal performance, as detailed in Table 2 (Zhu et al. 2022). Meanwhile, SVR, an extension of the Support Vector Machine for regression, is effective for handling nonlinear relationships and ensuring high prediction accuracy (Thomas et al. 2017). SVR finds an optimal hyperplane to bound predictions within a margin and is tuned using parameters like C, ε, and gamma, as shown in Table 3 (Astudillo, et al. 2020; Faris et al. 2018). XGBoost, a gradient-boosting algorithm, excels in predictive accuracy and efficiency, using regression trees to minimize prediction errors (Alibrahim and Ludwig 2021). It handles large datasets and complex relationships between input variables and the target variable (Marinov and Karapetyan 2019). Key hyperparameters such as "max_depth", "learning_rate", and "gamma" influence the model performance; their selected values for the current study are presented in Table 4 (Rehman et al. 2025; Prayogo et al. 2020; Nguyen et al. 2021).

Table 2 Hyperparameters and its ranges selected for optimization of RF

Table 3 Hyperparameters and its ranges selected for optimization of support vector regression

Table 4 Hyperparameters and its ranges selected for optimization of XGBoost

Moreover, K-fold cross validation (CV) is adopted in this study to develop models and avoid training bias caused by relying on a single data subset. In addition, hyperparameter tuning, which adjusts settings before training, is employed to enhance the modeling process (Probst et al. 2019; Sipper 2022). This involves comparing performance against previous configurations to identify the optimal setup (Elgeldawi et al. 2021). Two search-based algorithms are employed in this study to optimize machine learning models. The Random Search (RS) approach involves randomly selecting hyperparameter combinations from a search space to efficiently explore and identify optimal settings for a model (Ren, et al. 2022). RandomizedSearchCV (RSCV), the Python implementation of RS, samples fixed combinations of hyperparameters using the "param_distributions" and "n_iter" arguments and incorporates cross-validation for robust performance estimation. It stores results in "cv_results" and identifies the best-performing model through "best_estimator" and "best_params" (Sarajcev and Meglic 2022), as shown in Fig. 5. In contrast, Grid Search (GS) systematically explores all possible hyperparameter combinations using the "param_grid" argument in GridSearchCV (GSCV), with a focus on model performance and generalization (Alibrahim and Ludwig 2021; Marinov and Karapetyan 2019). GSCV begins with a wide search space and narrows it down based on literature insights to ensure optimal tuning and generalization (Belete and Huchaiah 2022), as illustrated in Fig. 6.

Fig. 5

The flow chart of RSCV

Fig. 6

The flow chart for GSCV

Meanwhile, this study utilized various model metrics to evaluate the performance of ML regression models by assessing their prediction accuracy on training and test datasets. For this purpose, the coefficient of determination (R²), Root Mean Squared Error (RMSE), Mean Absolute Error (MAE) and Mean Squared Error (MSE) are considered to gauge the statistical performance of the models as follows.

$$RMSE= \sqrt{\frac{{\sum }_{i=1}^{n}{({E}_{i}-{M}_{i})}^{2}}{n}}$$

(3)

$$\text{MAE}= \frac{{\sum }_{i=1}^{n}|{E}_{i}-{M}_{i}|}{n}$$

(4)

$$\text{MSE}= \frac{1}{n}{\sum }_{i=1}^{n}{({E}_{i}-{M}_{i})}^{2}$$

(5)

$${\text{R}}^{2}=1-\frac{{\sum }_{i=1}^{n}{({E}_{i}-{M}_{i})}^{2}}{{\sum }_{i=1}^{n}{({E}_{i}-{\overline{E} }_{i})}^{2}}$$

(6)

where ${E}_{i}$ is predicted value, ${M}_{i}$ is real value, ${\overline{E} }_{i}$ is the mean of predicted values and $n$ is total number of data points. In addition, learning curves are also used which assist in understanding how a model's performance changes with the size of the training data. They play a crucial role in various machine learning scenarios, helping with data acquisition, determining when to stop model training, and selecting the appropriate model (Mohr and Rijn 2201). Models are also evaluated using cross-validation scores for both training and testing, as well as assessing standard training and testing scores.

Results and discussion

Preprocessed data

The database after removing outliers consists of 447 data points and 8 parameters. The statistical summary of the cleaned data points is presented in Table 5. For this data, the D ranges from 300 to 400 mm, with a mean of 361.74 mm and a standard deviation of 48.65 mm. The Z₁, Z₂, and Z₃ exhibit mean values of 3.801 m, 6.815 m, and 0.346 m, respectively, with varying standard deviations, highlighting the variability in layer thicknesses. The parameter Z_p has a mean of 2.815 m, ranging from 1.95 to 3.4 m. The N_sh and the N_t show consistent averages of 10.965 and 7.177, respectively. The P_u varies widely, with values ranging from 407.2 to 1551, and a mean of 1005.519, reflecting significant variability across the dataset. Generally, outlier removal has led to slight changes in the mean, while the median and mode have remained mostly consistent. This is because outliers have a greater impact on the mean but less on the median and mode, which are more robust measures of central tendency.

Table 5 Statistical summary of cleaned data

For Z₁, the mean has decreased from 3.8238 to 3.8018, skewness has decreased from 0.5477 to 0.5034, resulting in a change in range of –0.97. For Z₂, the mean has increased from 6.579 to 6.8158, skewness has increased from −1.0327 to −0.3564, signifying a shift towards symmetry, with a change in range of −3.56. For Z₃, the mean has slightly increased from 0.3315 to 0.3462, skewness has reduced from 0.8428 to 0.7459, indicating a move towards a more symmetric distribution, with a change in the range of −0.47. For Z_p, the mean has slightly changed from 2.8038 to 2.8153, skewness has also slightly moved from −0.3423 to −0.3125, implying a small move towards symmetry, with a change in range of −0.73. For N_sh, the mean has increased from 10.7426 to 10.9660, skewness has moved from −0.2149 to −0.0408, indicating a more symmetric distribution post-outlier removal, with a change in range of −4.73. For N_t, the mean has increased from 7.0557 to 7.1776, and skewness has substantially changed from −1.9806 to almost 0.0071, implying a significant move towards symmetry, with a change in range of −2.33. Finally, for P_u, the mean has increased from 984.2014 to 1005.5187, skewness has moved from −0.0910 to −0.1725, indicating a slight move towards negative skewness, with no change in range.

In summary, the range has consistently decreased, indicating the successful removal of extreme values and resulting in a more compact dataset. The mode remains consistent across all inputs. The increase in the mean and the substantial shift in skewness from a negative value to one closer to zero suggest that the outliers were predominantly on the lower side, making the original distribution left-skewed. The slight reduction in both mean and skewness further indicates that the outliers had been slightly pulling the mean upward, making the distribution more positively skewed.

For each input parameter, minimum and maximum values are changed. After data normalization, the summary of rescaled data is shown in Table 6. The range of all input parameters was from 0–1 for normalized data with a standard deviation of 0.352 to 0.486. For the present study, the preprocessed data is divided into training, testing, and validation datasets. Specifically, 90% of data is used during hyperparameter tunning process, where cross-validation divides the data into multiple folds. In each fold, a portion of the data is used for training the model, while the remaining folds are used for testing. This iterative process ensures that this 90% of the data contributes to both training and testing the models. which is further divided through CV. Additionally, the remaining 10% of data is held out throughout the model development process and serves only as the external validation dataset to evaluate final model performance.

Table 6 Statistical summary of normalized data

Optimized hyperparameters

All the models undergo a search process to find the model with the highest mean scoring metric of five folds (i.e., K = 5). The best set of hyperparameters, which are shown in Tables 7, 8, and 9, along with corresponding mean train and test metrics, are presented in Fig. 7. There's a noticeable, though not substantial, difference between the mean training and test scores. This hints at slight overfitting, suggesting that the model might be capturing some noise present in the training data. The Fig. 7b, d, f and h contains the standard deviation between folds. A lower standard deviation, which is basically between model metrics during each fold split, indicates that the model's performance is consistent across different subsets of the dataset. The standard deviations of the metrics across the five folds (Fig. 7b, d, f, g, and h), particularly for the test scores, hint at the model's stability. R²_train (Fig. 7a) is higher than R²_test (Fig. 7c), which indicates there is some indication of the model capturing noise or exhibiting mild overfitting. As the model undergoes a cross-validation process, a scatter plot during each fold is shown in Fig. 8.

Table 7 Tuned/optimal hyperparameters for RF

Table 8 Tuned/optimal hyperparameters for SVR

Table 9 Tunned/optimal hyperparameters for XGBoost

Fig. 7

Mean training and test scores with standard deviations across five folds a R²_train; b R²_train(std); c RMSE_train; d RMSE_train(std); e R²_test; f R.²_test(std); g RMSE_test; h RMSE_test(std)

Fig. 8

The scatter plots for tuned models a RS-RF; b GS-RF; c RS-SVR; d GS-SVR; e RS-XGBoost; f GS-XGBoost

Performance evaluation

Evaluation of tuned models through learning curves

To gain insights into models that are presumably well-optimized for the dataset, learning curves are plotted. The variance (standard deviation) around both the training and validation scores provides an idea of how much variation exists in performance across different cross-validation folds.

For RF, the learning curve (Fig. 9) using RS with R² (Fig. 9a) for the training score curve, starting from 0.921 and ending around 0.949, indicates the model has a high goodness-of-fit on the training data. The slight dip and subsequent plateau suggest that with more data, the model's performance on the training set stabilizes. Regarding the validation curve, starting from 0.872 and ending around 0.929 shows that as the model is trained with more data, its ability to generalize improves, but there's still a gap between training and validation, suggesting slight overfitting. The same behavior can be observed in the curve with RMSE on the y-axis (Fig. 9b). For RF, learning curves depicting model fit with optimal hyperparameters obtained through GS (Fig. 10) with R² on the y-axis are displayed in Fig. 10a. These curves exhibit similar behavior to those derived from RS, indicating that RS performed adequately in identifying a set of hyperparameters comparable in performance to those discovered by GS. The slight upward trend observed at the end of the validation curve can be interpreted as a positive sign – suggesting that the model may benefit from additional training data, potentially leading to further improvement in its generalization capability. Employing both R² and RMSE for learning curves provides a comprehensive understanding of the model's behavior. While R² offers insights into the model's data-fitting ability, RMSE quantifies the magnitude of the error. The insights gleaned from these learning curves are valuable for comprehending the model's behavior and identifying potential areas for enhancement. The learning curves for SVR (Fig. 11) with optimal hyperparameters obtained from both RS (Fig. 11a and b) and GS (Fig. 11b and d) exhibit concurrent dips in the training and validation curves around 192 samples, suggesting that the model encountered challenges with that specific subset of the training data. This occurrence could be attributed to particular characteristics or noise within those samples. Furthermore, the downward trend observed in the validation curve towards the end suggests a potential slight overfitting of the model, particularly evident as the gap between the training and validation scores persists.

Fig. 9

Learning Curves for Random Forest tuned through Random Search a R²; b RMSE

Fig. 10

Learning Curves for Random Forest tunned through GS a R²; b RMSE

Fig. 11

Learning Curves for SVR tunned through RS a R²; b RMSE, and GS c R²; d RMSE

For XGBoost, the learning curve with optimal hyperparameter of RS (Fig. 12a and b), the training curve starts from 0.925 and steadily increases to stabilize at 0.954 and the validation curve starts from 0.876 and steadily increases to stabilize at 0.93. The consistent and small gap between the training and validation curves throughout suggests that the model is neither overfitting nor underfitting significantly. This indicates that the model generalizes well to unseen data. For the learning curve with optimal hyperparameters obtained from GS, as shown in Fig. 13a and b, the training curve initiates from a high point of 0.925, experiences a slight dip around 100 samples, though not significantly, and eventually stabilizes around 0.970. In contrast, the validation curve begins at 0.825, undergoes a minor decline around 225 samples, but ultimately stabilizes around 0.93. Although the gap between the training and validation curves diminishes over time, it persists around 0.04 by the end. While this difference is small, it is noticeable, indicating a mild overfitting as the model demonstrates a slightly superior fit on the training data compared to the validation data.

Fig. 12

Learning Curves for XGBoost tuned through RS a R²; b RMSE

Fig. 13

Learning Curves for XGBoost tuned through GS a R²; b RMSE

The model trained with RS (Fig. 12a and b) hyperparameters exhibits a narrower gap between the training and validation sets, suggesting it is slightly more robust than the one trained with GS (Fig. 13a and b) hyperparameters. Despite achieving a higher R² on the training set, the GS model displays a more prominent gap, indicating mild overfitting. The slight declines observed in the curves (at 100 samples for the training curve and 224 samples for the validation curve) in the GS (Fig. 13a and b) model might be attributed to specific data subsets, noise, or transient overfitting during training. Nevertheless, the R² scores are relatively high for both models, signifying good predictive capabilities.

Evaluation of tuned models on validation data

To evaluate and compare the performance of the tuned models, all models were validated on 10% of the data, which was set aside at the beginning of the process and not used in in model development. The results were then compared with those of untuned or base models (Fig. 14). Initially, XGBoost outperformed the other models with an R²_validation of 0.923 (Fig. 14a) and relatively lower error metrics (MAE_validation = 69.056, RMSE_validation = 98.091, MSE_validation=9621.96), showcasing its superior ability to capture the complexities in the data. RF and SVR also performed reasonably well but had higher errors, indicating room for improvement (Fig. 14). Hyperparameter tunned models significantly performed better than base models. For RS-RF, it led to substantial improvements, with GS-RF (Fig. 14). This highlighted the value of an exhaustive search in fine-tuning model parameters, particularly for complex tasks like predicting load-bearing capacity. Similarly, RS-SVR and GS-SVR showed notable gains from tuning, though it still lagged tunned RF and XGBoost, even after optimization.

Fig. 14

The comparison of base and tuned models on validation data a R²_validation; b MAE_validation; c RMSE_validation; d MSE_validation

XGBoost, however, benefited the most from tuning. The GSCV-tuned XGBoost model (GS-XGBoost) achieved the highest accuracy across all models, with an R² of 0.952 (Fig. 14a) and further reduced error metrics (MAE_validation = 60.88, RMSE_validation = 76.45, MSE_validation = 5844.96). This demonstrates that XGBoost, when fine-tuned, offers exceptional predictive power and reliability, making it the most robust model for this application. Overall, while all models improved with tuning, XGBoost stood out as the most effective tool, particularly when paired with GSCV for hyperparameter optimization. However, the differences in performance between the hyperparameter tuning methods are relatively small due to the higher number of iterations during RS. Overall, the results show the GS technique for hyperparameter tuning tends to yield slightly better results compared to RS, as evidenced by the higher R² and lower error metrics (i.e., MAE, RMSE, and MSE). Overall, this analysis can guide the selection of models and hyperparameter tuning methods for similar tasks, emphasizing the importance of model selection and hyperparameter optimization in achieving optimal performance in machine learning tasks. The graphs in Fig. 15 show the predicted values versus the actual values. Each algorithm, when tuned using RSCV, was compared against its counterpart tuned with GSCV.

Fig. 15

The scatter plots for tuned models on validation data a RS-RF, GS-RF; b RS-SVR, GS-SVR; c RS-XGBoost, GS-XGBoost

For the simultaneous visualization of multiple statistical measures, the predictions of each tuned model are evaluated on a Taylor diagram (Fig. 16). This analysis reveals that models tuned using the Grid Search technique, particularly the XGBoost and RF algorithms, demonstrate superior performance. They exhibit R² of 0.977 and 0.9769, and standard deviations of 338.552 and 336.988, respectively. This superiority is also evident in their lower RMSE values, indicating a strong linear relationship with actual values and minimal prediction errors. Additionally, the XGBoost model tuned through a RS technique shows promising results. The Grid Search—XGBoost (GS-XGBoost) model, in particular, stands out as the top performer, exhibiting the highest correlation and the lowest RMSE among all models evaluated. This suggests the remarkable capability of this model in capturing the underlying relationships in the data, making it a reliable choice for similar datasets or prediction tasks.

Fig. 16

Taylor diagram analysis of different models

The predictions by turned models on 10% validation data are presented in Table S2. All these tuned models are also evaluated by finding the distribution of errors, which is shown in Fig. 17 by revealing their performance through the spread of their error distributions. The overall graphs (Fig. 17) suggest that the models have their errors symmetrically distributed around zero, which indicates no significant bias in overestimating or underestimating. The tuned RF (Fig. 17a) and XGBoost (Fig. 17c) models show a tight clustering of errors when refined through GS, implying a potential refinement in prediction accuracy. The tuned SVR models (Fig. 17c) exhibit a broader spread in errors, indicating a range of predicted deviations. The visual overlap in error ranges between Random and Grid Search suggests a nuanced difference in their tuning effectiveness.

Fig. 17

The error rate percentage for the tuned models based on the normal distribution a RS-RF, GS-RF; b RS-SVR, GS-SVR; c RS-XGBoost, GS-XGBoost

Comparison with existing models

The prediction performance of the proposed tuned ML model, GS-XGBoost, was compared against those of other proposed models in the literature (Table 10). The performance indices of existing models were obtained based on data points used in the current study for model development ensuring comparability. It can be observed that the GA-tuned models significantly decreased the error of DLNN and RF. Also, metaheuristics approaches like PSO, and WAO also showed promising performance. Meanwhile, GS-XGBoost showed superior performance amongst all previously proposed models, with R² value of 0.952 outperforming alternatives such as GA-DLNN (R² = 0.923) and RF-PSO (R² = 0.924). The superior performance of GS-XGBoost can be explained based on different factors; firstly, the grid search optimization systematically explored the hyperparameter space, ensuring the best configuration for capturing the complex, nonlinear relationships inherent in the current dataset. Secondly, the XGBoost algorithm's inherent strengths, including gradient boosting with regularization, made it more effective at handling the high-dimensional, noisy, and variable geotechnical data compared to models like DLNN, which are more prone to overfitting, and RF, which may struggle with bias due to simpler ensembling (Belete and Huchaiah 2022). Additionally, GS offers a balanced exploration of the hyperparameter space, unlike stochastic methods such as PSO and WAO, which, while promising, may require larger computational budgets and may not always converge to the global optimum. This systematic approach, coupled with XGBoost's robustness and computational efficiency, allowed GS-XGBoost to achieve superior performance in both accuracy and generalization. Although differences in computational budgets and optimization strategies exist among these methods, the results highlight the effectiveness of the GS-XGBoost model proposed in this study, which leverages the combination of systematic tuning and the algorithmic strengths of XGBoost to deliver outstanding performance in static pile load test predictions.

Table 10 Existing and current model comparison 

Field Implications

Determining the pile loading bearing capacity is a significant geotechnical challenge as it directly impacts the stability and safety of constructions. Although pile loading tests are accurate and reliable across various design stages, their high costs and time requirements often deter designers from conducting field tests. Despite advancements in testing methodologies and analytical tools, this process remains challenging and often involves considerable uncertainty. The current study proposed an enhanced machine learning-based model for predicting pile bearing capacity by comprehensively applying optimization techniques and rigorously evaluating the models. To extend the practical implications of this research, the author developed a web-based application using the proposed model to assist engineers and designers in leveraging the best-performing model in real-world scenarios. The application accepts input parameters such as soil layer thickness, pile diameter, and SPT values along the shaft and tip, and computes the pile bearing capacity. This pile bearing capacity corresponds to a settlement of 10% of the pile diameter. The web application is accessible at: https://pilebearingcapacity.streamlit.app/. By using this tool, users can make informed decisions during the preliminary evaluation of piles based on readily known parameters from the field. Additionally, the application allows users to adjust various parameters and analyze the desirable soil properties for achieving specific pile bearing capacity thresholds. This enables engineers to explore different design scenarios, aiding in effective decision-making during the early stages of earthwork projects involving pile construction.

Moreover, it is important to note that the current model is trained on a specific range of data, as described earlier. Therefore, its application is limited to these ranges which are also indicated in the web-based application. Additionally, such models are suitable for the preliminary evaluation (Rehman et al. 2025; Prayogo et al. 2020); therefore, for detailed assessments, pile bearing capacity tests are recommended considering the nature of the project. To further enhance the applicability of the proposed models, future studies could expand the dataset used in this research. Although an extensive dataset has already been utilized, future enhancements could focus on broadening the input parameter ranges, incorporating additional parameters or material property information, and including more tests and bearing capacity data.

Conclusion

This study delved into the application of RS and GS optimized machine learning techniques for predictive modeling of pile bearing capacity. Leveraging an extensive dataset curated from literature, this study explored the efficacy of various machine learning algorithms including RF, SVM, and XGBoost. Through systematic experimentation, multiple models were generated and the performance of machine learning algorithms was refined by integrating RSCV and GSCV techniques. The best models for predicting pile-bearing capacity are identified using comprehensive evaluation criteria. The following are the main findings of this study.

• Both RS and GS-tuned machine learning models exhibit remarkable accuracy, boasting R² values exceeding 0.9 and RMSE values less than 100 kN across both testing and training datasets. Notably, GS demonstrates a marginally superior statistical performance compared to RS. Furthermore, our tuned models with RS and GS showcase commendable performance on the validation dataset, with GS consistently outperforming RS.

• In terms of machine learning model performance, XGBoost emerges as the top performer across training, testing, and validation datasets, followed closely by Random Forest and Support Vector Machine. This underscores the effectiveness of tree-based algorithms in effectively capturing the intricate geotechnical variability inherent within pile bearing data.

• The models proposed in this study serve as invaluable tools for predicting the preliminary evaluation of pile bearing capacity, enabling swift and cost-effective geotechnical characterization within an acceptable error margin. This research contributes to the advancement of predictive modeling in geotechnical engineering and underscores the significance of leveraging machine learning methodologies for informed decision-making in civil engineering applications.

This study also presents a web-based application for evaluating pile bearing capacity using the best-performing model. However, this model is applicable to specific input ranges and is recommended for preliminary evaluations. For detailed design, pile bearing capacity tests are advised. Future studies could focus on expanding the applicability of the proposed model by improving the dataset and utilizing advanced optimization techniques and error analysis across the subset for model training and testing.