Introduction

The design of pile foundation has drawn more attention than other type of foundation structures. The use of axial loaded pile is more frequent and is designed using equations of static equilibrium and other dynamic equations [38]. However, the lateral loaded piles are used in more difficult conditions, particularly in tall and offshore structures. The design of laterally loaded piles is more difficult and requires solution of nonlinear differential equations. The elastic analysis adopting Winkler soil model [38] is not suitable for the nonlinear soil behavior. Matlock and Reese [33] adopted elastic analysis using nonlinear lateral load capacity—deflection (p-y) curves. Portugal and Seco e Pinto [37] used nonlinear p-y curves and finite element method for prediction of the behavior of laterally loaded piles. The above two methods are more accurate and widely used. But, spatial variability of soil is inevitable. Thus, developing a sufficiently accurate site model for FEM analysis requires extensive site characterization effort and desired constitutive modeling of clayey soil is also very difficult, even with considerable laboratory testing. So methods based on field data [8, 27, 34] have become very much popular for the above study, particularly for the preliminary estimate of pile load capacity. These methods are based on pile load test case histories and involve statistically derived empirical equations for estimation of expected lateral load capacity.

Artificial intelligence (AI) techniques such as artificial neural networks (ANNs) and support vector machine (SVM) are considered as alternate statistical methods and are found to be more efficient compared to statistical methods [11, 13]. ANN method has been found to be efficient in predicting the pile load capacity in both cohesion- less soil and clayey soil compared to traditional empirical methods [2, 9, 25, 26, 31, 45]. The performance of SVM model was found to be better than that of ANN model for prediction of frictional resistance of pile in clay [41]. Similar studies have also been made for prediction of lateral load capacity of piles in clay using ANN [13]. Based on various statistical performance criteria, Das and Basudhar [13] observed that ANN model is better compared to Brom’s and Hansen’s method. Using the same data set, Pal and Deswal [36] developed Gaussian process regression (GPR) and SVM models. They observed that GPR model is better compared to SVM model. However, they have compared the GPR model with the ANN model of Das and Basudhar [13] only in terms of correlation coefficient (R) and root mean square error (RMSE). The R is a biased estimate [10] and it is difficult to assess the prediction of the model in terms of under prediction or over prediction on the basis of R value only. The RMSE explains the overall error of the dataset instead of the maximum deviation in the prediction of individual case.

The most important problem associated with efficient implementation of ANN is generalization for some complex problems. The developed model needs to be equally efficient for new data during testing or validation, which is called as generalization. There are different methods for generalization like early stopping and cross validation [6, 13]. The magnitude of weight is one of the reasons for poor generalization [5]. The methods like Bayesian regularization neural network (BRNN) [14] has been used to consider the magnitude of weights as the part of the error function. Another reason for the poor generalization is due to the optimization of error function of ANN. The error function associated with weights and sigmoid function is a highly non-linear optimization problem with many local minima [14]. As the characteristic of traditional nonlinear programming based optimization method are initial point dependent, the use of global optimization algorithms like genetic algorithm and simulated annealing are being widely used in training ANN model [3, 24, 35]. The training of the feed-forward neural network using differential evolution optimization is known as differential evolution neural network (DENN) [11, 28]. Das et al. [12] observed that performance of DENN is better than BRNN and traditionally used Levenberg–Marquardt neural network (LMNN) for the slope stability analysis. The ANN is termed as a ‘black box’ system unable to explain the input output relationships and in SVM error parameter ‘C’ and sensitive function ‘e’ are to be found out by trial and error. However, now it is possible to write down a model equation based on the trained ANN model [13, 14, 24] and SVM model [11, 16], still the developed model particularly SVM model is not comprehensive. The modified artificial intelligence techniques in the class of ‘grey box’ and ‘white box’ are now a day being popular [23]. The genetic programming (GP) is defined as next generation AI technique and also called as ‘grey box’ model [23] in which the mathematical structure of the model can be derived, allowing further information of the system behaviour. The GP and its variants have been applied to few difficult geotechnical engineering problems [4, 15, 19, 20, 29, 40, 46] with success. A modified statistical technique called multivariate adaptive regression spline (MARS) is popularized by Friedman [18] for solving regression-type problems. MARS is also called as ‘white box’ system of predictive model, as it is based on physical laws and underlying physical relationships of the system can be explained. The MARS technique is very popular in the area of data mining because it does not assume or impose any particular type or class of relationship (e.g., linear, logistic, etc.) between the predictor variables and the dependent (outcome) variables of interest. This makes MARS particularly suitable for problems with more number of variables. It has an increasing number of applications in many areas of economy, science and technology. However, its use in geotechnical engineering is very limited [42, 43]. It needs to compare efficacy of the present GP and MARS models vis-à-vis ANN, and other empirical models in terms of different statistical performance criteria.

In the present study prediction models for lateral load capacity of piles in clay under un-drained condition have been developed using GP, MARS and ANN (BRNN, DENN). Different statistical criteria like correlation coefficient (R), Nash-Sutcliff coefficient of efficiency (E) [14], absolute average error (AAE), maximum absolute error (MAE) and root mean square error (RMSE) are used to compare the GP and MARS models with ANN (DENN, BRNN) models and existing empirical models (Broms and Hansen’s). As R value alone is not a good indicator of prediction accuracy of a model, hence, a objective function (OBJ) as per Gandomi et al. [21] is used to compare the statistical performance of various models considered in the present study. The OBJ takes into account the change of R, RMSE and AAE together. A ranking system [1] using rank index (RI) has also been followed to compare the different models basing on four criteria: (i) the best fitness calculations (R and E) for predicted lateral load capacity (Q p ) and measured lateral load capacity (Q m ), (ii) arithmetic calculations (mean, µ and standard deviation, σ) of the ratio, Q p /Q m (iii) 50 and 90 % cumulative probabilities (P 50 and P 90) of the ratio, Q p /Q m. and (iv) the probability of pile load capacity within 20 % accuracy level in percentage using histogram and lognormal probability distribution of Q p /Q m .

Methodology

The ANN has been extensively used in geotechnical engineering. In the present study, the ANN models are trained with differential evolution and Bayesian regularization method and are defined as DENN and BRNN respectively. But the details of GP and MARS are discussed briefly as they are not very common to geotechnical engineering professionals.

Artificial Neural Networks (ANN)

In the present study, the ANN models are trained with differential evolution and Bayesian regularization method and are defined as DENN and BRNN respectively. The use of DENN and BRNN are limited in geotechnical engineering [12, 13, 14, 24]. A brief description about the Bayesian regularization and differential evolution neural network is presented here for completeness.

Bayesian Regularization Neural Network (BRNN)

In case of back propagation neural network (BPNN) the error function considered for minimization is the mean square error (MSE). This may lead to over-fitting due to unbounded values of the weights. The other method called as regularization, in which the performance function is changed by adding a term that consist of mean square error of weights and biases as given below

$$ MSEREG\, = \,\gamma \,MSE\, + \,(1 - \gamma )\,\,MSW $$
(1)

where MSE is the mean square error of the network, γ is the regularization parameter and

$$ MSW\, = \,\frac{1}{n}\sum\limits_{j\, = \,1}^{n} {w_{j}^{2} } $$
(2)

This performance function will cause the network to have smaller weights and biases thereby forcing networks less likely to be over-fit. The optimal regularization parameter γ is determined through Bayesian framework [17] as the low value of γ will not adequately fit the training data and high value of it may result in over-fitting. The number of network parameters (weights and biases) are being effectively used by the network can be found out by the above algorithm. The effective number of parameters remains the same irrespective of the total number of parameters in the network.

Differential Evolution Neural Network (DENN)

The differential evolution (DE) optimization is a population based heuristic global optimization method. Unlike other evolutionary optimization, in DE the vectors in current populations are randomly sampled and combined to create vectors for next generation with real valued crossover factor and mutation factor. The detail of DENN is available in Ilonen et al. [28].

Genetic Programming

Genetic Programming is a pattern recognition technique where the model is developed on the basis of adaptive learning over a number of cases of provided data, developed by Koza [30]. It mimics biological evolution of living organisms and makes use of principle of genetic algorithm (GA). In traditional regression analysis the user has to specify the structure of the model whereas in GP both structure and the parameters of the mathematical model are evolved automatically. It provides a solution in the form of tree structure or in the form of compact equation using the given dataset. A brief description about GP is presented for the completeness, but the details can be found in Koza [30].

GP model is composed of nodes, which resembles to a tree structure and thus, it is also known as GP tree. Nodes are the elements either from a functional set or terminal set. A functional set may include arithmetic operators (+, × , ÷, or −), mathematical functions (sin(.), cos(.), tanh(.) or ln(.)), Boolean operators (AND, OR, NOT etc.), logical expressions (IF, or THEN) or any other suitable functions defined by the user. The terminal set include variables (like x1, x2, x3, etc.) or constants (like 3, 5, 6, 9 etc.) or both. The functions and terminals are randomly chosen to form a GP tree with a root node and the branches extending from each function nodes to end in terminal nodes as shown in Fig. 1.

Fig. 1
figure 1

A GP tree expressing function: (5X1 + X2)2

Initially a set of GP trees, as per user defined population size, is randomly generated using various functions and terminals assigned by the user. The fitness criterion is calculated by the objective function and it determines the quality of the each individual in the population competing with the rest. At each generation a new population is created by selecting individuals as per the merit of their fitness from the initial population and then, implementing various evolutionary mechanisms like reproduction, crossover and mutation to the functions and terminals of the selected GP trees. The new population then replaces the existing population. This process is iterated until the termination criterion, which can be either a threshold fitness value or maximum number of generations, is satisfied. The best GP model, based on its fitness value that appeared in any generation, is selected as the result of genetic programming. A brief description on various evolutionary mechanisms in GP are presented below.

Initial Population

In the first step of genetic programming a number of GP trees are generated by randomly selecting user defined functions and terminals. These GP trees form initial population.

Reproduction

In the second stage of the GP, a proportion of the initial population is selected and copied to the next generation and this procedure is called reproduction. Roulette wheel selection, tournament selection, ranking selection etc. are the methods generally followed for the selection procedure.

Crossover

In crossover operation, two trees are selected randomly from the population in the mating pool. One node from each tree is selected randomly, the sub-trees under the selected nodes are swapped and two offsprings are generated as shown in Fig. 2.

Fig. 2
figure 2

A typical crossover mechanism in GP

Mutation

A GP tree is first selected randomly from the population in the mating pool and any node of the tree is replaced by any other node from the same function or terminal set as shown in Fig. 3. A function node can replace only a function node and the same principle is applicable for the terminal nodes.

Fig. 3
figure 3

A typical mutation mechanism in GP

The general form of proposed GP model can be presented as:

$$ Q_{p} = \sum\limits_{i = 1}^{n} {F\left[ {X,f\left( X \right),b{}_{i}} \right]} + b_{0} $$
(3)

where, F = the function created by the GP referred herein as pile load function, X = vector of input variables = {D, L, e, S u }, where D = diameter of pile, L = depth of pile embedment, e = eccentricity of load, S u  = un-drained shear strength of soil, b i is constant, f is the function defined by the user and n is the number of terms of target expression and b 0  = bias. The GP as per Searson et al. [44] is used and the present model is developed and implemented using Matlab [32].

Multivariate Adaptive Regression Spline (MARS)

MARS is an adaptive procedure because the selection of basis functions is data-based and specific to the problem at hand. This algorithm is a nonparametric regression procedure that makes no specific assumption about the underlying functional relationship between the dependent and independent variables. It is very useful for high dimensional problems. For this model an algorithm was proposed by Friedman [18] as a flexible approach to high dimensional nonparametric regression, based on a modified recursive partitioning methodology. MARS uses expansions in piecewise linear basis functions of the form Eq. (4)

$$ c^{ + } \left( {x, \tau } \right) = \left[ { + \left( {x - \tau } \right)} \right]_{ + }, \quad c^{\_} \left( {x, \tau } \right) = \left[ { - \left( {x - \tau } \right)} \right]_{ + } $$
(4)

where [q] = max{0, q} and τ is an univariate knot. Each function is piecewise linear, with a knot at the value τ, and it is called a reflected pair.

The points in Fig. 4 illustrate the data (xi, yi) (i = 1, 2, … N), composed by a p-dimensional input specification of the variable x and the corresponding 1-dimensional responses, which specify the variable y.

Fig. 4
figure 4

A basic element in the regression with MARS

Let us consider the following general model Eq. (5) on the relation between input and response:

$$ Y = f(X) + \varepsilon $$
(5)

where, Y is a response variable, X = (X1,X2, …, X Tn) is a vector of predictors and ε is an additive stochastic component, which is assumed to have zero mean and finite variance.

The goal is to construct reflected pairs for each input xj (j = 1, 2, …, p) with p-dimensional knots τi = (τi,1, τi,2, …, τi,p)T. Actually, we could even choose the knots τi,j more far away from the input values xi,j, if any such a position promises a better data fitting.

After these preparations, our set of basis functions is Eq. (6):

$$ \delta : = \{ (X_{j} - \tau )_{ + } ,(\tau - X_{j} )_{ + } |\tau \in \{ x_{1,j} ,x_{2,j} , \ldots ,\,x_{N,j} \} , \quad j \in \{ 1,\,2, \ldots ,p\} \} $$
(6)

If all of the input values are distinct, there are 2Np basis functions altogether. Thus, we can represent f (X) by a linear combination, which is successively built up by the set δ and with the intercept θ0, such that Eq. (6) takes the form

$$ Y = \theta_{0} + \sum\limits_{m = 1}^{M} {\theta_{m} \psi_{m} (X) + \varepsilon .} $$
(7)

Database and Preprocessing

In the present study the experimental database of Rao and Suresh Kumar [39] has been considered. Das and Basudhar [13] have developed ANN model and Pal and Deswal [36] have developed GPR and SVM models using the above database. The data consist of D, L, e, S u as the inputs and Q m as output. Out of the mentioned 38 data, 29 data are selected for training and remaining 09 data are used for testing the developed model as per Das and Basudhar [13]. The data are normalized in the range [0, 1] and [−1, 1] for MARS and ANN (DENN, BRNN) models respectively to avoid the dimensional effect of input parameters. In the GP modeling normalization or scaling of the data is not required.

Results and Discussion

In the present study each individual in the population consists of more than one gene and each gene is a traditional GP tree. Here, function set used include: +, × , ÷, −, sin(.), cos(.), tanh(.) and exp(.). As discussed earlier in GP procedure first a number of potential models are evolved at random. Each model is trained and tested using the training and testing cases respectively. The fitness of each model is determined by minimizing RMSE between the predicted (Q p ) and actual (Q m ) value of the output variable as the objective function,

$$ RMSE = f = \sqrt {\frac{{\sum\nolimits_{i = 1}^{n} {\left( {Q_{m} - Q_{p} } \right)^{2} } }}{n}} $$
(8)

where n = number of cases in the fitness group. If the errors calculated by using Eq. (8) for all the models in the existing population do not satisfy the termination criteria, the generation of new population continues till ‘best’ model is developed as per the earlier discussion. The ‘best’ Q p model was obtained with population size of 2,000 individuals and 150 generations with reproduction probability of 0.05, crossover probability of 0.85, mutation probability of 0.1 and with tournament selection. In GP model development it is important to make a tradeoff between accuracy in prediction of Q p and complexity of the model equation which is achieved by proper selection of number of genes and depth of GP tree. In this study optimum result was obtained with maximum number of genes as two and maximum depth of GP tree as four. The developed GP model can be described as Eq. (9) and shown below.

$$ Q_{p} = \exp \left[ {0.037(D - 6.35)} \right]\,\left[ \begin{array}{l} 0.032(L - 130)(S_{u} - 3.4) \\ \times \left[ {0.000035\,\,(L - 130)^{2} + \sin \left( {0.028(S_{u} - 3.4)} \right)} \right] \\ - 19.259\left( {0.02e - 3.625} \right)\left[ {0.037\left( {D - 6.35} \right) - 0.02e} \right] \\ \end{array} \right] + 81.307 $$
(9)

The ‘best’ MARS model has been developed with a six basis functions after several trials with different number of basis functions. Each set of basis functions was used to predict the pile load capacity (Q p ) and their correlation coefficient (R) was calculated. Figure 5 shows the plot of RMSE value versus number of basis functions considered for model generation, though the MARS model performance gets worst when very few number of basis function is used. However, as the number of basis function is increased, the complexity of model also increases; keeping this in mind six basis functions are adopted in the present study.

Fig. 5
figure 5

The performance MARS models on a wide range of number of basis functions

The coefficients of different basis functions produced for the developed MARS model, and the coefficient of intercept generated is presented in Table 1. Hence, model equations can be written using the obtained coefficients and basis functions as presented in Eq. (10) as follows:

$$ \begin{aligned} {\text{Qp }} & = { 68}. 7 5 8 { } + { 2}. 2 2 4 {\text{ h}}\left( {{\text{D }} - { 18}} \right) \, {-}{ 3}. 4 4 1 {\text{ h}}\left( { 18 - {\text{ D}}} \right) \, + \, 0. 9 5 4 {\text{h}}\left( {{\text{L }} - { 13}0} \right) \\ &\quad - { 2}. 9 2 1 {\text{ h}}\left( {{\text{e }} - \, 0} \right) \, + { 2}. 9 9 8 {\text{h}}\left( {{\text{S}}_{\text{u}} {-}{ 7}. 2} \right) \, - 1 1. 4 8 4 {\text{ h}}\left( { 7. 2{-}{\text{ S}}_{\text{u}} } \right) \\ \end{aligned} $$
(10)

where,

$$ {\text{h}}\left( {{\text{D }} - { 18}} \right) \, = { \hbox{max} }\left( {0,{\text{D }} - { 18}} \right) $$
(10a)
$$ {\text{h}}\left( { 1 8 { }{-}{\text{ D}}} \right) \, = { \hbox{max} }\left( {0, 1 8 { }{-}{\text{ D}}} \right) $$
(10b)
$$ {\text{h}}\left( {{\text{L }}{-}{ 13}0} \right) \, = { \hbox{max} }\left( {0,{\text{ L }}{-}{ 13}0} \right) $$
(10c)
$$ {\text{h}}\left( {{\text{e }}{-} \, 0} \right) \, = { \hbox{max} }\left( {0,{\text{ e }}{-} \, 0} \right) $$
(10d)
$$ {\text{h}}\left( {{\text{S}}_{\text{u}} {-}{ 7}. 2} \right) \, = { \hbox{max} }\left( {0,{\text{ S}}_{\text{u}} {-}{ 7}. 2} \right) $$
(10e)
$$ {\text{h}}\left( { 7. 2 { }{-}{\text{ S}}_{\text{u}} } \right) \, = { \hbox{max} }\left( {0,{ 7}. 2 { }{-}{\text{ S}}_{\text{u}} } \right) $$
(10f)
Table 1 The coefficients of basis functions of MARS model

Based on the DENN and BRNN analysis best models were developed with 3 and 2 hidden layer neurons respectively. Model equations for above two models can be written using the obtained weights and biases following Das and Basudhar [13, 14].

As it is important that the efficiency of model should be compared in terms of testing data than that with training data [14], in this study the comparisons of the methods are done on the basis of testing data only. Figure 6 shows the performance of predicted and observed values of lateral load capacity of piles for GP, MARS and ANN (DENN, BRNN) models. There is less scatter of data for the GP and MARS models compared to the other models. Table 2 shows the statistical performance in terms of R, E, AAE, MAE and RMSE for the GP and MARS model along with the results of ANN (DENN and BRNN), Broms and Hansen’s models for both training and testing data set. The developed GP, MARS and DENN models show good generalization in terms of close values of R and E for training and testing data. The OBJ for each of the developed models are evaluated and also presented in Table 2. Higher value of R and lower values of RMSE and AAE result in lower OBJ value, which indicates a more accurate model. Thus, it is evident from the Table 2 that the developed GP model is better than other models in terms of the calculated OBJ.

Fig. 6
figure 6

Comparisons of predicted load capacity of piles by different methods with measured load capacity

Table 2 Comparison of statistical performances of different models

While describing prediction of pile load capacity based on cone penetration test (CPT) [7] have emphasized that other statistical criteria should be used along with the correlation coefficient. Abu-Farsakh and Titi [1] and Das and Basudhar [13] have used the mean (μ) and standard deviation (σ) of ratio of predicted pile load capacity (Q p ) to the measured pile load capacity (Q m ) as important parameters in evaluating different models. The mean (μ) and standard deviation (σ) of Q p /Q m are important indicators of the accuracy and precision of the prediction method. Under ideal condition an accurate and precise method gives the mean value as 1.0 and the standard deviation to be 0. The μ value greater than 1.0 indicates over prediction and under prediction otherwise. In present study, the μ (1.006, 1.032) and σ (0.125, 0.141) of Q p /Q m for the MARS model is very close to those of GP [(1.007, 0.94), (0.090, 0.107)] and DENN [μ (1.018, 0.948), σ (0.106, 0.125)] for training and testing data. The values for BRNN (μ (1.042, 0.942), σ(0.143, 0.196) and other models as also presented in Table 3. The other criterion like cumulative probability of Q p /Q m [1, 13] should also be considered for the evaluation of performance of different models. The ratio Q p /Q m is arranged as per their values in an ascending order and the cumulative probability is calculated from the following equation:

$$ P = \frac{i}{n + 1} $$
(11)

where i is the order number given to the Q p /Q m ratio; n is the number of data points. If the computed value of 50 % cumulative probability (P 50) is less than unity, under prediction is implied; values greater than unity means over prediction. The ‘best’ model is corresponding to the P 50 value close to unity. The 90 % cumulative probability (P 90) reflects the variation in the ratio of Q p /Q m for the total observations. The model with P 90 for Q p /Q m close to 1.0 is a better model.

Table 3 Evaluation of performance of different prediction models considered in this study

Figure 7 shows the cumulative probability plots of Q p /Q m for different methods for both training and testing data. Based on the figure it can be seen that GP, MARS, DENN and BRNN models are very closely following each other. It can also be seen from Table 3 that P 50 values of MARS (1.004, 0.990), DENN (1.012,0.945), BRNN (1.005, 0.896) and GP (1.020, 0.885) models for training and testing data are comparable. whereas the Hansen method (P 50 = 0.542, 0.523) under predicts the pile load capacity and Broms method (P 50 = 1.112, 1.140) over-predicts the same. However, based on the P 90 value GP (1.096, 1.092) model is found to be close to MARS (1.178, 1.256) and DENN (1.156, 1.161) models and better than other models. The lognormal distributions of the Q p/Q m for different models are shown in Fig. 8. Based on the figure it can be seen that GP model is predicting lateral load capacity of the pile within 20 % accuracy level (i.e. Q p /Q m  = 0.8–1.2) better than MARS, DENN, BRNN and other statistical models as the shaded area under the lognormal distribution plot of GP model is more than those of the other models within the above limit.

Fig. 7
figure 7

Cumulative probability plots of Qp/Qm for different methods

Fig. 8
figure 8

Log normal distribution of Qp/Qm for different methods

As per the best fit calculations (R, E) (R1), arithmetic calculations of Q p /Q m (μ, σ) (R2), cumulative probability of Q p /Q m (P 50, P 90) (R3)and prediction of pile load capacity within 20 % accuracy level (R4), a ranking system is made among different models and presented in Table 3. The overall performance of the various models under present study is evaluated using RI as per Abu-Farsakh and Titi [1]. The RI is the sum of the ranks of different models as per the above four criteria (RI = R1 + R2 + R3 + R4). Lower the value of RI indicates better performance of the particular method. It can be seen from the Table 3 that GP model (RI = 5) is ‘best’ among various models used in the present study and is closely followed by MARS model (RI = 8) and other models [DENN (RI = 11), BRNN(RI = 16), Broom’s (RI = 20) and Hansen’s (RI = 24)].

The results of present developed models are also compared with the results of SVM and GPR models as given by Pal and Deswal [36]. However, the SVM and the GPR results are available for the testing data in terms of R and RMSE only. The R values of SVM and GPR models are 0.920 and 0.980 respectively. Similarly, the values of RMSE are 11.47 and 6.32 for SVM and GPR models respectively. Hence, the present GP (0.972, 8.194) and DENN (0.967, 8.549) models are found to be better than the SVM model as per R and RMSE values. The R value of GP (0.972) and MARS (0.980) models are comparable to GPR model, though GPR model is better than above two models in terms of RMSE value. However, due to absence of other data, performance of these two models based on other criteria as discussed in the above paragraph could not be made to make an elaborate comparison using RI.

Sensitivity Analysis

The sensitivity analysis is an important aspect of a developed model to find out important input parameters. In the present study sensitivity analysis was made for ANN (DENN, BRNN) models following [13]. For the developed GP model sensitivity analysis was made according to Gandomi et al. [22]. As per Gandomi et al. [22] the sensitivity (S i ) of each parameter, is expressed by Eq. 12a and 12b.

$$ S_{{_{i} }} = \frac{{N_{i} }}{{\sum\nolimits_{j = 1}^{n} {N_{j} } }} \times 100 $$
(12a)
$$ N_{i} = f_{\hbox{max} } \left( {x_{i} } \right) - f_{\hbox{min} } \left( {x_{i} } \right) $$
(12b)

where f max(x i ) and f min(x i ) are the maximum and minimum of the predicted output over the ith input domain, where the other variables are equal to their mean values. n is the number of variables. In the present study n = 4. Table 4 presents the results of above analyses. The sensitivity analysis for the MARS was done as per numbers of subsets (nsubsets), which is the number of subsets that include the variable, residual sum-of-squares (RSS) of the model, generalized cross validation (GCV) of the model [18] and presented in Table 5. As per ‘best’ model, GP, S u is the most important input parameter. Similar observation is also made by BRNN model (Garson algorithm and Connection weight approach). The other important inputs in descending order are D, e and L (Table 4).

Table 4 Sensitivity analysis of inputs as per different approaches
Table 5 Comparison and importance of various variables according to number of subsets (nsubsets), generalized cross validation (GCV) and residual sum-of-squares (RSS)

Conclusions

The following conclusions can be drawn from the above studies:

  1. (1)

    The proposed GP model is found to be effective and efficient than MARS, ANN (DENN, BRNN), SVM and other statistical models in predicting the lateral load capacity of piles in clay.

  2. (2)

    Using a ranking method based on different statistical criteria (statistical performances for predicted load capacity (Qp) and measured capacity (Qm), the mean and standard deviation of the ratio Qp/Qm, the cumulative probability for Qp/Qm. and prediction of load capacity within 20 % accuracy level) it has also been found that the developed GP model is more efficient compared to other AI and statistical models.

  3. (3)

    The developed model equation is found to more compact compared to the MARS and other AI models and can easily be used by the professionals with the help of a spreadsheet without going into the complexity of model development.

  4. (4)

    Based on sensitivity analysis undrained shear strength of soil is found to be the most important parameter followed by the diameter of pile, eccentricity and length of pile.