CN115935488A

CN115935488A - 一种基于监测数据的高拱坝应力预测模型构建方法

Info

Publication number: CN115935488A
Application number: CN202310024859.5A
Authority: CN
Inventors: 李波; 梁蔚; 余信江; 田亚岭; 郑谦; 杨胜梅; 张启灵; 胡超; 胡蕾; 李志�
Original assignee: Changjiang River Scientific Research Institute Changjiang Water Resources Commission
Current assignee: Changjiang River Scientific Research Institute Changjiang Water Resources Commission
Priority date: 2023-01-09
Filing date: 2023-01-09
Publication date: 2023-04-07

Abstract

本发明提供一种基于监测数据的高拱坝应力预测模型构建方法，包括：基于高拱坝无应力计和应变计组监测数据，建立LightGBM高拱坝应力预测模型；采用Tent混沌映射和折射反向学习策略增强麻雀初始种群多样性和质量，并基于自适应t分布变异使麻雀跳出局部最优位置，提高麻雀搜索算法的全局搜索能力；运用改进的麻雀搜索算法对LightGBM模型进行优化分析，确定最优超参数组合，并引入可解释黑箱模型的SHAP框架，建立一种可解释的高拱坝应力预测模型。本发明融合改进麻雀搜索算法、LightGBM和SHAP的高拱坝应力预测模型能准确预测高拱坝应力，识别出影响高拱坝应力的显著特征，为高拱坝健康监控和诊断提供决策依据。

Description

一种基于监测数据的高拱坝应力预测模型构建方法

技术领域

本发明涉及水工建筑物安全监测领域，具体是一种基于监测数据的高拱坝应力预测模型构建方法。

背景技术

高拱坝坝体在运行期承受较大的荷载，坝体可能出现超标应力。拉应力超标将引起坝体开裂，剪应力超标将导致坝体滑移，因此，坝体应力过大是引起坝体裂缝和整体稳定的主要原因。目前，常通过分析无应力计和应变计组监测数据，建立高拱坝应力预测模型，对高拱坝应力影响因素和高拱坝安全性态进行评估。

在对大坝应力应变监测数据分析并建立预测模型的众多方法中，一般采用主成分分析或偏最小二乘归法等线性分析方法建立统计模型。高拱坝运行期受到上下游水位、坝体温度、两岸山体谷幅收缩变形，以及混凝土徐变度等诸多因素的影响，高拱坝系统是一个复杂的非线性系统，常规的线性回归分析方法，无法建立准确的应力预测模型。一些解决复杂非线性问题的人工智能化方法，如人工神经网络、支持向量机、极限学习机等存在过学习，容易陷入局部最小等问题。因此，如何基于监测数据，建立一种预测能力高、可解释性强的高拱坝应力预测模型是需要亟待解决的问题。

发明内容

针对现有技术中存在的问题，本发明提供一种基于监测数据的高拱坝应力预测模型构建方法，该方法基于无应力计和应变计组监测数据分析，融合Tent混沌映射、折射反向学习策略和自适应t分布变异，提出一种改进麻雀搜索算法(ISSA)进行LightGBM高拱坝应力预测模型参数寻优，并引入SHAP框架增强模型的可解释性，该技术方案将ISSA、LightGBM与SHAP相结合建立的高拱坝应力预测模型具有很高的预测能力和可解释性。

本发明解决其技术问题所采取的技术方案为：

一种基于监测数据的高拱坝应力预测模型构建方法，包括下述步骤：

S1：应变计组应力计算：结合高拱坝无应力计和应变计组监测数据，计算高拱坝单轴应变；根据计算的高拱坝单轴应变，结合混凝土试验得到的弹性模量和徐变度公式，通过变形法计算高拱坝的应力；

S2：LightGBM高拱坝应力预测模型建立：将步骤S1计算所得高拱坝的应力和确定的高拱坝应力影响因素划分训练集和预测集，基于划分的训练集和预测集以及LightGBM算法建立高拱坝应力预测模型；

S3：种群初始化：设置麻雀搜索算法的相关参数和LightGBM超参数的范围，使用Tent混沌映射理论和折射反向学习初始化超参数的麻雀种群；

S4：最优超参数计算：利用步骤S2建立的LightGBM高拱坝应力预测模型对训练集进行训练，计算麻雀个体对预测集的适应度值；更新发现者、加入者和侦察者位置，计算麻雀种群适应度值；利用自适应t分布对所有麻雀位置进行变异，更新麻雀种群位置；判断算法是否达到最大迭代次数，不满足则返回计算麻雀个体的适应度值，满足则输出最优超参数；

S5：可解释的高拱坝应力预测模型建立：根据步骤S4计算的最优超参数，引入SHAP框架，建立可解释的LightGBM高拱坝应力预测模型，所述可解释的LightGBM高拱坝应力预测模型用于评估高拱坝应力影响因素的重要性。

进一步的，步骤S1中应变计组应力计算的具体步骤如下：

S11、结合高拱坝无应力计监测数据，建立无应力计统计模型方程：

ε₀＝a₀+a₁T+a₂t+a₃Ln(1+t)+a₄e^kt (1)

式中：ε₀为无应力计实测值，T为无应力计当前测点温度，t为距分析起始日期的时间长度，a₀,a₁,a₂,a₃,a₄为回归系数，k取-0.01；

S12、通过最小二乘法得到步骤S11中的回归系数，结合高拱坝应变计组监测数据，扣除混凝土自由体积变形得到外部荷载作用下的应力应变ε_n的计算公式如下：

式中：

为应变计组某一方向应变的实测值，

为应变计组实测温度的平均值；

S13、任一方向n的应力应变ε_n与正应变、切应变满足如下关系：

ε_n＝ε_xl²+ε_ym²+ε_zn²+γ_xylm+γ_yzmn+γ_zxnl (3)

式中：ε_x、ε_y、ε_z是x、y、z坐标轴向的正应变，γ_xy、γ_yz、γ_zx是xy、yz、zx坐标面上的切应变，l、m、n是方向余弦，l＝cos(n,x)、m＝cos(n,y)、n＝cos(n,z)；

解如下方程组将实测应变转换为正应变和切应变：

E_n＝Aε (4)

式中：E_n是各向应变计扣除自由体积变形后应变ε_n的向量，A是l、m、n构成的系数矩阵，ε是正应变和切应变的向量；

S14、应变计组单轴应变的计算式如下：

式中：ε^′ _x、ε^′ _y、ε_z ^′是x、y、z坐标轴向的单轴应变，γ_y ^′ _z、γ_x ^′ _z,、γ_x ^′ _y是xy、yz、zx坐标面上的单轴应变；

S15、将单轴应变过程线划分为一系列不等的时段，采用以下加载方式计算各时段应力增量：

式中：τ_i为第i个计算时段的终点龄期；

为第i个计算时段的中点龄期，

ε_i ^′为第i个计算时段终点龄期的单轴应变；

为

时刻混凝土的瞬时弹性模量；

表述以

为加荷龄期持续到τ_i的徐变度度；

表示以

龄期加荷单位应力持续到τ_i的总变形

的导数，称为τ_i时刻的持续弹性模量；Δσ(τ_i)为τ_i时刻的应力增量；

叠加得到任一时刻的应力：

式中：σ(τ_n)为任一时刻τ_n的拱向、径向或垂直向应力。

进一步的，步骤S2中LightGBM高拱坝应力预测模型建立的具体步骤如下：

S21、绘制步骤S1计算的高拱坝应力分别与上游水位、温度和谷幅变形的过程线，定性分析高拱坝应力影响因素的重要性；

S22、通过计算高拱坝应力与上游水位、温度和谷幅变形影响因素间的最大信息系数进行非线性相关性分析，并参考步骤S21的定性分析，将最大信息系数小于0.1的影响因素删掉，确定高拱坝应力y的影响因素x＝(x₁,x₂,…,x_n)；

S23、将高拱坝应力和步骤S22确定的应力影响因素划分为训练集和预测集；

S24、LightGBM算法的目的是找到某函数f(x)的近似值

使特定损失函数L(y,f(x))的期望值最小，如下所示：

最终加入分割后的目标函数为：

式中：g_i和h_i分别为损失函数的一阶与二阶梯度统计量，λ为L₂正则化的系数，I_L和I_R分别为左右分支的样本集，I＝I_L∪I_L为父样本集；

S25、基于步骤S23划分的训练集和预测集，以及步骤S24确定的LightGBM算法，通过垂直生长树，继续深层次的优化，建立高拱坝应力预测模型：

式中：x_i(i＝1,2,…,n)是高拱坝应力的影响因素，y为高拱坝应力。

进一步的，步骤S3中种群初始化的具体步骤如下：

S31、设置麻雀搜索算法的相关参数，包括种群个数M、迭代次数、预警值、发现者比例、侦察者比例，以及LightGBM超参数的范围，包括子叶个数num_leaves、子叶最小数据数min_data_in_leaf、最大深度max_depth和学习率learning_rate范围；

S32、对步骤S31设置的LightGBM超参数，采用Tent映射初始化麻雀种群，得到麻雀混沌种群的表达式为：

式中：x_ij为种群中第i只麻雀在第j维的位置；ub_j和lb_j分别为搜索空间第j维的最小值和最大值；y_j为第j维的Tent映射，其表达式为：

式中：0<α<1，取α＝0.7；y₁为-1到1的随机数；

S33、通过折射反向学习，提升步骤S32得到的麻雀混沌种群的个体质量，麻雀混沌种群折射反向位置的计算公式如下：

式中：

为x_ij的折射反向位置，k为透镜的缩放系数，其他参数见步骤S32中式(11)；

S34、合并步骤S32中的麻雀混沌种群x_ij和步骤S33中的折射反向种群

根据适应度值的升降进行排序，选取适应度值前M个麻雀个体作为LightGBM超参数初始种群。

进一步的，步骤S4中最优超参数计算的具体步骤如下：

S41、利用步骤S3确定的LightGBM超参数麻雀种群和步骤S2基于LightGBM建立的高拱坝应力预测模型对步骤S23划分的训练集进行训练，计算每只麻雀关于步骤S23划分的预测集的适应度值；适应度值的计算公式如下：

式中，y_i为实测值，

为预测值，

为均值，n为实测值个数，RMSE为均方根差；

S42、根据步骤S41计算的麻雀种群中个体适应度值，对所有麻雀个体适应度值进行排序，找到当前最佳和最差适应度值；

S43、更新发现者位置，公式为：

式中：

表示在t+1次迭代第i只麻雀的位置；

为常数，表示最大迭代次数；C_ξ∈(0,1)是一个随机数；R₂∈[0,1]，ST∈[0,1]分别表示预警值和安全值；C_Q为服从正态分布的随机数；L是一个一行多维的矩阵，且元素全为1；

S44、更新加入者位置，公式为：

式中：X_P为当前发现者所在的最佳位置；X_worst表示当前全局最差位置；A是一个各元素为1或-1的一行多维矩阵，A⁺＝A^T(AA^T)；

S45、更新侦察者位置，公式为：

式中：

是当前全局最佳位置；C_β是步长控制参数；C_K∈(0,1)是一个随机数；f_i是当前麻雀的适应度值，f_g和f_w是当前最佳适应度值和最差适应度值；C_ε是一个常数，用于避免分母为零，设置为10e-10；

S46、使用步骤S41中的式(14)计算更新后的麻雀种群适应度值；

S47：对麻雀个体的位置进行t分布变异的过程如下所示：

式中：x_i和

分别为变异前后第i只麻雀的位置，k为迭代次数，

为以迭代次数k作为自由度的t分布；服从自由度为k的t分布下的概率密度函数如下：

S48：使用步骤S41中的式(14)计算步骤S47中所有麻雀变异后的适应度值；与步骤S46确定的麻雀种群适应度值比较，重新排序，更新麻雀种群位置；

S49：判断算法是否达到最大迭代次数，不满足则返回步骤S41，满足则输出最优超参数。

进一步的，步骤S5中可解释的高拱坝应力预测模型建立的具体步骤如下：

S51：根据步骤S4计算的最优超参数，确定最终的LightGBM高拱坝应力预测模型；

S52：对步骤S51确定的LightGBM高拱坝应力预测模型，引入一种利用博弈论Shapley值的新型模型解释方法SHAP，Shapley值的计算如下：

式中：F表示所有特征的集合，S表示从F中去除第i个特征后所有特征子集，f_S∪{i}(x_S∪{i})表示包含第i个特征训练的模型，f_S(x_S)表示在没有第i个特征的情况下训练的模型，φ_i表示第i个特征的Shapley值；

可解释的高拱坝应力预测模型的计算式如下：

式中：g是解释模型；z^′∈{0,1+^M，当该特征存在时为1，否则为0；M是输入特征的个数。

本发明具有如下有益效果：

(1)本发明提出一种新的改进的麻雀搜索算法，该算法采用Tent混沌映射和折射反向学习策略增强麻雀初始种群多样性和质量，并基于自适应t分布变异使麻雀跳出局部最优位置，提高麻雀搜索算法的全局搜索能力；

(2)本发明运用改进的麻雀搜索算法对LightGBM模型进行优化分析，确定最优超参数组合，并引入可解释黑箱模型的SHAP框架，建立一种可解释的高拱坝应力预测模型；该模型可提升训练效率和预测精度，定量分析高拱坝应力影响因素的正负影响力。

附图说明

图1是本发明基于监测数据的高拱坝应力预测模型构建方法其中一个实施例的流程图；

图2是本发明实施例所使用的高拱坝应变计组整体平面布置图；

图3是本发明实施例高拱坝下游应变计组拱向应力空间分布；

图4是本发明实施例高拱坝下游典型应变计组S314-2的拱向应力与上游水位的过程线；

图5是本发明实施例高拱坝下游典型应变计组S314-2的拱向应力与温度测值平均值的过程线；

图6是本发明实施例高拱坝谷幅变形观测线的布置图；

图7是本发明实施例高拱坝下游典型应变计组S314-2的拱向应力与谷幅变形的过程线；

图8是本发明实施例ISSA-LightGBM、SSA-LightGBM和PSO-LightGBM模型的优化迭代过程；

图9是本发明实施例三种模型预测值和实测值的比较图；

图10是本发明实施例三种模型预测的残差图；

图11是本发明实施例ISSA-LightGBM模型的大坝应力各特征因素重要性分布图；

图12是本发明实施例解释模型Shapley值分布图。

具体实施方式

下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整地描述。

请参考图1，本发实施例提供一种基于监测数据的高拱坝应力预测模型构建方法，包括以下步骤：

S1：应变计组应力计算；

为分析、验证所提方法的有效性，选取溪洛渡高拱坝无应力计和应变计组监测数据进行分析；该高拱坝应变计组整体平面分布见图2；

由于该拱坝下游拱向压应力较大，主要对下游应变计组拱向最大压应力建立预测模型；2021年10月21日高拱坝下游应变计组拱向应力空间分布见图3所示；由图3可知，大坝下游拱向应力整体呈受压状态，高压应力区集中在下游河床坝段的中下部；当前拱向最大压应力为10.19MPa，位于14#坝段372m高程应变计组S314-2。

S2：LightGBM高拱坝应力预测模型建立；

根据坝工理论和力学知识，水压对坝体产生应力，由于高拱坝下游为水垫塘，水垫塘水位基本保持不变，因此，水压分量只需要考虑上游水压的影响；高拱坝下游典型应变计组S314-2的拱向应力与上游水位的过程线见图4；由图4可知，高拱坝下游拱向应力与上游水位整体呈负相关，拱向压应力有逐渐增大的趋势；

温度分量为坝体的变温引起的应力；高拱坝下游典型应变计组S314-2的拱向应力与温度测值平均值的过程线见图5；由图5可知，高拱坝下游拱向应力与温度具有明显的相关性，整体呈负相关；

高拱坝两岸山体的谷幅变形对高拱坝挤压引起坝体应力变化；大坝前后布置了7条谷幅测线(VDL01-VDR01，VDL02-VDR02，…，VDL07-VDR07)，其中坝前布置4条谷幅变形观测线，坝后布置3条谷幅观测线，谷幅变形观测线的布置见图6；高拱坝下游典型应变计组S314-2的拱向应力与谷幅变形的过程线见图7；由图7可知，两岸山体的谷幅变形未完全收敛，高拱坝下游拱向应力与谷幅变形整体呈正相关；

高拱坝应力主要与水压、温度、谷幅变形和时效等因素有关，由于谷幅变形也受时效因素的影响，为避免自变量因子的重复，高拱坝应力只需要考虑水压、温度和谷幅变形的影响；该高拱坝2014年3月6日全线浇筑到顶，大坝下游应变计组S314-2拱向应力预测模型选取2014年3月11日～2021年2月1日的监测数据为训练集，共240个样本，2021年2月15日～2021年10月21日的监测数据为预测集，共28个样本；基于LightGBM建立的高拱坝应力预测模型为：

式中：σ是高拱坝拱向应力，H是上游水位，T是应变计组实测温度的平均值，V_j是谷幅变形。

S3：种群初始化；

设置麻雀搜索算法(SSA)的种群个数为20、迭代次数为500、预警值为0.6、发现者比例为70％、侦察者比例为20％；为比较改进麻雀搜索算法(ISSA)的优越性，同时利用粒子群算法(PSO)对步骤S2建立的LightGBM高拱坝应力预测模型进行优化，PSO的学习因子都设置为0.2，加速度常数设置为2，惯性权重从0.9到0.4线性递减；LightGBM高拱坝应力预测模型超参数包括子叶个数num_leaves、子叶最小数据数min_data_in_leaf、最大深度max_depth和学习率learning_rate共4个参数，设置的范围见表1；使用Tent混沌映射理论和折射反向学习初始化超参数的麻雀种群。

表1 LightGBM模型超参数范围

超参数	范围
		num_leaves	[2,30]
min_data_in_leaf	[1,10]
		max_depth	[2,10]
learning_rate	[0.01,0.2]

S4：最优超参数计算；

在步骤S3确定的LightGBM超参数麻雀种群的基础上，运用ISSA、SSA和PSO对步骤S2建立的LightGBM高拱坝应力预测模型进行10次优化分析，计算10次的平均最优适应度值；

ISSA-LightGBM、SSA-LightGBM和PSO-LightGBM模型的优化迭代过程见图8；从图8可以看出，ISSA-LightGBM模型在收敛速度和精度方面没有出现过早收敛的问题；

三种模型预测值和实测值的比较见图9；从图9可以看出，ISSA-LightGBM模型的预测值与实测值基本一致，优于SSA-LightGBM和PSO-LightGBM；

图10为三种模型预测的残差；从图10可以看出，与其他预测模型相比，ISSA-LightGBM的预测残差均值较小，分布集中，表明该模型的预测效果最好；

表2给出了每个模型运行10次时预测集RMSE的最佳值、最差值、均值和标准差，结果表明，ISSA-LightGBM模型的预测精度优于其他两种模型；最终得到的三种LightGBM模型超参数的最优值如表3所示。

表2 不同模型预测集RMSE比较

模型	最佳值	最差值	均值	标准差
					ISSA-LightGBM	0.077	0.077	0.077	0.000
SSA-LightGBM	0.085	0.090	0.087	1.741e-3
					PSO-LightGBM	0.092	0.106	0.093	4.599e-3

表3 不同模型超参数的最优值

S5：可解释的高拱坝应力预测模型建立；

LightGBM可以输出各个特征(大坝应力影响因素)对预测结果的贡献程度，图11为ISSA-LightGBM模型的大坝应力各特征因素重要性分布；从图中可以看出，水压分量H和温度T的影响最显著，谷幅变形对模型影响最小；高拱坝应变计组S314-2位于下游河床坝段的中下部，该部位拱向应力受上游水位影响较小；该测点埋设高程高于水垫塘水位，应力受气温影响显著；谷幅收缩变形对下游坝体有明显的挤压作用，对高拱坝下游低高程部位拱向应力影响较大；综上分析可知，谷幅变形和温度对高拱坝应变计组S314-2拱向应力影响显著，上游水位对拱向应力影响较小，因此，通过LightGBM得到的特征因素重要性与实际不相符。

图12为高拱坝ISSA-LightGBM应力模型分析基础上，计算得到的所有特征的Shapley值分布图；图中每一行代表一个特征，横坐标为Shapley值；一个点代表一个样本，颜色越深说明特征本身数值越大，颜色越浅说明特征本身数值越小；从图中可以看出，谷幅变形V₃、V₁，温度T对模型的影响最显著，水压分量对模型影响最小，并且谷幅变形与大坝应力呈正相关，温度与大坝应力呈负相关，水压分量与大坝下游拱向应力的相关性较小，这些与大坝应力实际变化规律相符。

通过本发明的实施和实验分析，最终可以得到如下结论：

1、本发明提出的融合Tent混沌映射、折射反向学习策略和t分布的改进麻雀搜索算法增加了麻雀种群的搜索多样性和质量，提高了全局搜索能力。

2、本发明提出的提出的ISSA-LightGBM大坝应力预测模型，与SSA-LightGBM和PSO-LightGBM模型相比，具有更高的预测精度和泛化能力。

3、本发明使用SHAP可增强ISSA-LightGBM模型的解释性，评估大坝应力影响因素的重要性，识别出影响大坝应力的显著特征。

综上所述，本发明提出的基于监测数据的高拱坝应力预测模型构建方法是有效的，与传统的应力预测模型相比，具有较高的预测精度和可解释性，因此推荐在实际工程监测中进行推广应用。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何属于本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。