基于静力触探与机器学习的打入桩竖向承载力预测

摘要

分析了单桩承载力的主要影响因素，利用降维算法对静力触探数据进行预处理，使用树模型的重要性分析函数筛选主要特征，由贝叶斯优化算法确定超参数，提出了一种基于机器学习算法的打入桩竖向承载力预测模型。经过准度测试、蒙特卡洛模拟检验及与经验关系模型的对比，对预测模型的有效性和鲁棒性进行了评估。结果表明，支持向量机（SVM）模型、随机森林（RF）模型、极限梯度提升算法（XGBoost）模型均能较准确地预测打入桩竖向承载力，且准确性显著高于经验关系模型。根据综合分析结果，建议选择XGBoost模型用于打入桩竖向承载力的预测。

关键词

静力触探; 打入桩承载力; 机器学习; 降维算法; 极限梯度提升算法

打入桩具有成本低、沉降少、对周边环境影响小等优点，在基础工程中已得到广泛使用，其竖向承载力的准确确定对桩基的设计和施工至关重要。静载试验是工程中直接确定竖向承载力的常用方法，但其人力物力成本较高。如何精准经济地确定承载力一直是业内学者的重要研究课题^［1］。

静力触探（CPT）是一种高效的原位测试手段，由于CPT探头贯入过程与打入桩入土过程相似，国内外学者基于CPT数据和大量工程实践经验提出了许多预测打入桩单桩承载力的经验方法^［2-4］。经验方法简单便捷，但考虑的影响因素较少，包含较多无明确意义的修正系数，可靠性较差。鉴于此，国内外学者对CPT数据与桩承载力的关系开展了深入研究。李洪江等^［3］建立了孔压静力触探（CPTU）锥尖参数与水平荷载‒位移曲线（p‒y曲线）参数的综合对应关系，提出了一种高效确定桩基水平承载力的CPTU测试p‒y曲线法。李镜培等^［4］结合圆孔扩张理论建立了静压沉桩和CPTU贯入过程中土压力和超孔压的预测方法。Miao等^［5］提出一种基于CPTU数据预测静压桩单桩承载力的改进方法，通过改进的一阶二阶矩法得到的可靠性指标证明该方法比常规经验方法更加可靠。Moayedi等^［6］基于CPT数据对钻孔灌注桩荷载‒沉降响应进行研究，分别采用静态和动态神经网络建立了高精确度的人工智能模型，并给出了计算神经网络参数的简单计算公式。Singh等^［7］采用M5模型树、随机森林回归和神经网络建模方法，结合室内实验预测了群桩的斜荷载承载能力。Ebrahimian等^［8］采用进化多项式回归（EPR）提出了一种基于CPT数据的轴向桩承载力预测模型，该模型的训练使用了较为全面的数据库。Ghorbani等^［9］探讨了自适应神经模糊推理系统（ANFIS）在利用CPT数据预测桩极限轴向承载力方面的潜力。由于数据集的局限性，前人方法往往不能同时考虑桩的两种加载模式和泛用性，且在对CPT原始数据的处理上采用取均值法，此法会损失大量信息。

为建立广泛适用且可靠性高的打入桩极限承载力预测模型，本文基于浙大桩基数据库^［10］，首先对单桩竖向承载力的主要影响因素进行分析，再采用主成分分析降维算法处理CPT数据以解决经验方法中的信息损失问题。在采用树模型的重要性分析函数初步筛选主要特征后，利用3种主流机器学习算法建立预测模型，通过与经验关系模型的对比及对蒙特卡洛模拟结果的分析，对机器学习模型的有效性和鲁棒性进行验证。

1 单桩承载特性

1.1　单桩承载力影响因素

单桩的极限承载能力和很多因素有关，如桩长、桩径、桩周土的物理力学性质、入土深度等。此外，工程实践中发现桩的几何形状、材料类型，成桩过程中的成桩工艺、加荷速率，施工工艺、施工方法等也会影响承载力。这些影响因素和桩的承载力之间存在着复杂的非线性关系^［1］。打入桩可分为多种不同类型，如钢管桩、预制管桩、实心预制桩等，不同桩型桩壁厚度、土塞高度不同，其桩端承载力也存在一定差异。可根据桩体材料和铸壁厚度等桩身参数来表征上述不同桩型的差异。

本文从两方面考虑影响打入桩竖向承载力的主要因素，即桩周土的物理力学性质参数和桩身参数，前者包括地下水位深度、场地土物理力学性质（用CPT数据表征），后者包括桩身形状、桩体材料、桩长、桩外径、铸壁厚度、加载模式。

1.2　桩端影响区范围

基于CPT数据的经验方法采用不同的端阻影响区来确定锥尖阻力（$q_{\mathrm{c}}$）的均值，导致承载力计算结果差异较大。以美国某场地的抗压打入桩为例，分别采用《建筑桩基设计技术规范》^［2］（以下简称《桩基规范》）、《孔压静力触探测试技术规程》^［11］（以下简称《静探规程》）、LCPC法^［12］计算其所取$q_{\mathrm{c}}$均值和竖向承载力，桩外径为$d$，其端阻影响区分别为桩端平面$+\mathrm{4}d~-d$、$+\mathrm{4}d~-\mathrm{4}d$、$+\mathrm{8}d~-\mathrm{4}d$（“$+$”、“$-$”分别代表桩端平面往上和往下），如图1所示。

图1  不同经验方法的端阻影响区范围

Fig. 1  Range of end resistance influence zone of different empirical methods

根据图1，《桩基规范》、《静探规程》、LCPC法所取的$q_{\mathrm{c}}$均值分别为25.8、23.6、24.8 MPa，计算得桩端承载力分别为584.1、1 287.4、923.6 kN，竖向承载力$Q_{\mathrm{m}}$分别为1 795.8、3 379.3、2 515.4 kN。说明采用经验方法确定的单桩竖向承载力受到桩端阻力影响区范围的影响。

机器学习模型可以输入全部静力触探数据，并利用降维算法获取主成分以减少信息损失，可以最大限度利用现场试验数据，同时规避CPT深度不满足端阻影响区范围的问题。

2 桩基数据库与数据预处理

2.1　数据库

本文使用的数据来源于浙江大学搜集的全球52个场地及其116根打入桩（ZJU‒ICL）^［10］。根据1.1节中对桩基承载力影响因素的分析，初步选择地下水位深度、CPT数据、桩身形状、桩体材料、加载模式、桩长、桩外径、铸壁厚度共8个特征作为模型的输入参数，单桩极限承载力$Q_{\mathrm{m}}$作为模型的输出参数。

对桩身形状、桩体材料、加载模式等基数标量标签化处理以量化特征，如表1所示。

表1  部分特征的标签化处理

Tab. 1  Tagging of some features

特征	基数标量	标签值
桩身形状	圆桩	0
	方桩	1
桩体材料	混凝土桩	0
	钢桩	1
加载模式	抗压桩	0
	抗拔桩	1

其余输入参数中，地下水位深度、桩长、桩外径、铸壁厚度是比率标量；原始CPT数据是沿深度的矩阵，采用降维算法对其进行处理。

输出参数中，注意到No.40场地的某抗压桩$Q_{\mathrm{m}}$为34 680 kN，超出ZJU‒ICL数据库桩承载力平均值1 000%，对预测模型而言为异常值，舍去。

2.2　CPT数据的降维处理

2.2.1　锥尖阻力$q_{\mathrm{c}}$和侧壁摩阻力$f_{\mathrm{s}}$的横向降维

ZJU‒ICL数据库52个场地中有21个场地含有$q_{\mathrm{c}}$和$f_{\mathrm{s}}$数据，31个场地仅有$q_{\mathrm{c}}$数据。对这21个场地的$q_{\mathrm{c}}$和$f_{\mathrm{s}}$进行多重共线性分析，计算其皮尔逊（Pearson）相关系数，如图2所示。

图2  各场地双桥静探数据相关系数

Fig. 2  Correlation coefficient of CPT data at each site

图2中样本量最少的场地为No.18，含有71个样本，其自由度为69，显著性因子取$\alpha =\mathrm{0.05}$，根据线性相关显著性检验临界相关系数表可得最高临界相关系数为0.232。21个场地的双桥静探数据相关系数平均值为0.738，最大值为0.955，最小值为0.458，远高于对应临界相关系数值，认为总体相关关系显著。Eslami等^［13］认为将非独立的变量用于目标值的表征违反了数据表示的基本原则。从数理模型的角度看，解释变量之间存在高度相关关系也易使预测模型估计失真，故对$q_{\mathrm{c}}$和$f_{\mathrm{s}}$数据进行横向降维处理，形成融合数据$q_{\mathrm{f}}$。

主成分分析法（PCA）是解决多重共线性的主流数据降维算法，可以最大程度地保留数据信息。在线性变换过程中，计算数据协方差矩阵的特征向量，将其标准正交化，并按特征值的大小排序，前若干个特征值对应的数据变量方差最大，包含信息最多，其特征向量称为主成分，一般主成分的累积信息保留程度达85%以上即可^［14］。

以No.35场地为例，该场地CPT探测深度为20.6 m，探测精度为0.2 m，该组CPT数据包含103个样本，每组样本包含$q_{\mathrm{c}}$和$f_{\mathrm{s}}$，即

$${\mathit{X}}^{\mathrm{T}}=\left[{\mathit{x}}_{\mathrm{1}}^{\mathrm{T}},\cdots ,{\mathit{x}}_{\mathrm{103}}^{\mathrm{T}}\right]$$

（1）

式中：${\mathit{X}}^{\mathrm{T}}$为样本集，每一组CPT样本${\mathit{x}}_i\in {\mathbf{R}}^{\mathrm{2}}$。将$q_{\mathrm{c}}$和$f_{\mathrm{s}}$降成一维$q_{\mathrm{f}}$，设降维后数据为${\mathit{y}}_i$，标准正交基的集合为$\mathit{W}$，则${\mathit{y}}_i=\mathit{W}{\mathit{x}}_i$，其中${\mathit{y}}_i\in {\mathbf{R}}^{\mathrm{1}}$，$\mathit{W}$为

$$\mathit{W}=\left[{\mathit{w}}_{\mathrm{1}}^{\mathrm{T}},\cdots ,{\mathit{w}}_{\mathrm{103}}^{\mathrm{T}}\right]\in {\mathbf{R}}^{\mathrm{1}\times \mathrm{2}}$$

（2）

设样本的均值为$\overline{\mathit{x}}=\frac{{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{103}}}{\mathit{x}}_i}{\mathrm{103}}$，样本的投影的均值为$\mu =\frac{{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{103}}}{\mathit{w}}^{\mathrm{T}}{\mathit{x}}_i}{\mathrm{103}}$，则样本的协方差矩阵为

$$\mathrm{C}\mathrm{o}\mathrm{v}\left(\mathit{X}\right)=\frac{{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{103}}}({\mathit{x}}_i-\overline{\mathit{x}})({\mathit{x}}_i{-\overline{\mathit{x}})}^{\mathrm{T}}}{\mathrm{103}}$$

（3）

则样本投影后的方差为

$${\sigma }^{\mathrm{2}}=\frac{{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{103}}}{\left({\mathit{w}}^{\mathrm{T}}{\mathit{x}}_i-\mu \right)}^{\mathrm{2}}}{\mathrm{103}}={\mathit{w}}^{\mathrm{T}}\mathrm{C}\mathrm{o}\mathrm{v}\left(\mathit{X}\right)\mathit{w}$$

（4）

要使得样本投影后的方差最大，则对于正交基$\mathit{w}$要求$\widehat{\mathit{w}}=\mathrm{a}\mathrm{r}\mathrm{g}\underset{\mathit{w}}{\mathrm{m}\mathrm{a}\mathrm{x}}{\mathit{w}}^{\mathrm{T}}\mathrm{C}\mathrm{o}\mathrm{v}\left(\mathit{X}\right)\mathit{w}$，加上约束条件${\mathit{w}}^{\mathrm{T}}\mathit{w}=\mathrm{1}$，故可通过拉格朗日乘数法构造一个目标函数：

$$\mathit{L}\left(\mathit{w},\lambda \right)={\mathit{w}}^{\mathrm{T}}\mathrm{C}\mathrm{o}\mathrm{v}\left(\mathit{X}\right)\mathit{w}+\lambda \left(\mathrm{1}-{\mathit{w}}^{\mathrm{T}}\mathit{w}\right)$$

（5）

式中：$\mathit{L}\left(\mathit{w},\lambda \right)$为PCA降维时的目标函数，$\lambda$为拉格朗日乘数。

对目标函数求偏导且令偏导为0，可得

$$\mathrm{C}\mathrm{o}\mathrm{v}\left(\mathit{X}\right)\mathit{w}=\lambda \mathit{w}$$

（6）

观察式（6）可知$\lambda$是$\mathrm{C}\mathrm{o}\mathrm{v}\left(\mathit{X}\right)$的特征值；$\mathit{w}$是对应的特征向量，将所有的特征值排序，选择最大的若干个特征值对应的特征向量，即得到$\mathit{W}$。

特征间相关度越高，PCA降维后信息损失越小。经过计算，各场地$q_{\mathrm{c}}$和$f_{\mathrm{s}}$数据横向降至一维后信息保留程度在99.55%至99.99%之间，如图3所示。

图3  CPT降维的信息保留程度

Fig. 3  Information retention of CPT dimension reduction

2.2.2　沿深度的纵向降维

将21个场地的$q_{\mathrm{f}}$值与31个场地的$q_{\mathrm{c}}$值沿深度方向用PCA算法降至低维，得到信息高度浓缩的CPT数据$c_{\mathrm{p}}$，其信息保留程度图如图3所示。

从图3中可见，CPT数据维数低于4时，信息保留程度将低于85%，考虑到休斯现象及PCA信息保留要求，取纵向降维到4维。

以No.35场地为例，其CPT数据在横纵向降维后的结果如表2所示，其中$q_{\mathrm{f}}$和$c_{\mathrm{p}}$量纲为1。

表2  ${\mathit{q}}_{\mathbf{c}}$和${\mathit{f}}_{\mathbf{s}}$横向和纵向降维结果示例

Tab. 2  Examples of horizontal and vertical dimension reduction results of ${\mathit{q}}_{\mathbf{c}}$ and ${\mathit{f}}_{\mathbf{s}}$

深度/ m	$q_{\mathrm{c}}$/ MPa	$f_{\mathrm{s}}$/ MPa	$q_{\mathrm{f}}$	$c_{\mathrm{p}}$
				$c_{\mathrm{p}1}$	$c_{\mathrm{p}2}$	$c_{\mathrm{p}3}$	$c_{\mathrm{p}4}$
0.20	3.40	0.02	-13.85	-74.02	-16.21	-22.27	0.47
0.40	5.61	0.04	-11.64
$⋮$	$⋮$	$⋮$	$⋮$
20.40	34.98	0.31	17.73
20.60	38.79	0.41	21.54

2.3　特征工程

为降低模型复杂度并提高泛用性，使用树模型的特征重要性分析函数对输入参数进行初步筛选。以在树模型中作为划分属性的次数为重要性计算指标，11个输入参数的重要性分数如图4所示。

图4  模型特征重要性分数

Fig. 4  Importance scores of model features

由图4可知，作为基数标量的3个特征中，加载模式对预测模型的重要性远超桩体材料和桩身形状，故仅保留加载模式；作为比率标量的4个特征中，桩长有着最高的重要性分数，其次是桩外径、地下水位深度，铸壁厚度的重要性可忽略不计，舍去；CPT数据降维后的四维数据中，$c_{\mathrm{p}\mathrm{1}}$作为第一维保留着绝大部分的CPT数据集信息，且对模型的贡献度超过其他三维数据的总和，故仅保留$c_{\mathrm{p}\mathrm{1}}$。经过特征选择，暂时保留桩长、加载模式、桩外径、$c_{\mathrm{p}\mathrm{1}}$、地下水位深度共5个特征，后续通过蒙特卡洛模拟对是否保留地下水位深度做出进一步判断。

3 算法模型

本文采用3种主流的机器学习算法对单桩承载力$Q_{\mathrm{m}}$进行预测计算，分别是梯度下降算法的代表支持向量机（SVM）、梯度提升算法的代表随机森林（RF）、梯度提升算法的新秀XGBoost。

3.1　支持向量机模型

支持向量机（SVM）的基本模型是定义在特征空间上的间隔最大的线性分类器，通过使用核函数，支持向量机回归（SVR）可改变原数据维度以在新空间中进行线性回归^［15］。

对于$\{(x_{\mathrm{1}},y_{\mathrm{1}}),(x_{\mathrm{2}},y_{\mathrm{2}}),\cdots ,(x_m,y_m)\}\left(x_i,y_i\in \mathbf{R}\right)$这样的训练样本，希望通过训练学习得到一个回归模型$f\left({\mathit{x}}_i\right)={\mathit{w}}^{\mathrm{T}}{\mathit{x}}_i+b$以使得$f\left({\mathit{x}}_i\right)$与${\mathit{y}}_i$尽可能接近，${\mathit{w}}^{\mathrm{T}}$为权值向量，$b$为偏置值。SVR能容忍$f\left({\mathit{x}}_i\right)$与${\mathit{y}}_i$之间最多有$\epsilon$的偏差，故SVR问题可形式化为如下优化问题：

$$\underset{\mathit{w},b}{\mathrm{m}\mathrm{i}\mathrm{n}}\frac{\mathrm{1}}{\mathrm{2}}{‖\mathit{w}‖}^{\mathrm{2}}+C{\displaystyle \sum _{i=\mathrm{1}}^m}{l}_{\epsilon }\left(f\left({\mathit{x}}_i\right)-{\mathit{y}}_i\right)$$

（7）

式中：$C$为正则化常数；$m$为样本数；$l_{\epsilon }$为不敏感损失函数，如式（8）所示。

$$l_{\epsilon }\left(z\right)=\left\{\begin{array}{c}\mathrm{0},\left|z\right|\le \epsilon \\ \left|z\right|-\epsilon ,\left|z\right|>\epsilon \end{array}\right.$$

（8）

引入松弛变量${\mathit{\xi }}_i\mathrm{、}{\mathit{\xi }}_i^{\mathit{*}}$，优化问题转化为

$$\underset{\mathit{w},b,{\mathit{\xi }}_i,\widehat{{\mathit{\xi }}_{\mathrm{i}}}}{\mathrm{m}\mathrm{i}\mathrm{n}}\frac{\mathrm{1}}{\mathrm{2}}{‖\mathit{w}‖}^{\mathrm{2}}+C{\displaystyle \sum _{i=\mathrm{1}}^m}\left({\mathit{\xi }}_i+{\mathit{\xi }}_i^{\mathbf{*}}\right)$$

（9）

$$\mathrm{s}.\mathrm{t}.f\left({\mathit{x}}_i\right)-{\mathit{y}}_i\le \epsilon +{\mathit{\xi }}_i$$$${\mathit{y}}_i-f\left({\mathit{x}}_i\right)\le \epsilon +{\mathit{\xi }}_i^{\mathbf{*}}$$$${\mathit{\xi }}_i\ge \mathrm{0},{\mathit{\xi }}_i^{\mathbf{*}}\ge \mathrm{0},i=\mathrm{1,2},\cdots ,m$$

使用拉格朗日乘子法求得SVR的解如下所示：

$$f\left(\mathit{x}\right)={\displaystyle \sum _{i=\mathrm{1}}^m}\left({\mathit{\alpha }}_i^{\mathbf{*}}-{\mathit{\alpha }}_i\right){{\mathit{x}}_i}^{\mathrm{T}}{\mathit{x}}_i+b$$

（10）

式中：$f\left(\mathit{x}\right)$为打入桩竖向承载力；${\mathit{\alpha }}_i^{\mathbf{*}}$为${\mathit{\alpha }}_i$的伴随矩阵，能使得$\left({\mathit{\alpha }}_i^{\mathit{*}}-{\mathit{\alpha }}_i\right)$不为0的样本即为支持向量。

特征缩放对于梯度下降算法至关重要，对于SVM常采取Z-score标准化方法。经处理的数据符合均值为0，标准差为1的正态分布，公式如下：

$$x^{\mathrm{*}}=\frac{x-\mu }{\sigma }$$

（11）

式中：$x$为标准化前的参数；$\mu$为相应参数的均值；$\sigma$为相应参数的标准差。

标准化的具体做法为：将训练集的特征参数和$Q_{\mathrm{m}}$均标准化后再投入模型训练，然后输入测试集的特征参数到模型以预测测试集的单桩承载力。值得注意的是，对于测试集单桩承载力的反标准化，使用的是训练集中单桩承载力的均值和标准差。

3.2　随机森林模型

随机森林（RF）是采用树模型的梯度提升算法中的代表，它通过组合多个单一决策树，实现了弱评估器的互补，提高了集成算法的整体性能^［16］。

给定$h_{\mathrm{1}}\left(x\right),h_{\mathrm{2}}\left(x\right),\cdots ,h_k(x)$为组成RF的$k$棵决策树，$\mathit{X}$、$\mathit{Y}$为两个随机向量，边缘函数（marginal function）为

$$\begin{array}{l}\mathrm{m}\mathrm{g}\left(\mathit{X},\mathit{Y}\right)=\mathrm{a}\mathrm{v}\mathrm{k}\left(I\left(h_k\left(\mathit{X}\right)=\mathit{Y}\right)\right)-\\ \underset{j\ne Y}{\mathrm{m}\mathrm{a}\mathrm{x}}\mathrm{a}\mathrm{v}\mathrm{k}\left(I\left(h_k\left(\mathit{X}\right)=j\right)\right)\end{array}$$

（12）

式中：$\mathrm{m}\mathrm{g}\left(\mathit{X},\mathit{Y}\right)$为$\mathit{X}$和$\mathit{Y}$的边缘函数；$I()$代表事件发生与0，1二值函数之间的对应关系；分类正确的向量用$\mathit{Y}$表示；分类不正确的向量用$j$表示；$\mathrm{a}\mathrm{v}\mathrm{k}()$代表取平均值。

边缘函数在RF算法中的概念为

$$\mathrm{m}\mathrm{r}\left(\mathit{X},\mathit{Y}\right)=P\left(h_k\left(\mathit{X}\right)=\mathit{Y}\right)-\underset{j\ne Y}{\mathrm{m}\mathrm{a}\mathrm{x}}P\left(h_k\left(\mathit{X}\right)=j\right)$$

（13）

式中：$\mathrm{m}\mathrm{r}\left(\mathit{X},\mathit{Y}\right)$为RF算法中$\mathit{X}$和$\mathit{Y}$的边缘函数；预测正确的分类概率用$P\left(h_k\left(\mathit{X}\right)=y\right)$来表示。

作为一个弱评估器，每棵决策树都与数据集中的一个原始数据集和一个未被抽取的样本组成的数据集互相对应。$Q_k\left(x\right)$表示未被抽取的数据集数目，$Q_k\left(x,y_j\right)$为对于输入的随机向量$x$在$Q_k\left(x\right)$中投票的分类类别为$y_j$的比例，则

$$Q_k\left(x,y_j\right)=\frac{\sum _{k}I\left(h_k\left(x\right)=y_j,\left(x,y\right)\in Q_k\left(x\right)\right)}{\sum _{k}I\left(h_k\left(x\right),\left(x,y\right)\in Q_k\left(x\right)\right)}$$

（14）

RF模型的强度可用边缘函数的期望表征，即

$$s=E\left(P\left(h_k\left(\mathit{X}\right)=\mathit{Y}\right)\right)-\underset{j\ne Y}{\mathrm{m}\mathrm{a}\mathrm{x}}P\left(h_k\left(\mathit{X}\right)=j\right)$$

（15）

式中：$s$为量化的RF算法的强度。

将公式（14）得出的$Q_k\left(x,y_j\right)$代入到公式（15）中，则得到RF集成树模型的强度估计为

$$s=\frac{\mathrm{1}}{n}{\displaystyle \sum _{i=\mathrm{1}}^n}\left(Q\left(x_i,y\right)-\underset{j\ne Y}{\mathrm{m}\mathrm{a}\mathrm{x}}Q\left(x_i,j\right)\right)$$

（16）

对于作为梯度提升算法的RF和XGBoost，则不需要对数据做归一化处理。

3.3　XGBoost模型

XGBoost即极限梯度提升算法，它在普通树模型的基础上改进算法以提高精度，充分调用CPU的多线程并行以提升树达到自身的计算极限^［17］。

基于梯度增强决策树算法，XGBoost通过最小化目标函数至期望范围，预测值计算公式如下：

$${\widehat{y}}_i=\phi \left({\mathit{X}}_{\mathit{i}}\right)={\displaystyle \sum _{k=\mathrm{1}}^K}{f}_k\left({\mathit{X}}_{\mathit{i}}\right),f_k\in F,i\in n$$

（17）

式中：${\widehat{y}}_i$为预测值；${\mathit{X}}_{\mathit{i}}$为输入变量；$i$为均值样本；$f_k$为第$k$个弱评估函数；$n$为样本个数；$F$为所有分类回归树的集合。

目标函数计算公式如下：

$$L\left(\phi \right)={\displaystyle \sum _{i=\mathrm{1}}^n}l\left({\widehat{y}}_i,y_i\right)+{\displaystyle \sum _{k=\mathrm{1}}^K}\Omega \left(f_k\right)$$

（18）

式中：${\displaystyle \sum _{i=\mathrm{1}}^n}l\left({\widehat{y}}_i,y_i\right)$为预测值与实测值的偏差；${\displaystyle \sum _{k=\mathrm{1}}^K}\Omega \left(f_k\right)$为限制算法复杂度的正则化项。

最小化目标函数，通过正则化项平滑算法的学习权重后，最终得到目标函数的最优解如下：

$${\tilde{L}}^{\left(t\right)}\left(q\right)=-\frac{\mathrm{1}}{\mathrm{2}}{\displaystyle \sum _{j=\mathrm{1}}^T}\frac{{\left(\sum _{i\in I_j}{g}_i\right)}^{\mathrm{2}}}{\sum _{i\in I_j}{h}_i+\lambda }+\gamma T$$

（19）

式中：$q$为树结构；$t$为迭代次数；$T$为叶节点的数目；$I_j$为叶子$j$的实例集；$\lambda$和$\gamma$为正则化系数；$g_i$和$h_i$分别为损失函数的一阶导数和二阶导数。

4 预测模型的建立与性能检验

4.1　模型的建立与训练

根据3.3节的特征工程结果，选择桩长、加载模式、桩外径、$c_{\mathrm{p}\mathrm{1}}$、地下水位深度共5个特征作为模型的输入参数，单桩静载试验结果作为模型的输出参数，分别采用SVM、RF、XGBoost 3种机器学习算法建立预测模型。将ZJU‒ICL数据库的115组数据以7∶3的比例划分为训练集和测试集，考虑到数据集的不同划分会对模型的训练和测试造成较大影响，故将随机状态初始化为固定值。

在超参数的设置上，考虑到本文问题为非凸问题，采用贝叶斯优化算法针对训练集进行超参数寻优。选择三项式核函数构建SVM模型，对RF模型锁定random_state为0。调参结果如表3所示。

表3  各模型超参数的贝叶斯寻优结果

Tab. 3  Bayesian optimization results of hyperparameters of each model

机器学习模型	超参数	最优值
SVM	C	4.41
	coef0	3.59
	gamma	0.02
RF	n_estimators	10
	max_depth	15
	min_samples_split	2
	min_samples_leaf	1
	max_features	1
XGBoost	num_round	111
	eta	0.43
	subsample	0.57
	max_depth	7
	colsample_bytree	0.41
	colsample_bylevel	0.82
	colsample_bynode	0.11

对采用最优超参数的各模型进行训练，将模型在训练集上的预测结果与静载试验结果进行对比，并对预测值和实测值进行线性拟合，如图5所示。

图5  各预测模型预测值与实测值比较（训练集）

Fig. 5  Comparison of predicted values and measured values of each prediction model (training set)

由图5可见，3个模型在训练集上都有较好的学习表现。XGBoost模型的学习效果最好，对训练集的拟合优度$R^{\mathrm{2}}$达到了0.999；其次是RF模型和SVM模型，$R^{\mathrm{2}}$分别为0.957和0.935。线性拟合结果也说明XGBoost对训练集的学习最为充分。

4.2　模型的测试与评估

将训练好的模型应用于测试集，预测值与实测值对比结果如图6所示。

图6  各预测模型预测值与实测值比较（测试集）

Fig. 6  Comparison of predicted values and measured values of each prediction model (testing set)

从图6中可见，XGBoost模型在测试集上有着最好的预测效果，决定系数$R^{\mathrm{2}}$达到了0.965；其次是SVM模型和RF模型，分别为0.930和0.906。可见XGBoost模型泛化能力最强，同为梯度提升算法的RF则出现了较严重的过拟合，SVM模型表现介于二者之间。

对各预测模型在数据集上的预测效果使用决定系数$R^{\mathrm{2}}$、均方根误差$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}$、平均绝对误差$\mathrm{M}\mathrm{A}\mathrm{E}$、平均绝对百分比误差$\mathrm{M}\mathrm{A}\mathrm{P}\mathrm{E}$等4种评估指标进行评价，以$R^{\mathrm{2}}$为主要指标降序排序，如表4所示。

表4  各预测模型评估指标

Tab. 4  Evaluation indexes of each prediction model

数据集	预测模型	$R^2$	$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}$/ $\mathrm{k}\mathrm{N}$	$\mathrm{M}\mathrm{A}\mathrm{E}$/ $\mathrm{k}\mathrm{N}$	$\mathrm{M}\mathrm{A}\mathrm{P}\mathrm{E}$/ $\mathrm{\%}$
训练集	XGBoost	0.999	125.9	59.9	4.2
	RF	0.957	685.0	355.0	20.1
	SVM	0.935	840.2	595.2	43.1
测试集	XGBoost	0.965	780.3	541.2	24.0
	SVM	0.930	1 100.5	770.0	30.2
	RF	0.906	1 270.7	779.1	31.2

由表4可以看出，无论在训练集还是测试集上，XGBoost模型相对其他模型均有更优的评估指标。值得注意的是，虽然SVM模型在测试集上的$R^{\mathrm{2}}$和$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}$显著优于RF模型，但$\mathrm{M}\mathrm{A}\mathrm{E}$和$\mathrm{M}\mathrm{A}\mathrm{P}\mathrm{E}$却与RF模型接近。这是因为RF模型对高水平承载力预测效果较差，导致较少的“大偏差值”严重降低了$R^{\mathrm{2}}$和提高了$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}$，而“大偏差值”对$\mathrm{M}\mathrm{A}\mathrm{E}$和$\mathrm{M}\mathrm{A}\mathrm{P}\mathrm{E}$的影响较小。各模型的真实性能需要进一步检验。

4.3　模型的鲁棒性分析

蒙特卡洛模拟本质上是一种使用随机数来解决计算问题的方法，常用于检验模型在一定参数摄动下的鲁棒性^［18］。利用划分数据集的随机数生成 5 000个蒙特卡洛模拟，让各模型连续学习随机划分出的训练集并预测对应的测试集。蒙特卡洛模拟结果的评估指标概率密度分布如图7所示。

图7  蒙特卡洛模拟评估指标概率密度

Fig. 7  Probability density of Monte Carlo simulation evaluation index

由图7可见，XGBoost模型的$R^{\mathrm{2}}$分布最接近1，3个模型$R^{\mathrm{2}}$分布峰值分别在0.90、0.85、0.82左右。XGBoost模型的$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}$分布总体均值最小，其次是RF和SVM模型。$\mathrm{M}\mathrm{A}\mathrm{E}$分布中，XGBoost模型与RF模型峰值对应$\mathrm{M}\mathrm{A}\mathrm{E}$接近，但XGBoost模型峰值更高，SVM模型则表现稍差。结合4.2节中模型的测试结果，可知3个模型均能较准确地预测打入桩竖向承载力。XGBoost模型总体性能最强且鲁棒性最好，RF模型总体性能较好但存在对高水平承载力预测效果较差的问题，SVM模型总体性能相对较差。

利用蒙特卡洛模拟对是否需要地下水位深度这项特征进行检验。以表现最佳的XGBoost模型为例，仅将地下水位深度从输入特征中剔除，其余参数不变，5 000次蒙特卡洛模拟结果的$R^{\mathrm{2}}$概率密度分布如图8所示。

图8  蒙特卡洛模拟${\mathit{R}}^{\mathrm{2}}$概率密度

Fig. 8  R² probability density of Monte Carlo simulation

由图8可见，剔除地下水位深度这项特征后，XGBoost模型的$R^{\mathrm{2}}$分布峰值仍保持在0.90左右，且高水平$R^{\mathrm{2}}$的概率密度变大，可见剔除该特征能一定程度提升模型的鲁棒性和泛用性。最终保留桩长、加载模式、桩外径、CPT降维数据$c_{\mathrm{p}\mathrm{1}}$共4项特征。

4.4　模型的有效性分析

在ZJU‒ICL数据库52个场地中，共有12个场地存在静力触探深度小于该场地中打入桩桩长的情况，经过插值处理后，能运用经验公式的打入桩为87根。为了进一步评价本文模型的性能，选取《建筑桩基设计技术规范》^［2］、《孔压静力触探测试技术规程》^［11］、LCPC法^［12］和机器学习模型对87根桩的竖向承载力进行计算。同时为了检验CPT降维数据$c_{\mathrm{p}\mathrm{1}}$的有效性，采用静探规程中所取的CPT数据均值替换$c_{\mathrm{p}\mathrm{1}}$作为模型输入特征进行对比预测。各模型在测试集上的结果取对数对比如图9所示。

图9  机器学习模型与经验关系模型对比

Fig. 9  Comparison of machine learning models and empirical relation models

由图9可见，经验模型预测结果与真值线偏差较大。这是因为经验模型考虑影响因素不充分，对数据挖掘深度不够，比如在对CPT原始数据的处理上采用取均值法，此法会损失大量信息。而且本文数据来源于全球各地，各经验模型受地区局限性影响较大。计算各模型结果的评估指标，如表5所示。

表5  各预测模型预测值与实测值距离的评估指标

Tab. 5  Evaluation indexes of distance between predicted value and measured value of each prediction model

预测方法	预测模型	$R^2$	$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{E}$/ $\mathrm{k}\mathrm{N}$	$\mathrm{M}\mathrm{A}\mathrm{E}$/ $\mathrm{k}\mathrm{N}$	$\mathrm{M}\mathrm{A}\mathrm{P}\mathrm{E}$/ $\mathrm{\%}$
机器学习模型 （采用$c_{\mathrm{p}1}$）	XGBoost	0.966	883.7	566.7	23.9
	SVM	0.923	1 341.1	926.1	61.7
	RF	0.907	1 470.0	781.9	26.2
机器学习模型 （采用均值）	XGBoost	0.941	1 174.3	717.8	25.7
	RF	0.892	1 587.1	838.3	26.3
	SVM	0.866	1 762.8	1 192.6	81.3
经验关系模型	桩基规范	0.764	2 344.6	1 470.0	65.6
	LCPC	0.749	2 415.7	1 144.1	35.0
	静探规程	0.733	2 493.8	1 511.1	83.7

由表5和图9可见，机器学习模型对打入桩单桩竖向承载力的预测准确性显著高于经验关系模型。输入特征中，采用CPT降维数据$c_{\mathrm{p}\mathrm{1}}$的结果优于采用经验公式所取均值。机器学习模型中，XGBoost模型表现最好，预测值拟合优度$R^{\mathrm{2}}$达0.966。3种经验关系模型表现接近，《建筑桩基设计技术规范》^［2］中经验公式的预测效果相对较好，$R^{\mathrm{2}}$为0.764；其次分别为LCPC法^［12］、《孔压静力触探测试技术规程》^［11］，其$R^{\mathrm{2}}$均在0.75以下。

5 结论

本文采用PCA降维算法处理CPT数据，基于SVM、RF、XGBoost算法构建打入桩竖向承载力的预测模型，对模型的有效性和鲁棒性进行分析。主要结论如下：

（1）对单桩承载力的主要影响因素进行分析，经过特征工程筛选和蒙特卡洛模拟验证，最终确定影响打入桩竖向承载力的4个主要因素：桩长、加载模式、桩外径、桩周土参数（CPT数据表征）。

（2）采用PCA降维算法处理CPT数据以减少信息损失，同时规避CPT深度不满足端阻影响区范围的问题。模型的有效性分析结果表明，此处理方法切实可行。

（3）模型测试及蒙特卡洛模拟结果表明，SVM、RF、XGBoost模型均能较准确地预测打入桩竖向承载力。XGBoost模型总体性能最强且鲁棒性最好，RF模型总体性能较好但存在对高水平承载力预测效果较差的问题，SVM模型总体性能相对较差。

（4）机器学习模型的预测准确性显著高于经验关系模型。经验关系模型中，《建筑桩基设计技术规范》^［2］对ZJU‒ICL数据库中的打入桩竖向承载力有较好的预测效果，可作为工程实践中采用经验公式时的优选。

作者贡献声明

莫品强：数据采集，方案设计，指导论文修改。

林浩东：研究构思，算法实现，稿件撰写。

胡 静：稿件修订，稿件校核。

高 柳：协助完成论文内容。

庄培芝：参与论文的修改。