ARTICLE DETAIL

资讯详情

深耕编程入门与网站建设的一线实战洞察。

基于随机森林的桥梁易损性分析:PGA-位移预测与曲线绘制

基于随机森林的桥梁易损性分析:PGA-位移预测与曲线绘制 简介桥梁地震易损性分析的Python随机森林实现面向地震工程、桥梁安全评估与机器学习交叉领域的工程师、研究人员及学生旨在通过数据驱动方式评估不同地震强度下桥梁损伤风险并可绘制易损性曲线。包内以1个docx文档封装体积仅15KB包含可直接运行的Python代码和配套使用解释。内容覆盖数据加载与预处理、训练集与测试集划分、随机森林回归模型构建与预测以及基于位移阈值计算超过概率等关键流程同时针对实际工程中多特征建模与概率精度提升给出提醒。从导入库到曲线绘制均有步骤说明适合具备一定Python编程和机器学习基础的读者快速上手。目前已有75人浏览学习可帮助缩短代码调试时间快速完成从数据到易损性曲线的实验流程为桥梁抗震性能研究提供参考。1. 为什么是随机森林桥梁易损性分析的非参数回归选择拿到一批不同 PGA 下的桥墩位移响应数据第一反应是线性回归直接拟合一条直线但地震作用下的桥梁响应几乎不满足线性假设。桥墩屈服后位移增速明显变快结构进入塑性阶段PGA 与位移之间普遍存在拐点和平台段此时随机森林作为非参数回归器的优势就体现出来了——它不要求预设函数形式树模型通过特征空间切分得到分段常数近似能天然追踪这类非线性变化。另一个现实因素是工程数据通常样本量不大几十到几百条记录很常见随机森林对噪声和缺失值没有神经网络敏感默认参数跑出来的结果一般不至于离谱。这篇内容面向有 Python 和 sklearn 基础的工程师与研究人员侧重数据组织、模型参数和曲线绘制三个环节的具体做法。2. 数据准备与特征工程PGA之外的桥梁响应建模2.1 从原始记录到建模表字段怎么留在项目里最花时间的往往不是模型而是数据表整理。原始地震响应记录里通常包含多条测点、多个工况直接拿来做模型会混入不该有的信息。常见做法是维护一张宽表每一行代表一次地震动输入下的桥梁响应样本至少包含字段示例建模型时是否保留地震动记录编号EQ-014可作索引不进特征PGA (g)0.35保留为特征墩顶位移 (m)0.42作为目标变量损伤状态标签DS1/DS2/DS3可作为分类目标或辅助验证墩高 (m)12.5多特征建模时才保留场地类别II多特征建模时才保留如果手头数据只有 PGA 和位移两列也不影响流程启动只是模型可解释的空间会小一些。读取后用data.info()和data.describe()看一眼字段完整度和量纲这是地震工程数据的习惯动作因为位移的单位可能是毫米、厘米或米量纲不统一会让阈值设定环节直接崩掉。2.2 数据检查与缺失值处理shuffle 前必须做的三件事import pandas as pd import numpy as np data pd.read_csv(bridge_data.csv) print(data.head()) print(data.isnull().sum()) # 工程数据里常见问题是 PGA 列出现 0 或负值先过滤异常记录 data data[data[PGA] 0] data data.dropna(subset[PGA, Displacement]) # 统一单位这里假设位移列以米为单位如果数据是 mm 则除以 1000 # data[Displacement] data[Displacement] / 1000.0 X data[PGA].values.reshape(-1, 1) y data[Displacement].valuesreshape(-1, 1)是把一维数组变成列向量sklearn 的特征矩阵要求二维结构行数自动推断列数固定为 1。这一步初看只是在迁就接口实际上是在明确一个语义当前模型只接受一个输入特征 PGA。如果后续要加入墩高和场地类别就要从 DataFrame 里一次性选出多列构造 X而不是继续用单列的写法。对于缺失值随机森林本身不直接支持 NaN所以要么删除要么用均值或中位数填充。在桥梁响应数据里位移列为空的记录往往对应传感器失效或采样异常直接删除通常是比填充更稳妥的选择因为这些样本的标签本身并不可信。data.dropna(subset[PGA, Displacement])只针对建模必需的两列做删除其他字段的缺失可以后面再处理避免过早丢失整行数据。2.3 单特征与多特征什么时候不该只用 PGAPGA 是最常用的地震动强度指标但不是唯一选项。Sa(T1)结构基本周期对应的谱加速度在很多易损性研究里比 PGA 更有效因为它考虑了结构自振频率与地震动频谱的耦合。如果原始数据里只有 PGA模型给出的位移预测就已经隐含了“数据里只有强震幅值信息”这个限制曲线形状会偏保守或偏乐观取决于桥梁周期落在哪个频段。当数据表里存在墩高、场地类别、支座类型等字段时把它们一起放进特征矩阵会显著提升位移预测的精度。比如两个相同 PGA 工况一个在 II 类场地、一个在 IV 类场地墩顶位移可能差两倍仅靠 PGA 无法区分这种差异。随机森林对类别特征的处理方式是直接做数值编码用pd.get_dummies()或OrdinalEncoder都行树模型的切分不受特征缩放影响这是它比 SVM 或 KNN 更省事的地方。# 多特征建模示例PGA 为连续特征场地类别做独热编码 feature_cols [PGA, Pier_Height] X_multi data[feature_cols].copy() X_multi pd.concat( [X_multi, pd.get_dummies(data[Site_Class], prefixSite)], axis1 ) print(X_multi.head())pd.get_dummies会把 II、III、IV 这类场地类别展开成多列 0/1 标志位比如Site_II、Site_III、Site_IV。这样做的好处是树模型可以在切分时分别利用每个类别的信息缺点是类别过多时特征矩阵膨胀在本就小的工程数据集上容易引入噪声。类别少于 5 种时用独热编码没问题超过 5 种可以考虑用OrdinalEncoder做序数编码来压低维度。建模前建议做一次相关性检查data.corr()输出里如果两个特征相关系数超过 0.8保留其中一个即可避免冗余特征稀释树分裂的随机性。墩高和截面惯性矩常常高度相关两者同时入模并不会带来信息增量却会让特征重要性排在后面的那一个被错误解读成“不重要”。3. 随机森林回归器构建与参数调优从默认值到可解释结果3.1 训练测试划分随机种子与样本量约束from sklearn.model_selection import train_test_split X_train, X_test, y_train, y_test train_test_split( X, y, test_size0.2, random_state42 )test_size0.2意味着如果样本总量只有 50 条测试集只有 10 条评价指标波动会很大。工程数据常见做法是先看样本量再定划分比例样本量低于 100 条时test_size建议设为 0.15 到 0.2 之间random_state固定下来方便复现。random_state42本身没有魔法它只是保证每次运行都切出同一份训练集和测试集这在项目评审里很重要——不同次运行结果不一致别人很难复核。如果数据表里带有桥梁编号划分时要考虑是否按桥梁分组。train_test_split默认是逐行随机划分同一个桥梁的记录可能同时出现在训练集和测试集里这会让模型在测试集上表现虚高因为分布特征已经被训练阶段见过。更严谨的做法是用GroupShuffleSplit按桥梁编号分组代价是测试集样本量和代表性都会打折扣具体取舍要看项目目标——发论文建议分组做传感器布设方案可以接受逐行划分。3.2 RandomForestRegressor 参数语义与工程取值rf_model RandomForestRegressor( n_estimators200, max_depthNone, min_samples_split3, min_samples_leaf1, oob_scoreTrue, random_state42 ) rf_model.fit(X_train, y_train)n_estimators是树的数量。100 到 300 之间是较常见的工程区间树太少会增大单树影响、预测方差大超过 300 收益快速衰减计算时间却线性增长。max_depthNone表示让树自由生长到叶子纯化在样本量小时不容易出现过拟合因为每棵树只看到随机抽样的子集和特征子集。min_samples_split3是说内部节点至少要有 3 个样本才继续分裂这个值对噪声多的位移数据挺有用默认 2 会让树对个别异常点过度敏感。min_samples_leaf1保持默认即可位移目标是连续值没必要强制叶子包含多个样本。oob_scoreTrue会输出袋外数据上的 R²用于快速检验模型没见过的数据表现。参数默认值工程建议区间说明n_estimators100100-300树数量越大越稳定max_depthNoneNone 或 10-20限制深度可抑制过拟合min_samples_split22-5样本少时调大min_samples_leaf11-3连续目标保持 1 即可oob_scoreFalseTrue快速看袋外精度上面的参数并不是最优解而是在 200 条左右的桥梁响应数据上先跑通流程的保守配置。再往下调参通常用GridSearchCV做小范围搜索对这类数据我一般只搜索max_depth和min_samples_split两个参数n_estimators不参与搜索直接取一个稍大的固定值原因是树数量对最终精度的影响远小于分裂策略搜索它能得到的收益很有限。3.3 训练后的第一件事OOB 分数与残差检查y_pred rf_model.predict(X_test) print(OOB R^2:, rf_model.oob_score_) print(Test R^2:, rf_model.score(X_test, y_test)) residuals y_test - y_pred print(残差均值:, np.mean(residuals)) print(残差标准差:, np.std(residuals))OOB 分数和测试集 R² 的差异值得留意。两者差距超过 0.15 通常说明训练集和测试集分布不一致典型原因是测试集里高 PGA 样本比例过少模型在未见过的强度区间只能外推。残差均值应该接近 0如果系统性偏正或偏负可能是位移的单位在中途被修改过或者样本里混入了不同量测位置的数据。R² 在易损性分析里的参考意义有限更重要的评估是看预测位移在损伤阈值附近的分类准确率这个放到后面一章展开。4. 易损性曲线绘制从位移预测到超阈值概率4.1 散点与曲线直接 plot 预测值为什么是错的示例代码里有一行plt.plot(X_test, y_pred, colorred)这在测试集顺序没有按 PGA 排序时会给出一条来回穿越的折线看起来像噪声而不是曲线。正确做法有两种一是在绘图前对X_test和y_pred一起按 PGA 排序二是生成一段连续 PGA 网格用模型逐个预测再连线。第二种做法更通用也方便后续叠加多个阈值的概率曲线。import matplotlib.pyplot as plt # 先生成连续的 PGA 网格覆盖训练数据的取值范围 pga_grid np.linspace( X[PGA].min(), X[PGA].max(), 200 ).reshape(-1, 1) y_grid_pred rf_model.predict(pga_grid) # 以训练集为背景散点以网格预测为红色实线 plt.figure(figsize(8, 5)) plt.scatter(X[PGA], y, alpha0.4, labelObserved) plt.plot(pga_grid, y_grid_pred, colorred, linewidth2, labelRF prediction) plt.xlabel(PGA (g)) plt.ylabel(Displacement (m)) plt.title(Displacement vs PGA with Random Forest) plt.legend() plt.grid(alpha0.3) plt.show()np.linspace(start, stop, 200)在 PGA 最小值到最大值之间均匀生成 200 个点reshape(-1, 1)再降维成 sklearn 需要的二维特征矩阵。网格密度取 200 个点足够画出一条平滑的曲线再多对视觉没有增益。散点透明度alpha0.4是为了让重叠点不至于糊成一片这点在数据量大时尤其明显。如果希望位移预测曲线在尾部不出现下弯要检查 PGA 最大值附近是否有异常样本。随机森林是分段常数模型它在训练数据范围之外的预测是最近叶子的均值天然不能外插所以曲线两端通常出现平台段这是算法特性不是代码 bug。4.2 超过某一位移阈值的概率硬分类近似不满足概率光滑性示例代码里的易损性概率计算用(y_pred_threshold threshold).astype(int)把每个 PGA 网格点的预测位移硬编码成 0 或 1。这在语义上是一个确定性判断画出来是一条阶跃曲线会在阈值对应位置突变工程上不太好拿来做决策。更稳妥的做法是利用随机森林内置的预测分布或者退一步做分类问题。from sklearn.ensemble import RandomForestClassifier threshold 0.5 clf RandomForestClassifier( n_estimators200, max_depth5, min_samples_leaf5, class_weightbalanced, random_state42 ) clf.fit(X_train, (y_train threshold).astype(int)) pga_grid np.linspace(X[PGA].min(), X[PGA].max(), 200).reshape(-1, 1) prob_exceed clf.predict_proba(pga_grid)[:, 1] plt.plot(pga_grid, prob_exceed, labelfP(D {threshold} m)) plt.xlabel(PGA (g)) plt.ylabel(Probability of exceedance) plt.ylim(0, 1.05) plt.legend() plt.grid(alpha0.3) plt.show()概率分类器直接输出连续的 0 到 1 之间的概率估计比硬阈值判断平滑得多。class_weightbalanced的作用是让少数类样本获得更高权重位移超过阈值的事件在低 PGA 区间几乎不会发生正负样本天然不平衡不加这个参数分类器会倾向于把一切都预测为“未超限”。max_depth5限制树的复杂度因为分类目标的噪声比回归目标更大压缩深度能减少训练集上的过拟合。用分类概率代替回归预测与阈值的比较本质上是改变了问题定义回归模型的超阈值概率是从连续预测中推导出来的分类模型直接逼近条件概率 P(位移阈值 | PGA)。这两者在高度非线性的响应模式上会给出不同曲线若工程报告指定了某一种方法不要随意切换。当样本量大到每个 PGA bin 都有足够数量时回归加阈值判断和分类概率应该收敛到相近结果如果差距很大优先怀疑训练集划分的随机性。4.3 多阈值绘图一条易损性曲线不够一组曲线才有工程意义thresholds [0.2, 0.4, 0.6, 0.8] fig, ax plt.subplots(figsize(9, 5)) for thr in thresholds: clf RandomForestClassifier( n_estimators200, max_depth5, class_weightbalanced, random_state42 ) clf.fit(X_train, (y_train thr).astype(int)) prob clf.predict_proba(pga_grid)[:, 1] ax.plot(pga_grid, prob, labelfD {thr} m) ax.set_xlabel(PGA (g)) ax.set_ylabel(Probability of exceedance) ax.set_title(Bridge seismic fragility curves by threshold) ax.legend() ax.grid(alpha0.3) plt.tight_layout() plt.show()四条曲线分别对应不同损伤状态阈值。0.2m 对应轻微损伤、0.8m 接近严重损伤这样的划分需要结合桥梁设计位移角限值标定。绘制时要注意不同阈值曲线在低 PGA 端不应交叉如果出现交叉大概率是某个阈值档位正负样本太少分类器在低 PGA 段经历了明显过拟合需要考虑降低该阈值档位的树深度或增加样本量。5. 验证与扩展从单条曲线到可用的工程结论5.1 阈值敏感性一组阈值比一条曲线信息量大易损性曲线的工程价值在于回答“PGA 达到多少时某损伤状态开始出现”。单条曲线能给出中位值50% 概率对应的 PGA但决策者更需要知道不同阈值下的曲线簇是否分层清晰。把 0.2、0.4、0.6、0.8 四个阈值画在同一张图里能直接读出中等损伤与严重损伤之间的 PGA 间隔如果间隔过小说明结构在某一强度区间内迅速从轻微损伤滑向严重损伤这本身是设计上需要注意的信号。做阈值敏感性分析时可以每次都改变阈值并重训练分类器观察中位 PGA 是否随阈值单调递增。用scipy.interpolate.interp1d从曲线上反解 0.5 概率对应的 PGA 值比直接用坐标数组找离 0.5 最近的索引更准确避免网格间距导致读数偏差。若出现中位 PGA 下降的非单调情况常见原因是分类器在某个阈值下遇到了类别极端不平衡此时检查该档位正样本数量低于 15 个时建议采用 SMOTE 或直接下调min_samples_leaf来增加决策边界稳定性。5.2 与有限元易损性结果的对比检验随机森林给出的易损性曲线本质上是数据驱动拟合缺乏物理机制约束。拿它与 OpenSees 或 Abaqus 的非线性时程分析结果做叠加对比是验证流程的标准动作如果两条曲线在中位 PGA 处的偏差超过 20%先检查输入数据的 PGA 分布是否覆盖了结构屈服的敏感区间再检查训练集的损伤状态标签划分是否符合有限元模型定义的损伤等级。数据驱动的曲线不需要和有限元完全重合但趋势和量级应一致否则模型训练好也不能进入工程报告。对比检验时还有个容易忽略的细节有限元分析通常给出的是能力谱曲线纵轴往往是位移角或延性系数而随机森林模型直接用的是墩顶绝对位移。两者之间差一个墩高的换算系数对比前先把单位统一到同一物理量上。曲线簇绘制完成后可以输出一张汇总表列出每个损伤阈值对应的中位 PGA 与 16%/84% 分位数这个格式在审查时最好用能直接进入报告正文的表格区。本文还有配套的精品资源点击获取
返回列表