Optimization for the Assessment of Spudcan Peak Resistance in Clay–Sand–Clay Deposits

Clay–sand–clay deposits are commonly encountered in the offshore field. For spudcan installation in this soil stratigraphy, the potential for punch-through exists, with the peak penetration resistance formed within the interbedded sand layer. Therefore, a careful assessment of the penetration resistance profile has to be performed. Based on the recently proposed failure-stress-dependent model, this paper presents a modified predictive model for estimating the peak resistance. The modified model incorporates the bearing capacity depth factor and the protruded soil plug in the bottom clay layer into the formulation. It is proven that the modified predictive model provides improved deterministic estimations for the peak resistances measured in centrifuge tests. Based on the modified predictive model, a parameter optimization technique is utilized to optimize the prediction of peak resistance using penetration resistances observed beforehand. A detailed application procedure is proposed and applied to the centrifuge tests accumulated from existing publications, with further improvement on the predictions demonstrated. The proposed parameter optimization procedure combined with the modified predictive model provides an approach to perform real-time optimization for assessing spudcan peak resistance in clay–sand–clay deposits.


Introduction
The increasing demand for hydrocarbon resources and renewable energy is driving offshore explorations into regions where highly layered seabed conditions are prevalent [1,2]. On such soil stratigraphy, the installation of mobile jack-up rigs involves more complexity and uncertainty. This is particularly the case for the assessment of punch-through incident, which refers to a sudden and large penetration of the spudcan foundation attached to the bottom of the jack-up leg due to a negative or near zero gradient of penetration resistance profile.
Clay-sand-clay stratigraphy is recognized as one of the soil profiles commonly encountered in the offshore field with a significant potential for punch-through failure [1,3]. The typical penetration resistance profile of spudcan in such deposits is shown in Figure 1, which is presented as the total penetration resistance q versus the penetration depth d measured from the seabed surface to the lowest point of the largest cross-section of spudcan. As the spudcan penetrates in the top clay layer, either a single-layer response (i.e., the penetration response in single-layer clay, without being affected by the underlying layers) or a squeezing response occurs depending on the distance between the spudcan base and the sand layer [1,4,5]. Once the squeezing response is triggered at the depth of d sq , the penetration resistance will increase sharply until reaching the peak value of q peak in the interbedded sand layer at the depth of d peak . The potential for punch-through accidents should be assessed by comparing the designed preload with q peak . Since the preload is applied through ballast water, which can only be dumped gradually, the spudcan may experience uncontrolled penetration once the preload exceeds q peak . Therefore, accurate prediction of the peak resistance is a vital component for spudcan penetration in clay-sandclay deposits. Recently, deterministic design methods have been proposed by Ullah et al. [4] and Zheng et al. [5] based on centrifuge model tests and large-deformation finite-element (LDFE) analyses. The peak resistance of spudcan in the interbedded sand layer is calculated by modifying the failure-stress-dependent predictive model established by Lee et al. [6] and Hu et al. [7] for a spudcan on sand-over-clay deposits. Both deterministic predictive methods provide reasonable predictions for the peak resistances from centrifuge tests and LDFE analyses. Nevertheless, as will be detailed later, the depth where the peak resistance occurs is underestimated by these design methods. This deficiency leads to uncertainties in the prediction of spudcan peak resistance and should be rectified by modifying the predictive model. Therefore, further improvement is proposed in this study to obtain a more realistic representation of the true soil failure mechanism at punch-through. Although deterministic design methods have been satisfactorily verified against centrifuge tests and numerical analyses, their practical applications can still lead to considerable divergences between predicted and recorded results, even for relatively simple single-layer and double-layer soil stratigraphies [8,9]. This is, among others, mainly due to the uncertainties associated with (a) the error of the predictive model itself and (b) the determination of model parameters as input for the calculation [10][11][12]. The former is induced as a result of simplifications of the physical problem. The latter can be more significant for practical applications than physical modeling tests, as much more uncertainties are involved in the investigation and interpretation of in-situ soil properties. Therefore, a technique that allows fast optimization of the prediction based on penetration resistances measured real-time will definitely favor the assessment of punch-through. Probabilistic methods provide a feasible way to perform real-time optimization for the prediction based on the monitored data. The percentile curves or probabilistic contours of the spudcan load-penetration response can be predicted through Monte Carlo simulations by account-ing for the uncertainties and probabilistic distribution of the inputs in the predictive model [11,[13][14][15]. Bayes' theorem can then be adopted to update these predictions according to the measured load-penetration data [16,17]. The parameter optimization technique (POT) is another alternative option, which is also formulated based on Bayes' theorem, but the prediction is updated by optimizing the model parameters according to the measured data. Most recently, POT was applied to the real-time assessment of spudcan installation in sand over clay deposits [18]. The resistances measured just below the sand surface were used to optimize the strength parameters of the clay layer, which were later input into the predictive model for an optimized prediction of peak resistance. Calibrations against a series of centrifuge tests and two case histories have demonstrated promising performance of the approach proposed by Jiang et al. [18]. However, the applicability of POT for more complex soil stratigraphies, e.g., clay-sand-clay deposits, has yet to be examined. Moreover, the current POT scheme optimizes the predicted peak resistance only after the spudcan penetrates into the sand layer, which does not allow sufficient time for the rig operator to adjust the preload strategy if a punch-through incident is deemed to occur. A more robust POT scheme is hence required.
In this paper, a deterministic predictive model for the calculation of spudcan peak resistance in clay-sand-clay deposits is presented first, which is modified from that proposed by Ullah et al. [4] to rectify the underestimated peak resistance depth. Accordingly, design formulas are proposed to account for the effects of the bearing capacity depth factor and the protruded soil plug in the bottom clay layer in the calculation of peak resistance. After that, a detailed procedure is proposed for the application of POT in clay-sand-clay deposits. By the proposed approach, the prediction of peak resistance in the interbedded sand layer is optimized based on the resistances measured in the overlying clay layer. The feasibility and performance of the proposed POT are investigated based on a series of existing centrifuge tests.

Soil Parameters
A schematic diagram of a spudcan of diameter D penetrating in clay-sand-clay deposits is illustrated in Figure 1. The sand layer with effective unit weight γ s and thickness H s , is interbedded by clay layers. The top clay layer has an effective unit weight of γ ct , undrained shear strength of s uts at the mudline, strength gradient of k t and thickness of H ct . The bottom clay layer has an effective unit weight of γ cb , undrained shear strength of s ubs at the sand-clay interface, strength gradient of k b and nominally infinite depth.

Modified Predictive Model
The predictive model proposed by Ullah et al. [4] has been validated against the centrifuge tests reported by Ullah et al. [19] and Hossain [1]. Based on the geometries and soil parameters summarized in Table 1, the peak resistances can mostly be predicted with an error bounded by ±20%, as shown in Figure 2 (solid grey circles) in terms of the ratio r between predicted and measured peak resistances. For the conciseness of the narration, the mathematic framework of Ullah et al.'s [4] method is provided in Appendix A. The predictive model proposed in this study is modified from Ullah et al. [4] so that the soil failure mechanism is modeled more realistically. The modifications include four aspects, which are detailed as follows.  Figure 2. Comparison between measured and predicted peak resistances.

Distribution Factor D F
Zheng et al. [20] reported that the peak resistances were significantly overestimated for sand overlying clay with high undrained shear strength using the distribution factors estimated by Equation (A3). They have hence collected centrifuge test data from existing publications, carried out additional centrifuge tests and LDFE analyses, and proposed a new design formula for D F , which is expressed as The new formula covers a more comprehensive range of soil properties and geometries, and provides an improved performance for the prediction of peak resistance in sand-over-clay deposits by Hu et al.'s [7] method, based on which Ullah et al.'s [4] model was developed. Therefore, it is also adopted in the modified predictive model of this study.

Peak Resistance Depth d peak
By plotting the measured d peak from all the related centrifuge tests reported by Ullah et al. [19] and Hossain [1] against the predictions from Equation (A5) of Appendix A, as shown in Figure 3, it is found that the formula tends to underestimate d peak . The measured values of d peak from centrifuge tests can be better fitted by For H ct = 0, Equation (2) reduces to d peak = 0.12H s for sand-over-clay deposits. The comparison between the measured and predicted d peak for Equation (2) is also plotted in Figure 3, which shows improvement in the estimations.

Soil Resistances from the Sand Plug
With Equation (2) being used and with the assumption that the thickness of the trapped clay frustum is H c = 0.07H ct following Ullah et al. [4,19], the effective sand layer height measured from the top of the sand plug to the original sand-clay interface, H eff , is changed from 0.88H s to By substituting Equation (3) into Equation (A1) with the rest of the design formulas, the same as Ullah et al. [4], the resulting peak resistance ratios are shown as crosses in Figure 2. It is seen that most of the peak resistances are overly underestimated, with the largest error of −44.5%. This is because, for d peak from Equation (2), the spudcan has penetrated into the sand layer by a depth that the sand-bottom clay interface can no longer be regarded as undeformed. Instead, this deformed interface, which forms the bottom end of the soil plug, protrudes into the bottom clay layer. This is evidenced by Figure 11 of Ullah et al. [19] and Figures 11 and 17a of Hossain [1], which demonstrate the digital images taken from half-model centrifuge tests. Additional resistance is provided by the soil plug below the original sand-bottom clay interface, which contributes an important portion of the peak resistance.
Based on the above failure mechanism, the predictive model by Ullah et al. [4] is modified to incorporate the effects of the soil plug below the original sand-bottom clay interface, as illustrated in Figure 4a. According to the equilibrium of the thin disc element, Equation (A1) can still be used to calculate q peak . However, the term N c0 s ubs , which is equal to the vertical stress within the sand frustum at z = H eff , has to be modified to consider the base and frictional resistances provided by the sand plug with a height of H plug3 in the bottom clay layer. As such, the vertical equilibrium is established according to Figure 4b for the sand frustum in the bottom clay layer, with the assumption that the shear stress on the periphery is constant and equal to the average undrained shear strength s u,plug3 of the clay around the frustum which can be transformed to In Figure 4, D 1 , D 2 , and H plug3 are the upper base diameter, lower base diameter, and height of the sand frustum in the bottom clay layer, respectively; and q plug is the vertical resistance acting on the lower base of the frustum due to soil shearing. To determine these parameters, Test T3SP, as reported by Ullah et al. [19], was reproduced using the LDFE approach in Zheng et al. [5]. Similar to the finding by Zheng et al. [5], a reasonable agreement between the measured and computed penetration resistance profiles is obtained, as shown in Figure 5a. As for the soil plug heights, which are compared in Figure 5b, the total plug height of about 0.05H ct + 0.9H s is found from the numerical analysis, close to 0.07H ct + 0.9H s reported and assumed by Ullah et al. [4,19]. Therefore, a trapped clay height of H c = 0.07H ct and a sand plug height, which is the summation of H eff and H plug3 , of 0.9H s are also assumed in this study, and hence the parameters in Equation (5) can be calculated as: The plug end resistance q plug can be calculated by multiplying the bearing capacity factor by the corresponding undrained shear strength. The design formula recommended by ISO [21], which considers the embedment depth and soil non-homogeneity after Skempton [22] and Young et al. [23], respectively, is adopted to calculate q plug as The effects of the bearing capacity depth factor and protruded soil plug in the bottom clay layer are captured through Equation (5) in combination with Equations (6)- (10). The modified predictive model hence estimates the peak resistance according to Equations (A1) and (A2) with D F from Equation (1), H eff from Equation (3) and N c0 s ubs from Equation (5).

Non-Linear Load-Penetration Response before Peak
The squeezing response was simplified as a straight line by Hu et al. [7] and Ullah et al. [4]. This is acceptable in the deterministic models as the penetration resistance prior to the peak is non-crucial to the assessment of punch-through. However, this curve has to be described appropriately for the application of POT in optimizing the prediction of peak resistance [18]. A non-linear curve is thus proposed to model the load-penetration response between d = 0.9H ct and d peak , as shown in Figure 6. According to the centrifuge test data reported by Ullah et al. [19] and Hossain [1], it is found that an exponential function may be used to approximate the non-linear response as where q s is the penetration resistance at d = 0.9H ct and A is a parameter that is best fitted as A = 1.04. The comparison between Equation (11) and the centrifuge data is plotted in Figure 7. For deterministic predictions, the value of q s can be estimated based on a linear interpolation between the resistances at d sq and d peak .

Application of POT for Spudcan in Clay-Sand-Clay Deposits
The methodology of POT is based on the data assimilation theory, which is a particular application of Bayes' theorem [24]. It provides a real-time optimization for the prediction by "adjusting" the model parameters based on the observational information of the variables (e.g., penetration resistance in this study). The technique has been widely incorporated in the forecast models of earth systems [25][26][27][28][29][30][31][32]. The steps for the implementation of POT are summarized in Appendix B.

POT Scheme for q peak in Clay-Sand-Clay
Optimized estimation of q peak in sand over clay deposits was proposed by Jiang et al. [18], in which POT was performed based on the penetration resistances observed at depths of 0-0.025H s below the sand surface. A straight line was used to simplify the load-penetration response above d peak . However, this may not be appropriate for clay-sand-clay deposits since the load-penetration response before peak resistance is much more non-linear, with the gradient reducing as the spudcan approaches d peak . Therefore, a more reasonable description, e.g., Equation (11), should be adopted. Moreover, the observation locations selected by Jiang et al. [18] are below the sand surface and close to d peak , while an earlier optimization of the predicted peak resistance is preferred, which allows the rig operator to make a judgment as early as possible about whether a punch-through incident will occur or not for spudcans under the designated preload.
For the deterministic model of spudcan in clay-sand-clay deposits introduced in Section 2.2, it is suggested in this study that the observation locations above the top clay-sand interface are adopted so that the prediction of q peak can be optimized before punch-through occurs. Here, observation locations at d = 0.925, 0.95, 0.975, and 1.0H ct are selected to verify the performance of the proposed POT in Section 4. With q s measured at d = 0.9H ct and the penetration resistance profile estimated from Equation (11) of the modified predictive model, POT can then be performed at the selected observation locations. As justified by Jiang et al. [18], among the model parameters that are determined with relatively significant uncertainties, the most dominant parameter for q peak in sand overlying clay is the undrained shear strength at the sand-clay interface. A similar phenomenon was found in the trial calculations for clay-sand-clay deposits. Therefore, only s ubs is adopted as the parameter to be optimized in order to rectify all the uncertainties associated with the prediction of q peak .

Procedure of Application
For a modern jack-up rig with the ability to monitor applied vertical load and corresponding spudcan penetration depth during the installation stage, the modified predictive model can be combined with the POT scheme to provide real-time optimization for the prediction of spudcan peak resistance in clay-sand-clay deposits. The detailed procedure is proposed and illustrated in Figure 8, which is summarized as follows: (1) Carry out an assessment of punch-through before jack-up installation by making a deterministic prediction of q peak according to the modified predictive model, i.e., Equations (A1), (A2), (3) and (5)-(10) with D F from Equation (1). Note that D F from Equation (1) is only applicable for 0.16 ≤ H s /D ≤ 1.0. If the designated preload is much lower than the predicted q peak , e.g., less than 0.75q peak , POT may not be necessary as the spudcan is highly likely to rest in the top clay layer. However, caution should still be exercised, especially when the spudcan continues to approach the top clay-sand interface with the increase of ballast water.  Due to the promising convergent rate of POT, the above procedure can be implemented within seconds through, e.g., object-oriented programming languages. This enables realtime optimization in practice for the prediction of the spudcan peak resistance in claysand-clay deposits.

Performance of Modified Predictive Model
The peak resistance ratios resulting from the deterministic predictions of the modified predictive model are plotted as triangles in Figure 2 for the centrifuge tests in Table 1. At first glance, both the modified predictive model and the original one by Ullah et al. [4] show similar performances, resulting in prediction errors mostly within ±20%. However, the statistics summarized in Table 2 indicate an improved performance of the modified predictive model, with the SD and mean absolute error (MAE) reducing from 0.126 to 0.094 and 9.7 to 8.3%, respectively. The peak resistance ratios resulting from the modified predictive model also fell into a narrower range of 0.77-1.101, compared with the range of 0.732-1.267 from Ullah et al.'s [4] approach. This is because the modified predictive model is a more realistic representation of the true soil failure mechanism at peak resistance. Additionally, it also provides the basis for the application of POT in clay-sand-clay deposits.

Performance of POT
Following the procedure developed in Section 3.2, the proposed POT is used to optimize the prediction of peak resistance for spudcan penetration in clay-sand-clay deposits. Although there are a number of case histories reported by the joint industry project InSafeJIP [3], the recorded load-penetration data were either terminated before the peak resistance or incomplete for the application of POT. Therefore, only centrifuge tests listed in Table 1 are considered, which are all the centrifuge tests available in the existing publications. As mentioned previously, the observation locations in each test are selected at d/H ct = 0.925, 0.95, 0.975, and 1.0.

Overall Performance
After probabilistic analyses with POT at all observation locations, the absolute errors of the peak resistance ratios are plotted in Figure 9 to compare with results from the corresponding deterministic predictions. Improved predictions of q peak are observed for most cases, especially for those with relatively large errors from deterministic predictions. However, a relatively obvious increase in the absolute error is also observed in three centrifuge tests, i.e., Tests SPb6, SPb8, and SPb12, with |r − 1| increasing from 3.7, 3.8, and 1% to 19.7, 9.1, and 9.7% respectively. This may be because the modified predictive model is still somewhat biased due to simplifications (e.g., the idealized shape of the soil plug and the simplified representation of the load-penetration squeezing response), despite the modifications proposed in this study. Nevertheless, the absolute errors mostly drop below 10% after the optimization of POT. The overall performance of POT is quantified in Table 3, where the statistics for the peak resistance ratios after the optimization at all observation locations are listed. Compared with the data in Table 2, the overall performance is further improved, with the SD and MAE of the peak resistance ratios reducing from 0.094 and 8.3% to 0.075 and 5.9%, respectively. Table 3. Statistics for peak resistance ratios after parameter optimization.

Predictive Model
Design Formulas Peak Resistance Ratio, r MAE * Min.
Max. SD

Discussions
Four centrifuge tests, i.e., Tests SPc16, SPb16, T6SP, and SPb6 reported by Ullah et al. [19], are taken as examples to illustrate the optimization process and effectiveness. For each case, the optimized prediction of the penetration resistance profile at each observation location are compared with the centrifuge test result and deterministic prediction in Figure 10. The variations of the ensemble spread (i.e., SD), ensemble mean and experimental value for s ubs , and the measured and predicted penetration resistances at different observation locations are listed in Tables 4-7. There are some interesting points to note in Figure 10 and Tables 5-7, as discussed below.
While the proposed approach does not work perfectly for all cases, Tests SPc16 and SPb16 are representative cases in terms of the variations of optimized resistance profile and parameter ensemble with the observation location. As shown in Figure 10a,b, the optimized penetration resistance profile evolves towards the measured profile along with the spudcan penetration after the first observation location, providing greater confidence in the assessment of punch-through. The best performance is observed for Test T6SP (Figure 10c), where the absolute error reduces from 23.0% to 2.4% after the optimization at all observation locations. This may reflect the relatively large experimental error in the measurement of s ubs for Test T6SP. For these three cases where POT leads to an improved prediction, the statistics in Tables 4-6 demonstrate a consistent trend: after the first observation location, the variation of s ubs is evolutionary (by less than 5% after each optimization), and the lower degree of variation is associated with better performance; the optimized q at each observation location is close to the measured value, with errors less than 1%.
However, it should be noted that for cases where the POT has resulted in an improved prediction of q peak , the optimized resistance profile does not necessarily evolve towards the measured data all the way as the spudcan approaches the d peak . For example, in Test SPb16 (Figure 10b), the optimization of q peak at d = 0.925H ct results in a profile that is even further away from the measured data, although the following optimized profiles gradually get closer to the measured data as the spudcan penetrates further and leads to a smaller predictive error after all observation locations. Test T6SP is another example (Figure 10c), where the discrepancy between the optimized and measured penetration resistance profiles gets larger after each observation location. Nevertheless, the optimized penetration resistance profile gets very close to the measured data immediately at the first observation location (i.e., d = 0.925H ct ), and its changes in the following observation locations are quite negligible (less than 2%).
The worst performance of POT is observed for Test SPb6 (Figure 10d). Table 7 shows that at the first observation location of d = 0.925H ct , there is a sharp change in s ubs from 24.6 to 7.7 kPa and the predictive error of q peak changes from 3.7% to −31.6%, although the optimized prediction of q becomes almost identical to the observed value with the updated s ubs . At subsequent observation locations, the optimized resistance profiles evolve towards the measured data, but the variation in s ubs is still significant (by~20% after each optimization), and the gap between the optimized and observed values of q also grows with the penetration depth, with a difference of about 9% (461.9 vs. 507.8 kPa) at the last observation location. This is in contrast to the cases that show an improved prediction with POT, where the variation in s ubs is evolutionary and the optimized and observed values of q only have a minimal difference. The main reason may be because the simplified non-linear curve is not a good representation for the squeezing response of Test SPb6.
The different patterns of variations discussed above lead to the following implications for the practical application of the proposed POT approach: (a) an evolutionary variation of the optimized s ubs (e.g., by <5%) and a close agreement between the optimized and measured q values usually indicate an improved prediction of q peak ; and (b) extreme cautions should be taken during the optimization process if the variation of the optimized s ubs is significant (e.g., by >5%) and there is an increasing error between the optimized and measured q values, which, although rare, may indicate a worsened prediction by POT.

Conclusions
This paper describes a modified predictive model for assessing the peak resistance of spudcan in clay-sand-clay deposits. Model modifications mainly include the incorporation of the effects of the bearing capacity depth factor and the protruded soil plug in the bottom clay layer. The modified model provides a more realistic representation for the true soil failure mechanism at punch-through, and the verification against existing centrifuge tests shows an improved performance compared with the original model. Based on the modified predictive model, the parameter optimization technique has been employed to optimize the prediction of peak resistance in clay-sand-clay deposits using the penetration resistances observed at locations above the top clay-sand interface. A detailed procedure for applying the parameter optimization technique has been summarized and then applied to the centrifuge tests. Further improvement on the prediction of peak resistance has been demonstrated compared with the deterministic predictions from the modified predictive model, which confirms the effectiveness of the proposed approach. Discussions have been made on the optimization process in terms of the variations of the parameter ensemble and the optimized penetration resistances. Based on the discussions, some guidelines for the application of the proposed approach have been suggested.

Appendix A. The Mathematic Framework of Predictive Model for Ullah et al.'s Method
In the design approach developed by Ullah et al. [4], the penetration resistance profile of a spudcan in clay-sand-clay deposits is simplified as four segments of straight lines, as illustrated in Figure 1. The approach has been calibrated against centrifuge tests. For spudcan in the top clay layer, the penetration resistances at depths above d sq are estimated using the bearing capacity factors reported by Houlsby and Martin [33] for single-layer clay. The penetration resistance profile in the bottom clay layer is calculated with a deep bearing capacity factor, which is expressed as a function of the thicknesses of the overlying clay and sand layers.
The peak resistance focused on in this study is formulated based on a punching shear mechanism (see Figure 4 of Ullah et al. [4]). Therefore, q peak mainly comprises the contributions from the frictional resistance and the end bearing capacity of the layered soil plug (or soil frustum) that consists of trapped clay and sand. According to the equilibrium of an element disc in the sand frustum, a general expression for q peak can be given as where N c0 is the bearing capacity factor of clay at the sand-bottom clay interface (i.e., z = H eff ); H eff is the effective sand layer height measured from the top of the sand plug to the original sand-clay interface; q 0 is the surcharge on the soil surface; V f is the volume of spudcan embedded by the soil; and ψ is the dilation angle of sand. The parameter E* is expressed as where φ* is the reduced friction angle calculated from the operative friction and dilation angles following Drescher and Detournay [34], while the operative friction and dilation angles are determined through iterative calculations based on the strength-dilatancy re-lationships modified from those proposed by Bolton [35]. D F is the distribution factor following Hu et al. [7], which was proposed for spudcan penetration in sand-over-clay deposits and expressed as a function of H s /D In the first term of Equation (A1), N c0 s ubs calculates the resistance at the sand-clay interface due to soil shearing in the bottom clay layer. Ullah et al. [4] suggested estimating N c0 following the lower bound solution of a flat rough-based surface footing (i.e., no depth factor being considered) reported by Houlsby and Martin [33], which is expressed as N c0 s ubs = (6.34 + 0.56k b D/s ubs )s ubs (A4) The last bracketed term of Equation (A1) calculates the contributions from the frictional resistance and weight of the trapped clay frustum, and the weight of the backfill soil. The height and average undrained shear strength of the trapped clay frustum are H c and s u,Hc , respectively, and the height of the backfill soil is H BF . Based on the observations from 11 centrifuge tests on half-model spudcan and flat-based foundations, Ullah et al. [4,19] reported H eff = 0.88H s and H c = 0.07H ct . Therefore, the depth of peak resistance is determined as By assuming s u,Hc = the average strength of the top clay layer and H BF = 0.5H ct , Ullah et al. [4] predicts q peak using Equations (A1) and (A2) with D F from Equation (A3), H eff = 0.88H s and N c0 s ubs from Equation (A4).

Appendix B. Steps for Implementation of POT
The implementation of POT can be summarized as follows: (1) Generate a number of random values that follow a standard normal distribution for each of the model parameters to be optimized, which serve as the parameter ensemble. (2) Carry out probabilistic predictions based on the predictive model using the parameter ensemble generated in the last step, and the predicted results form the model ensemble. (3) Calculate the observational increment ∆V i O for each member of the model ensemble based on the relationship between the observed data and the model predictions at the observation location, i.e., where V represents the observable variable at the observation location; σ V is the standard deviation of V; the superscript i denotes the ith ensemble member; the subscripts O and M refers to "observation" and "model", respectively; and the overbar denotes the ensemble mean. where P denotes the parameter to be optimized, and hence ∆P i is the increment of the ith member of the parameter ensemble, while cov(P M , V M ) refers to the covariance between the parameter ensemble and model ensemble. (5) The parameter ensemble is updated by adding ∆P i to the corresponding ensemble member. With the mean value of the updated parameter ensemble, an optimized prediction can be obtained.