贝叶斯网络故障诊断实战：基于pgmpy的受电弓健康管理

简介：面向西安地铁2号线车辆受电弓无法升弓故障诊断需求，这份资源提供基于贝叶斯网络的完整实现方案与可运行Python代码，适合具备一定编程基础的地铁车辆维护工程师、故障诊断研究人员以及对贝叶斯网络应用感兴趣的开发者。资源为单个docx文档，约42KB，内容既有论文核心思路梳理，也有利用pgmpy库完成网络结构构建、EM算法参数学习、缺失数据处理与诊断推理的逐步讲解，并覆盖故障树向贝叶斯网络转换、条件概率表设定、模拟故障数据生成及后验概率查询等细节。目前已有65人浏览学习。学习者可以结合理论与实践，调整网络结构和参数以适配其他受电弓或车辆子系统故障场景，进而为检修团队快速定位故障原因、优化检修流程提供直接可复用的技术参考。

1. 贝叶斯网络为什么适合地铁受电弓故障诊断

地铁车辆每天跑完一整天，受电弓碳滑板磨到限位、气囊漏气导致弓网接触压力不足、升弓机构卡滞导致拉弧——这些故障在运行数据里往往不是单一传感器超限，而是多个参数同时偏离正常范围。传统的阈值告警会报告一排红灯，维修人员面对几十条报警反而定位不到根因。贝叶斯网络的思路刚好反过来：先把“故障原因”和“可观测信号”之间的因果依赖画成一张有向图，再用条件概率描述每个节点受父节点影响的程度。一旦收到新的传感器证据，网络就能算出“当前最可能的故障是什么”，输出一个带概率的排序结果，而不是孤立的一堆报警。这篇文章从受电弓的故障机理开始，带你一步步完成建模、参数学习、推理和工程化落地，全部核心代码用 Python + pgmpy 实现，可直接复制修改。适合做车辆运维诊断、PHM 方向的一线工程师，也适合刚接触贝叶斯网络的开发者在真实场景里快速上手。

2. 受电弓故障机理与贝叶斯网络建模：从故障树到有向无环图

2.1 先从故障树出发，找到受电弓的可观测证据

受电弓的常见故障模式，在不同线路上略有差异，但基本可以分成三类：气路故障（气囊漏气、减压阀失效）、机械结构故障（铰链卡滞、升弓弹簧断裂）、电气接触故障（碳滑板磨损、弓网接触压力不足导致的拉弧）。做贝叶斯网络之前，我习惯先列出每个故障原因对应的“可观测量”，因为网络里每个叶子节点必须对应一个能拿到数据的传感器或检测项。

常见故障模式与观测变量的对应关系如下表：

这一步的关键在于：不要把所有测点都直接连到故障原因上，中间通常还有“中间状态”。比如“接触压力不足”既是气囊漏气的结果，又是拉弧的原因。把它单独拿出来作为中间节点，能让网络更贴合物理过程，也让后验概率的解释更自然。

2.2 故障树转化为贝叶斯网络时的三类节点

故障树里的“顶事件”“中间事件”“底事件”，对应到贝叶斯网络里就是三类节点：

• 根节点（无父节点）：代表独立的故障原因，比如“气囊漏气”“碳滑板磨损”“机构卡滞”，每个根节点有自己的先验概率。
• 中间节点：代表物理状态的中间变化，比如“接触压力不足”“弓网接触电阻异常”，作用是连接根因和更下游的症状。
• 叶子节点（无子节点）：代表可以直接测到的量，比如“气路压力低”“拉弧次数高”“升弓时间长”。

边（有向箭头）的方向必须遵循因果方向：从“原因”指向“结果”。不能因为“气路压力低时通常能观察到气囊漏气”就把箭头反过来。贝叶斯网络推理虽然允许反向传播证据，但结构上必须是物理因果结构。

以常见做法为例，受电弓诊断网络的邻接关系可以设置成四组因果链：

• 气囊漏气 → 气路压力低 + 接触压力不足
• 接触压力不足 → 拉弧次数高 + 牵引电流波动
• 碳滑板磨损 → 滑板厚度偏小 + 接触电阻偏大
• 升弓机构卡滞 → 升弓时间长 + 接触压力波动

2.3 用 pgmpy 搭建受电弓诊断网络的最小框架

理论部分落到代码上，我用 pgmpy 的BayesianNetwork来定义网络结构。以下代码可以直接复制运行，前提是环境里已经装好 pgmpy：

from pgmpy.models import BayesianNetwork from pgmpy.factors.discrete import TabularCPD # 定义网络结构：每个元组表示一条有向边，方向为 父节点 -> 子节点 model = BayesianNetwork([ ('airbag_leak', 'air_pressure'), # 气囊漏气 -> 气路压力低 ('airbag_leak', 'contact_force'), # 气囊漏气 -> 接触压力不足 ('contact_force', 'arcing'), # 接触压力不足 -> 拉弧次数高 ('slider_wear', 'slider_thickness'), # 碳滑板磨损 -> 滑板厚度偏小 ('slider_wear', 'contact_resistance'), # 碳滑板磨损 -> 接触电阻偏大 ('mechanism_stuck', 'raise_time'), # 机构卡滞 -> 升弓时间长 ('mechanism_stuck', 'contact_force') # 机构卡滞 -> 接触压力波动 ]) # 检查结构是否合法：不允许有环，每个节点必须有定义好的状态 print("网络节点:", model.nodes()) print("网络边:", model.edges())

逻辑说明：BayesianNetwork构造函数的参数是一个包含所有有向边的列表，每条边写成(父节点, 子节点)的元组。节点名建议直接用带下划线的简写，避免在后面的代码里反复写中文变量名。构造完成后调用model.check_model()可以检查图结构是否是合法的有向无环图，如果建边时不小心形成了环，这一步会直接报错。

参数说明：这里引入了 7 个节点，其中airbag_leak、slider_wear、mechanism_stuck是根节点，contact_force是中间节点，air_pressure、arcing、slider_thickness、contact_resistance、raise_time是叶节点。实际项目中节点数往往会更多，比如会增加“环境湿度”“接触网硬点”等外部因素，但结构搭建的方式完全一样。

2.4 条件概率表（CPT）的两种来源：专家经验与历史数据

贝叶斯网络除了结构，还必须给每个节点配备条件概率表。根节点需要先验概率，非根节点需要给定父节点状态下各取值的概率。两种主流做法：

第一种是纯专家经验。比如根据某线路历史统计，气囊漏气的先验概率是 0.02，碳滑板磨损是 0.05。CPT 里的数值可以来自故障树分析中的概率计算，也可以来自检修记录里故障频率的统计。这种方式适合“故障样本少、机理清晰”的场景，新线路没有历史数据时基本只能靠它。

第二种是用历史数据学习。把运行记录和最终的故障定责结果对齐，形成一张“传感器证据 + 故障标签”的数据表，然后调用 pgmpy 的fit()方法从数据里估计条件概率。这种方式适合线路已经运营了足够长时间、积累了大量故障样本的情况。通常我会先用专家经验建一个初始的贝叶斯网络上线跑，跑上三个月后，再用积累的数据重新训练 CPT，形成“先专家后数据”的演进路径。

手动设置 CPT 的代码片段如下：

# 根节点先验概率：airbag_leak 的先验为 P(no)=0.98, P(yes)=0.02 cpd_leak = TabularCPD('airbag_leak', 2, [[0.98, 0.02]]) # 中间节点 contact_force：依赖 airbag_leak 和 mechanism_stuck 两个父节点 # 状态顺序为 [no, yes]，父节点取值顺序从左到右为 # (airbag_leak=no, mechanism_stuck=no), (no, yes), (yes, no), (yes, yes) cpd_force = TabularCPD( 'contact_force', 2, [[0.96, 0.55, 0.30, 0.05], [0.04, 0.45, 0.70, 0.95]], evidence=['airbag_leak', 'mechanism_stuck'], evidence_card=[2, 2] )

逻辑说明：TabularCPD第一个参数是节点名，第二个参数是该节点的状态数量，第三个参数是概率矩阵。矩阵的行对应contact_force的每个取值（第一行是 no，第二行是 yes），矩阵的列对应父节点取值的所有组合，顺序由evidence_card=[2, 2]决定。这里“气囊漏气且机构卡滞同时存在”时接触压力正常的概率只有 0.05，对应矩阵最后一列，能反映多故障并发时的叠加效果。

注意：每一列的概率之和必须等于 1，否则check_model()校验不通过。参数里的 0.04、0.45、0.70、0.95 这些数值直接影响后验推理结果的灵敏度，不建议拍脑袋乱填，至少要跟故障树里的底层事件概率对上。

3. 诊断系统的核心代码实现：参数学习、推理与特征处理

3.1 数据准备：把传感器记录转成贝叶斯网络能吃的离散样本

贝叶斯网络在 pgmpy 里的标准输入是一个 pandas DataFrame，每一行是一条样本，每一列是一个节点变量，取值必须是离散的字符串或整数。但是受电弓的传感器数据大多是连续量，比如气路压力是 460 kPa、升弓时间是 6.8 秒。直接丢进网络之前必须先做离散化。

习惯做法是给每个连续量设定阈值或分箱区间：

import pandas as pd def preprocess_sensor_data(raw_df): df = raw_df.copy() # 气路压力：低于 450 kPa 视为 low，否则 normal df['air_pressure'] = df['air_pressure'].apply( lambda v: 'low' if v < 450 else 'normal' ) # 滑板剩余厚度：小于 22mm 时视为 thin df['slider_thickness'] = df['slider_thickness'].apply( lambda v: 'thin' if v < 22 else 'normal' ) # 升弓时间：大于 9 秒视为 slow df['raise_time'] = df['raise_time'].apply( lambda v: 'slow' if v > 9 else 'normal' ) # 拉弧次数：日累计超过 15 次视为 high df['arcing'] = df['arcing'].apply( lambda v: 'high' if v > 15 else 'normal' ) # 接触电阻和接触压力按相关系数分箱，这里用中位数做二值化 df['contact_resistance'] = pd.qcut( df['contact_resistance'], q=2, labels=['small', 'large'] ) df['contact_force'] = pd.qcut( df['contact_force'], q=2, labels=['low', 'normal'] ) return df

逻辑说明：上限阈值（比如气路压力 450 kPa、升弓时间 9 秒）一般来自受电弓检修规程，可以在部署时配置成 JSON 参数文件，方便不同线路按自己的标准调整。pd.qcut的分箱方式是按样本分布等频切分，这适合没有明确物理阈值的量，有些项目也会改成固定阈值cut，两者各有好处。注意离散化之后必须保证每个节点在训练集和推理集里出现的取值类别完全一致，否则fit()会报 “ValueError: Data must be a DataFrame” 或状态数量不匹配的异常。

3.2 参数学习：有数据时用最大似然估计，数据少时用贝叶斯估计

当历史数据里有明确的故障标签时，就可以跳过手动设置 CPT 的繁琐步骤，直接学到整张概率表。下面这段代码把刚才定义好的空模型model（只包含结构，不包含 CPD）和样本数据train_df放在一起，调用fit()：

from pgmpy.estimators import MaximumLikelihoodEstimator, BayesianEstimator # 方法一：最大似然估计，样本量大时推荐 model_mle = model.copy() model_mle.fit(train_df, estimator=MaximumLikelihoodEstimator) # 方法二：贝叶斯估计，加入 BDeu 先验，避免零概率 model_bdeu = model.copy() model_bdeu.fit( train_df, estimator=BayesianEstimator, prior_type='BDeu', equivalent_sample_size=10 ) # 检查学习到的 CPT print(model_bdeu.get_cpds('contact_force'))

逻辑说明：MaximumLikelihoodEstimator的做法是直接统计每个父节点组合下子节点取值的频率。这种方法在样本充足时无偏，但是当某个父节点组合下没有任何样本时，CPT 里会出现 0 概率，后续推理一旦落到这个分支就会直接得出不可能的结论。BayesianEstimator加了一个equivalent_sample_size=10的 BDeu 先验，相当于给每个父节点组合都预置了 10 条虚拟样本，从数学上完全消除零概率问题。在受电弓故障样本很少的场景下，我一般直接用第二种方法。

参数说明：equivalent_sample_size表示虚拟样本量的强度，设得越大，先验对最终概率的影响越大；设得太小，平滑效果不明显。一般取 5 到 20 之间即可。如果现场故障样本特别稀疏，可以调大这个值，让专家先验继续主导一小段时间。

3.3 故障推理：变元消去算法与证据更新

网络结构和 CPT 都齐了之后，核心工作就是推理：给定一组观测到的传感器证据，计算每个故障原因的后验概率。pgmpy 的VariableElimination提供精确推理，适合几十个节点以内的网络。

from pgmpy.inference import VariableElimination infer = VariableElimination(model_bdeu) # 假设当前列车检测到：气路压力低、拉弧次数高、滑板厚度偏小 evidence = { 'air_pressure': 'low', 'arcing': 'high', 'slider_thickness': 'thin' } # 分别查询三个根节点的后验概率 result_leak = infer.query(variables=['airbag_leak'], evidence=evidence) result_wear = infer.query(variables=['slider_wear'], evidence=evidence) result_stuck = infer.query(variables=['mechanism_stuck'], evidence=evidence) print("气囊漏气概率:", result_leak) print("碳滑板磨损概率:", result_wear) print("机构卡滞概率:", result_stuck)

逻辑说明：query()的第一个参数是想要查询的目标节点列表，第二个参数是证据字典。推理引擎会利用贝叶斯网络的条件独立性把联合分布拆成多个因子相乘，再通过消元逐步边缘化掉无关变量，从而高效算出后验分布。在这个例子里，“拉弧次数高”并不是直接证据对应的根因，但它会通过contact_force这一中间节点间接提高“气囊漏气”和“机构卡滞”的后验概率，这正是贝叶斯网络相对阈值判据的核心优势：能利用间接证据互相印证。

参数说明：evidence字典里的 key 必须和网络里设置的节点名完全一致，value 必须是该节点状态列表里的合法取值。如果跑出KeyError，优先检查预处理函数里的取值是否需要做字符串归一化。如果跑出ValueError: Can not use the given evidence，说明状态类别和 CPT 里的类别对不上，需要查看model_bdeu.get_cpds('节点名').state_names来确认真实的取值集合。

3.4 连续量节点处理的另一种思路：高斯贝叶斯网络

二值化离散化的优点是简单可解释，缺点是对阈值敏感。如果怕阈值选不好导致推理结果不稳定，可以改用LinearGaussianBayesianNetwork，让连续传感器变量直接作为节点，然后让离散故障节点作为父节点、连续变量作为子节点，为每个父节点状态组合拟合一个高斯分布。不过这个路线在工程上比纯离散网络麻烦，因为它的 CPT 变成了均值方差参数表，推理用的是联合树算法，结果也不如离散网络直观。我的建议是：第一阶段先用离散化方案上线，等确认了稳定阈值再考虑连续高斯方案。

4. 实际故障诊断流程与验证：从后验概率到检修决策

4.1 完整的诊断判断逻辑：同时输出多条可能原因

单次推理只给出一个后验概率还不够，现场维修需要的是一个可解释的结论列表。实际系统里我不会只依赖某一个概率阈值，而是先给根节点排序，再结合一个置信度阈值做筛选：

def diagnose(model, evidence, top_k=3, threshold=0.3): infer = VariableElimination(model) root_nodes = ['airbag_leak', 'slider_wear', 'mechanism_stuck'] result = {} for node in root_nodes: q = infer.query(variables=[node], evidence=evidence) # 取出 'yes' 状态对应的后验概率 state_index = q.state_names[node].index('yes') result[node] = q.values[state_index] # 按后验概率降序排列 sorted_result = sorted(result.items(), key=lambda kv: kv[1], reverse=True) # 过滤低于阈值的项 final = [(name, prob) for name, prob in sorted_result if prob >= threshold] return final[:top_k] # 示例：某次运行数据 diagnosis = diagnose( model_bdeu, {'air_pressure': 'low', 'arcing': 'high', 'raise_time': 'slow'} ) for node, prob in diagnosis: print(f"故障原因: {node}, 后验概率: {prob:.2%}")

逻辑说明：这段代码把三个根节点分别查一遍后验概率，然后排序和过滤，形成真正的故障诊断输出。state_index = q.state_names[node].index('yes')这一步很关键，因为 pgmpy 返回的values数组顺序可能和state_names列表一致，但显式地用索引取值更稳妥，避免代码里硬编码第几个下标。处理后，输出结果就可以直接推送到检修终端或生成工单。

参数说明：threshold=0.3是诊断灵敏度的调节旋钮。阈值调低，会报告更多潜在故障，召回率高但容易误报；阈值调高，只报概率最高的一两个故障，精确率高但可能漏掉低概率并发故障。地铁运维场景里，我会区分“预警”和“告警”两档：0.3 以上只推送至日检列表，0.6 以上直接生成故障工单。

4.2 用一个完整案例把全流程串起来

拿一列地铁某日早高峰的传感器数据举例。原始记录值如下：

• 气路压力值：428 kPa → 离散化为low
• 日累计拉弧次数：27 次 → 离散化为high
• 滑板剩余厚度：24.5 mm → 离散化为normal
• 升弓时间：10.2 秒 → 离散化为slow

用上面编号的网络推理，得到的后验概率排序可能为：升弓机构卡滞约 0.63，气囊漏气约 0.28，碳滑板磨损约 0.05。这个结果说明“升弓时间长+接触压力波动”这两个证据强指向机构卡滞，而滑板厚度正常又大幅降低了磨损的后验概率。把这份输出交给检修人员，相当于直接把边界清晰的原因定位从几十条报警里提出来。

4.3 用历史故障样本验证诊断准确率

系统上线前必须回答一个问题：诊断准确率到底多高？验证方法和普通分类模型类似，但要注意贝叶斯网络的“标签”是根节点的真实故障状态。用留一法交叉验证的思路写一个评估循环：

from sklearn.metrics import accuracy_score, precision_recall_fscore_support def evaluate(model, data, labels, root_node): infer = VariableElimination(model) preds = [] for _, row in data.iterrows(): # 把这一行的所有观测变量转换为证据 evidence = { col: row[col] for col in data.columns } q = infer.query(variables=[root_node], evidence=evidence) state_index = q.state_names[root_node].index('yes') preds.append(q.values[state_index] >= 0.5) acc = accuracy_score(labels, preds) precision, recall, f1, _ = precision_recall_fscore_support( labels, preds, average='binary' ) return acc, precision, recall, f1 acc, prec, rec, f1 = evaluate( model_bdeu, test_df[observable_columns], test_df['slider_wear_true'] ) print(f"Accuracy: {acc:.3f}, Precision: {prec:.3f}, Recall: {rec:.3f}, F1: {f1:.3f}")

逻辑说明：对每条测试样本，把除根节点外的所有可观测列都作为证据传入推理引擎，然后看yes状态的后验概率是否超过 0.5。注意这里的data.columns必须只包含模型里的观测变量，不能把slider_wear_true这类标记列也放进去，否则等于把答案提前告诉了推理引擎。

参数说明：0.5 只是默认阈值，验证时可以循环换成 0.3、0.4、0.6 等数值，画出不同阈值下的精确率和召回率曲线。这一步能帮助业务方决定线上系统把预警线放在哪里。如果某条线路的故障率本来就低，看整体 Accuracy 意义不大，更应该关注 Minority 类别的 Recall。

5. 部署经验与调试技巧：样本不平衡、零概率与代码运行环境

5.1 样本不平衡：真实故障太少，模型学不到“异常模式”

地铁受电弓的真实故障率非常低，气囊漏气可能一年也没几次，直接拿全量样本做参数学习，学到的先验概率会极度偏向正常。这会导致一个典型问题：即使气路压力已经低到报警线，由于先验概率正常类太高，后验概率依然不高，诊断被“正常”主导。

解决办法有三个层次。第一，采样时做类平衡，训练贝叶斯网络时把故障样本复制若干倍，或者用 SMOTE 这类算法合成故障样本，但合成样本要谨慎，避免制造出不存在的父节点组合。第二，调大 BDeu 的equivalent_sample_size会影响先验强度，但对样本量极端的类，更有效的是直接手动设置根节点先验，把现场统计到的故障频率写进 CPT。第三，不要在整车上训练一个“通用受电弓模型”，而是按线路分组训练，因为不同线路的接触网状态差异很大，合并训练会把线路特征稀释掉。

5.2 零概率灾难：某个证据组合在训练集里从没出现过

当证据组合复杂到一定程度，训练数据覆盖不到某些父节点取值组合是常见现象。比如airbag_leak=yes且mechanism_stuck=yes的并发故障，在历史样本里一次都没出现，最大似然估计就会把contact_force=normal的概率学到 0。line 推理时一旦证据命中这个分支，概率结果直接失真。

解决办法是前面提到的贝叶斯估计加prior_type='BDeu'。此外还有一个更细的操作：对所有 CPT 里为 0 的单元格做拉普拉斯平滑，在拟合完成后再次检查：

def smooth_zero_cpds(model, alpha=0.01): for cpd in model.get_cpds(): values = cpd.values.copy() values[values == 0] = alpha # 重新归一化，确保每列概率和为 1 normalized = values / values.sum(axis=0) cpd.values = normalized return model

逻辑说明：这段代码在模型推理之前把所有零概率替换成一个极小值，然后按列重新归一化。alpha=0.01适合样本量较大的情况，如果样本量小可以取 0.05 到 0.1。这个技巧在“专家搭建的 CPT + 少量数据修正”混合模式下尤其好用，因为专家也可能漏填某个组合。

5.3 用 BIF 文件持久化模型，方便跨环境复用

模型训练好之后，不能每次都从原始数据重新训练。pgmpy 支持把整个贝叶斯网络（包括结构、CPT、节点状态）写成标准 BIF 格式，这样在列车诊断主机、地面服务器、云平台之间迁移时就只要拷贝文件：

from pgmpy.readwrite import BIFReader, BIFWriter # 保存 BIFWriter(model_bdeu).write_bif('pantograph_diagnosis.bif') # 加载 model = BIFReader('pantograph_diagnosis.bif').get_model()

逻辑说明：BIF 文件是纯文本的 XML 格式，可以通过版本管理工具统一追踪每条线路的网络版本。现场升级模型时，只需要把新的 BIF 文件下发到列车的边缘计算单元，重启推理服务即可生效。注意保存前确认所有 CPT 都已完成平滑，否则加载后再去改会比较绕。

5.4 关于代码运行环境的一个实操提示

这套代码不挑运行环境，最常见的执行方式有三种：在 WSL Ubuntu 的终端里直接跑 Python 脚本，在 VS Code 里新建 Jupyter Notebook 一步步调试，或者把诊断代码封装成 FastAPI 服务部署在服务器上。我建议在本地调试阶段直接用python diagnosis.py运行，把上面提到的函数按顺序组织在同一个文件里，入口处加上if __name__ == '__main__':调用一次完整推理，方便验证数据流是否畅通。把print("模型加载完成，等待传感器数据...")放在入口位置，测试时能清晰看到脚本执行到哪一步停了。

最后想留一个实用的验证技巧：当你在现场拿到一批新的运行数据，不要急着把整个模型重新训练一遍，而是先用旧模型对单条异常样本做推理，把根节点的后验概率分布画成柱状图，对比新增样本前后排序是否发生变化。如果某条证据加入后概率排序基本不变，说明这条证据对区分故障原因几乎没有贡献，后续可以把它从观测变量里逐步剔除，简化网络的维护成本。

本文还有配套的精品资源，点击获取