MATLAB实现基于STFT-RF短时傅里叶变换(STFT)结合随机森林(RF)进行故障诊断分类预测的详细项目实例

请注意此篇内容只是一个项目介绍 更多详细内容可直接联系博主本人 

 或者访问对应标题的完整博客或者文档下载页面(含完整的程序,GUI设计和代码详解)

在现代工业生产与装备运维场景中,旋转机械、电机、轴承、齿轮箱以及各类电力设备的运行状态越来越复杂,运行工况呈现出多工况、高负载、强干扰、强非线性等特点。传统依靠人工经验巡检和简单门限报警的方式,很难在设备早期故障阶段及时识别潜在隐患,往往等到故障扩大、振动剧烈或停机之后才被发现,造成停机损失、维修成本增加,甚至引发安全事故。这就推动了一整套基于信号处理和机器学习的故障诊断方法的快速发展,其中振动信号、电流信号、声发射信号等时序数据被广泛作为状态监测的主要数据来源。

在众多信号分析工具中,短时傅里叶变换是一种兼顾时域和频域信息的经典时频分析方法。面对机械故障信号中频率成分随时间变化的特点,仅仅依靠传统的傅里叶变换只能获得整体频谱,无法刻画出故障特征在时间轴上的演化过程。短时傅里叶变换通过在时间轴上引入滑动窗口,将长时间信号分段分析,从而得到时间-频率二维分布,使信号在不同时间片上的频谱变化得以清晰呈现。在轴承滚动体剥落、齿轮断齿、电机偏心等故障情形中,冲击成分、调制特征常常是间歇性出现或随工况变化而改变,短时傅里叶变换可以在二维时频平面上将其强化,便于提取更具有判别力的特征。

另一方面,机器学习与模式识别技术已经在故障诊断中广泛应用,常见方法包括支持向量机、k近邻、神经网络以及集成学习方法等。随机森林作为集成学习的代表之一,通过构建大量决策树并进行投票,具有良好的泛化性能和抗过拟合能力。相较于单一决策树,随机森林利用样本随机抽取与特征随机选择的机制,在保持模型表达能力的同时显著提升稳定性。此外,随机森林不依赖对数据分布的强假设,能够处理高维特征和非线性边界,非常适合用于复杂机械故障数据的分类任务。

将短时傅里叶变换与随机森林结合,可以形成一条较为完整的“信号处理+特征提取+机器学习分类”的故障诊断流程。短时傅里叶变换负责将原始的时序信号映射为时频表示,通过能量统计、频带能量分布、时频纹理量化等方式得到结构化特征,再将这些特征作为输入喂给随机森林训练分类模型。短时傅里叶变换提供了充分的物理意义和可解释性,能够对应机械系统的固有频率、故障特征频率、倍频和边带;随机森林则为这些特征提供了有效的模式识别能力,使得不同故障类型在特征空间中能够被准确区分。

基于MATLAB R2025b进行实现,具有良好的工程落地价值。MATLAB在信号处理、时频分析和机器学习方面集成了成熟函数库,提供了stft、spectrogram、trainRFClassifier等工具(在本项目中使用fitcensemble构建随机森林),支持对数据进行灵活的预处理、可视化和结果分析。同时MATLAB的脚本化特征和矩阵运算优势,使得从原始振动数据读取、短时傅里叶变换计算、特征构建到模型训练和测试可以在同一环境中顺畅完成。配合R2025b版本中对图形界面函数、可视化对象的更新要求,可以构建兼顾稳定性和可维护性的诊断程序。

此类基于短时傅里叶变换与随机森林的故障诊断项目,既具有明显的工程应用价值,又具备研究意义。一方面可以应用于实际轴承、齿轮、泵、风机、机床主轴、电机等关键设备的状态监测,提升预测性维护水平,减少计划外停机;另一方面也为时频分析与机器学习结合提供了一个可扩展框架,可以进一步与经验模态分解、小波变换、卷积神经网络等方法融合,形成多源多模态智能诊断系统。

在很多公开数据集中,已经可以看到振动信号的采集和标注,例如不同故障类型、不同故障程度、不同载荷工况等标签信息,为构建监督学习模型提供了良好的基础。在此基础上,通过合理设计短时傅里叶变换的窗口大小、重叠率和频率分辨率,以及精心筛选的特征集(如特定频带能量、谱熵、时频矩等),可以显著强化不同故障类别之间的可分性。再通过随机森林进行分类训练,并利用交叉验证和混淆矩阵进行效果评估,可形成一条相对完整、可复现、可扩展的MATLAB实现路径。

此外,短时傅里叶变换的时频图本身具有直观可视化特性,适合作为工程师理解故障机理的辅助工具。通过在MATLAB中绘制时频图,可以观察故障发展过程中频谱结构的变化,验证特征提取的合理性。同时,随机森林提供的特征重要性评估功能也能够帮助分析哪些时频特征在分类中起到关键作用,从而指导传感器布置方案优化、采样参数设置和特征工程策略调整,进一步提高故障诊断系统的整体性能。

综上,基于短时傅里叶变换与随机森林的故障诊断项目,立足于振动信号的时频分析,融合了统计学习的分类能力,兼具可解释性和实际效果,适合在MATLAB R2025b环境中进行全面实现与推广。通过本项目,可以搭建一套从数据采集、信号预处理、特征提取到分类预测与评估的完整流程,为工业设备智能运维提供技术支撑,也为后续深度学习迁移和多算法融合奠定基础。

项目目标与意义

故障可视化与时频特征挖掘

本项目首先关注的目标是实现机械故障信号在时间和频率两个维度上的可视化与特征挖掘。传统时域特征如均方根、峰值因子、峭度等虽然能在一定程度上反映故障信息,但难以全面刻画频率成分随时间变化的动态过程。短时傅里叶变换通过滑动时间窗,对每个时间片进行局部傅里叶变换,生成二维时频幅值矩阵,使得故障冲击、调制、周期性失谐等复杂现象能够在图像中直观体现。通过合理设计窗函数长度、重叠比例和频率分辨率,项目旨在得到具有足够清晰度的时频图,以便后续提取多种统计量和纹理量作为特征。通过这一目标的实现,可以帮助工程技术人员从频域、时频域角度理解故障机理,洞察哪些频带、哪些时间片段对故障分类最为关键,同时也为机器学习算法提供更丰富、更具判别力的输入信息,从而弥补单纯时域分析的不足。

构建鲁棒的随机森林故障分类模型

第二个核心目标是基于短时傅里叶变换提取的特征,构建一个鲁棒、泛化性能强的随机森林故障分类模型。随机森林通过集成多棵决策树,每棵树在样本和特征层面引入随机性,能够有效降低对噪声和异常值的敏感程度,使得模型即使在存在一定测量误差和环境干扰时依然保持较高分类精度。项目将利用MATLAB R2025b中的统计与机器学习工具,针对不同故障类别和不同工况状态进行模型训练和测试,评估随机森林的分类性能。通过调节树的数量、最大深度、分裂特征数等参数,使模型在准确率、训练时间和预测速度之间取得平衡。这一目标的实现,可以为实际工业场景提供一个可直接部署的故障智能识别模块,用于在线监控或者离线诊断,提升维护决策的科学性和及时性。

搭建端到端的MATLAB故障诊断实现流程

第三个目标是利用MATLAB R2025b搭建一条端到端的故障诊断实现流程,从原始数据读取、短时傅里叶变换计算、特征工程、随机森林训练与验证,到结果可视化和性能评估,形成完整的、可复现的工程实现案例。这一流程不仅涵盖信号处理与机器学习部分,还考虑新版本环境中的函数用法规范和图形接口变化,保证脚本兼容性与可维护性。通过将每一环节进行模块化设计,便于后续在不同数据集、不同设备类型、不同故障模式上快速迁移和扩展。此外,流程中还会包括样本划分、交叉验证、混淆矩阵展示等评估机制,帮助使用者全面理解模型表现。这一目标的实现,使故障诊断不再停留在理论推演,而是在具体工程环境中具备直接落地的能力,减少研发投入和试错成本。

提供可扩展的研究与应用基础平台

第四个目标是将基于短时傅里叶变换与随机森林的故障诊断项目打造为一个可扩展的研究与应用基础平台。通过清晰的特征提取接口和模型训练接口,后续可以方便地接入其他信号处理方法,例如小波变换、经验模态分解、变分模态分解等;在分类器层面也可以拓展到支持向量机、梯度提升树、卷积神经网络等更复杂模型。同时,随机森林提供的特征重要性评估结果,也可以用于指导特征降维、传感器布置优化等进一步研究。项目在MATLAB环境中的一体化实现,使研究人员能够快速比较不同算法组合的效果,为学术论文撰写、算法改进验证提供实验平台。在工程应用方面,基于本项目的框架可以逐步演化为在线监控系统中的核心模块,通过与数据采集硬件和上层决策系统对接,形成完整的智能运维解决方案。因此,这一目标的意义在于搭建一个兼具研究和工程推广价值的起点,为后续多方扩展预留足够空间。

项目挑战及解决方案

非平稳振动信号的时频分辨率折中

机械故障振动信号具有显著的非平稳特性,频率成分随时间变化,且冲击、调制和周期性特征往往以短时突变形式出现。短时傅里叶变换虽然能够同时在时间和频率两个维度进行分析,但窗长选择带来分辨率折中:窗长较短时,时间分辨率较高,可以更清晰地定位故障冲击发生的时间位置,但频率分辨率较低,难以区分相邻频率成分;窗长较长时,频率分辨率提高,便于识别特征频率和边带结构,但时间分辨率下降,短时突发特征可能被平均掩盖。这种相互制约关系在实际项目中构成重要挑战。解决方案上,可以从多角度入手。首先,根据典型机械旋转频率、采样频率和故障特征周期,给出一个经验窗长范围,确保至少覆盖若干个转动周期,同时不过度拉长。其次,可以在项目中尝试多组窗长与重叠率组合,对同一数据集分别计算时频图并提取特征,然后采用验证集对比不同参数下随机森林的分类性能,以数据驱动方式选择最优配置。还可以针对不同工况区间,使用不同窗参数进行分段分析,使得在不同转速、不同负载下都能保持较好分辨率。通过这些手段减轻非平稳性对短时傅里叶变换效果的影响,为特征提取打下良好基础。

高维时频特征的冗余与过拟合风险

短时傅里叶变换产生的是一个时间×频率的二维矩阵,直接将所有时频点作为特征输入分类器,会导致特征维度极高,不仅增加计算负担,还容易引入大量冗余信息和噪声,使模型出现过拟合,泛化能力下降。尤其在样本数量有限的条件下,高维特征空间带来的“维度灾难”问题尤为突出。解决这一挑战,需要在特征工程环节进行精细设计。项目中采用特定频带能量统计、主频峰值、边带能量比、谱熵、时频矩等聚合型特征,而不是直接使用原始时频矩阵。通过对每一时频图进行频带划分和统计汇总,将二维大矩阵压缩为几十维到几百维的特征向量,既保留主要故障信息,又避免维度过高。随机森林本身具备一定的特征选择能力,通过随机选取部分特征进行分裂,可以在一定程度上抑制冗余特征影响。此外,还可以利用特征重要性指标,对特征进行筛选和排序,将贡献较小或与噪声相关的特征剔除或降权,进一步降低过拟合风险。结合交叉验证和独立测试集评估模型性能,可以在保证分类精度的前提下获得更紧凑的特征子集。

MATLAB R2025b环境下的实现规范与兼容问题

第三个挑战来自实现环境本身。MATLAB R2025b在图形界面、机器学习函数、可视化组件等方面引入了一些新的规范与限制。例如GUI开发中不再推荐使用uifigure系列组件,项目中需要使用传统figure结合uicontrol进行界面布局;ConfusionMatrixChart对象不再作为Axes的子级,颜色映射需要通过对figure或axes调用colormap,并明确使用turbo等有效颜色图;在机器学习部分,某些回归函数参数调整策略发生变化,虽然本项目主要使用分类模型,但在调试和扩展过程中也需要关注fitrlinear、fitrnet等函数的参数限制。针对这些变化,项目在设计时采用更原生、稳定的接口:在界面部分仅利用figure、subplot、uicontrol来构建简单可视化界面,避免使用在新版本中行为变化较大的高级组件;在可视化混淆矩阵时,采用confusionchart创建独立对象,并通过colormap(gcf,turbo)或colormap(ax,turbo)的方式设置颜色映射,避免将ConfusionMatrixChart对象直接传入colormap。针对随机森林模型训练,则使用fitcensemble结合模板树模板来构建bagged trees,避免依赖在版本间变化较大的接口。通过在项目初期充分了解R2025b的相关限制,并在代码中遵守这些规范,可以显著减少运行错误和兼容问题,为项目在不同环境中的部署提供保障。

项目模型架构

信号获取与预处理模块

整体模型架构的起点是信号获取与预处理模块,该模块负责从实验平台或工业现场采集到的振动、转速、电流等原始信号中提取出适合短时傅里叶变换分析的时间序列。通常使用加速度传感器安装在轴承座、机壳或者支撑结构上,以固定采样频率采集振动信号。在数据进入分析链之前,需要进行多步预处理。首先进行直流分量去除,消除传感器偏置和安装导致的偏移,使信号围绕零均值波动。其次根据采样频率和后续分析目标,进行带通滤波,滤除与目标设备运行无关的低频漂移和高频噪声,保留设备固有频率范围、故障特征频率附近的频段。再者,对数据进行归一化或标准化,将不同工况、不同测点的信号幅值调整到可比尺度,以避免大幅值信号在特征统计中占据过大权重。对于长时间连续监测数据,还需要进行分段,将长序列切分成固定长度的片段,每个片段对应一个样本窗口,以便后续批量进行短时傅里叶变换。预处理模块的设计要兼顾实时性和信息保留,既不能过度滤波导致故障信息丢失,又要降低噪声干扰,保证时频分析的清晰度和随机森林分类的稳定性。

短时傅里叶变换时频分析模块

在信号预处理完成后,进入短时傅里叶变换时频分析模块,这是将一维时间序列映射为二维时频表示的核心环节。短时傅里叶变换在数学上通过滑动窗函数w(t)对信号x(t)进行加窗,再对每个时间位置τ上的窗口内信号执行傅里叶变换,得到X(τ, f),其幅值或功率谱矩阵构成时频图。实现时需要选择合适的窗函数类型(如汉宁窗、汉明窗等)、窗长和相邻窗之间的重叠长度。窗函数类型影响旁瓣泄漏和主瓣宽度,窗长则决定时间与频率分辨率之间的折中。通过stft或spectrogram函数,可以方便地计算出时频矩阵以及对应的时间和频率轴。对于每个故障样本片段,获得一个二维的时频幅值或功率分布图。在项目中还可以对时频矩阵进行对数幅值变换,以压缩动态范围,使得弱小故障特征在图像中更加突出。时频分析模块的输出不仅用于后续特征提取,也可以直接用作可视化,用于观察不同故障类型在时频平面上的典型模式差异,为特征选择和模型设计提供直观参考。

时频特征提取与构造模块

时频特征提取与构造模块负责将高维的时频矩阵转换为结构化、低维度的特征向量,使其适合输入随机森林分类器。在这一模块中,时频矩阵被视为一个二维能量分布场,通过多种统计方式进行压缩与抽象。典型特征包括:特定频带能量总和,将频率轴划分为若干子带,统计每个子带在整个时间范围内的能量积分;能量中心频率与带宽,通过计算加权平均频率和能量分布方差,刻画能量集中位置与扩展程度;谱熵,通过对归一化频谱分布计算熵,衡量频谱复杂度和“平坦程度”;时频矩,通过对时频图进行不同阶矩运算,使其对能量在时间和频率方向的分布形状产生敏感反应。此外,还可以提取最大幅值对应的频率轨迹统计量、多峰分布特征、边带结构能量比等更专业的指标。提取完成后,为每一时频图生成一个特征向量,并与对应的故障类别标签配对。特征提取模块需要在保留主要区分信息的同时尽量控制维度规模,并保证特征对噪声具有一定鲁棒性,以免在实际工况下性能大幅波动。通过合适的特征设计,可以显著提升随机森林在高噪声、高工况变化环境中的分类表现。

随机森林分类与模型训练模块

随机森林分类与模型训练模块是整个架构的机器学习核心部分。该模块接收来自特征提取模块的特征矩阵和对应的类别标签,使用集成学习方法构建故障分类模型。随机森林作为多棵决策树的集合,每棵树在训练时使用自助采样法从原始训练集抽取样本,且在每个节点分裂时只在随机选取的一部分特征中寻找最佳分裂特征,从而引入多重随机性,降低不同树之间的相关性。最终预测阶段,各棵树对输入样本进行独立分类,整体输出采用多数投票结果。项目中使用MATLAB R2025b中的fitcensemble函数,以bagging方式构建基于决策树的随机森林,通过模板树设置最大分裂数、叶节点最小样本数等参数,以控制树的复杂度和训练时间。此外,随机森林能提供特征重要性指标,评估每个特征在整体分类决策中的贡献度,方便后续特征优化。训练过程中采用训练集和验证集划分,以及交叉验证评估模型的泛化性能。通过不断调整树数量、树深度等参数,使得模型在分类精度、稳定性和计算成本之间达到合理平衡。训练完成后,该模块输出可供在线或离线预测使用的故障分类模型。

诊断结果评估与可视化模块

诊断结果评估与可视化模块负责在模型训练与测试后,对分类效果进行量化评价和直观展示。常用的评价指标包括整体分类准确率、各类别的精确率和召回率、F1值等;同时通过绘制混淆矩阵,可以清楚地观察不同故障类别之间的误判情况,例如某些轻微故障是否容易被判为正常或相邻故障等级。MATLAB中可以使用confusionchart绘制混淆矩阵对象,通过设置颜色映射为turbo,提高可读性。此外,还可以对随机森林输出的特征重要性进行可视化,展示贡献度排名靠前的特征,帮助理解模型决策依据。对于时频分析部分,可以选取代表性样本绘制时频图,配合分类结果展示某些关键时频模式与对应故障类型之间的对应关系。在整体评估过程中,还可以通过ROC曲线、错误分类样本的时频图分析等方式,对模型误差来源进行深入排查,指导后续特征工程和模型结构调整。评估模块的存在,使诊断系统不只是输出一个简单的分类结果,而是提供一整套可解释、可分析的诊断证据,为工程师和研究人员提供决策支撑与改进方向。

项目模型描述及代码示例

数据读取与预处理示例
dataFile = 'bearing_data.mat'; % 指定包含振动信号与标签的数据文件名称,用于后续加载样本数据
load(dataFile,'vib_signal','fs','label'); % 从MAT文件中加载振动信号矩阵vib_signal、采样频率fs以及对应的故障类别标签label
vib_signal = vib_signal - mean(vib_signal,1); % 对每个通道信号按列减去均值,去除直流偏置使信号围绕零均值波动
[b_band,a_band] = butter(4,[500 8000]/(fs/2),'bandpass'); % 设计四阶巴特沃斯带通滤波器,通带设定在500到8000Hz以保留主要机械特征频段
vib_filt = filtfilt(b_band,a_band,vib_signal); % 使用零相位前后向滤波对振动信号进行带通处理,避免相位畸变并抑制边缘效应
numSamples = size(vib_filt,1); % 获取滤波后振动信号的总采样点数,用于后续分段处理
segmentLength = round(0.2fs); % 将每个样本窗口长度设为0.2秒,对应的点数为采样频率乘以0.2
overlapLength = round(0.1fs); % 设置相邻样本之间的重叠长度为0.1秒,增加样本数量并提高时频分析平滑性
stepLength = segmentLength - overlapLength; % 计算每次滑动窗口起点之间的步长,用于按固定步长切分整段数据
startIdx = 1:stepLength:(numSamples-segmentLength+1); % 生成每个样本窗口的起始索引集合,保证最后一个窗口完整不越界
numSegments = numel(startIdx); % 计算总的样本窗口数量,用于预分配后续特征矩阵与标签向量大小
segData = zeros(segmentLength,numSegments); % 预先分配样本窗口矩阵,每一列存放一个时间片段,提升内存使用效率
for k = 1:numSegments % 遍历所有窗口索引,逐个提取时间片段构成样本集合
idxRange = startIdx(k):(startIdx(k)+segmentLength-1); % 计算当前窗口对应的样本索引范围,从起始点到起始点加长度减一
segData(:,k) = vib_filt(idxRange,1); % 将第一通道滤波信号在指定范围内的数据拷入样本矩阵的第k列作为当前样本
end
segData = segData./max(abs(segData),[],1); % 对每个样本窗口按列归一化,将最大绝对值缩放为1,减小不同窗口间幅值差异影响
STFT时频计算示例
winLength = round(0.02fs); % 设置短时傅里叶变换的分析窗长度为0.02秒,使时间分辨率与频率分辨率达到折中
stftWin = hamming(winLength,'periodic'); % 构造汉明窗作为分析窗,采用周期性形式以减少频谱泄漏与旁瓣影响
stftOverlap = round(0.5winLength); % 设置相邻短时窗之间的重叠长度为窗长的一半,提高时频图的平滑性和时间采样密度
nfft = 2^nextpow2(winLength); % 将FFT点数设为不小于窗长的二次幂,提升FFT计算效率并获得较细的频率分辨率
sampleIdx = 10; % 选择一个样本索引用于展示短时傅里叶变换结果,这里示例性选取第10个样本窗口
x_seg = segData(:,sampleIdx); % 从样本矩阵中取出对应的时间片段信号,用于时频分析演示
[stftSpec,f_axis,t_axis] = stft(x_seg,fs,'Window',stftWin,'OverlapLength',stftOverlap,'FFTLength',nfft); % 调用短时傅里叶变换函数计算复频谱矩阵及对应的频率和时间轴
stftMag = abs(stftSpec); % 对短时傅里叶变换的复数结果取幅值,将实部和虚部合成为幅度谱用于能量分析
stftPow = stftMag.^2; % 将幅值平方得到功率谱密度估计,用于强调能量分布并便于后续统计特征提取
figure; % 新建一个图形窗口,用于绘制该样本的时频能量分布图
imagesc(t_axis,f_axis,10*log10(stftPow)); % 使用imagesc函数将时频功率谱的对数值作为伪彩色图展示,纵轴为频率横轴为时间
axis xy; % 调整坐标轴方向,使频率轴从下到上增加,时间轴从左到右增加,符合常规时频图习惯
xlabel('Time / s'); % 设置横轴标签为时间单位秒,便于理解短时窗的位置分布
ylabel('Frequency / Hz'); % 设置纵轴标签为频率单位赫兹,标明频率分布范围信息
title('STFT Time-Frequency Representation of One Segment'); % 设置图像标题说明该图展示的是某一样本窗口的短时傅里叶时频图
colormap(gcf,turbo); % 将当前图窗的颜色映射设置为turbo,提高时频图的对比度和视觉识别度
时频特征提取示例
numFreqBands = 8; % 将频率轴划分为8个子带,用于统计每个子带的能量特征以构造频带能量向量
maxFreq = fs/2; % 计算奈奎斯特频率作为频轴上限,用于确定频带边界的划分范围
freqEdges = linspace(0,maxFreq,numFreqBands+1); % 在从0到奈奎斯特频率的区间内均匀划分频带边界,共有numFreqBands个频段
numFeatures = numFreqBands + 3; % 特征总数设置为频带能量数加3个全局统计特征,用于构成综合特征向量示例
features = zeros(numSegments,numFeatures); % 预分配特征矩阵,每一行对应一个样本,每一列对应一个特征分量
for k = 1:numSegments % 遍历每个样本窗口,针对每个时间片段计算时频特征向量
x_k = segData(:,k); % 取出第k个样本的时间信号窗口,作为短时傅里叶变换的输入信号
[S_k,f_k,t_k] = stft(x_k,fs,'Window',stftWin,'OverlapLength',stftOverlap,'FFTLength',nfft); % 对该样本进行短时傅里叶变换获得复谱矩阵以及对应的频率和时间轴
P_k = abs(S_k).^2; % 对复谱矩阵取幅值平方得到功率谱矩阵,作为能量分布基础数据
totalEnergy = sum(P_k(:)); % 计算当前样本时频图的总能量,将整个功率矩阵所有元素求和  
if totalEnergy == 0 % 检查总能量是否为零,以防止后续归一化出现除零情况  
    totalEnergy = 1; % 若总能量为零则强制设为1,保证分母非零同时不改变后续归一化结果的有效性  
end  
P_norm = P_k./totalEnergy; % 将功率谱矩阵按总能量归一化,使整体能量为1便于计算谱熵等规范化特征  
bandEnergy = zeros(1,numFreqBands); % 预分配频带能量数组,用于存放每个频段内的能量总和  
for b = 1:numFreqBands % 遍历每个频带,根据频率索引对能量进行聚合统计  
    f_low = freqEdges(b); % 当前频带的下边界频率,用于构建频率范围条件  
    f_high = freqEdges(b+1); % 当前频带的上边界频率,用于构建频率范围条件  
    idxBand = (f_k>=f_low) & (f_k<f_high); % 根据频率轴选择属于当前频带的频率索引布尔向量  
    bandEnergy(b) = sum(P_k(idxBand,:),'all'); % 对选定频带对应的功率谱元素求和,得到该频带能量总和  
end  
bandEnergyNorm = bandEnergy./sum(bandEnergy); % 将频带能量归一化,使所有频段能量之和为1,消除不同样本总能量差异影响  
f_center = sum(f_k.*sum(P_norm,2)); % 通过频率轴加权求和计算能量中心频率,反映主要能量分布的重心位置  
f_var = sum((f_k - f_center).^2 .* sum(P_norm,2)); % 计算能量分布相对于中心频率的方差,刻画频谱带宽与扩展程度  
p_vec = P_norm(:); % 将归一化功率谱矩阵按列展开成向量,用于计算谱熵等统计量  
p_vec = p_vec + eps; % 为避免对数运算中的零值问题,给每个元素加上机器精度级别的微小正数  
specEntropy = -sum(p_vec.*log(p_vec)); % 计算谱熵作为频谱复杂度度量,能量越分散谱熵值越高,越平坦  
features(k,1:numFreqBands) = bandEnergyNorm; % 将归一化频带能量特征存入特征矩阵前numFreqBands个位置  
features(k,numFreqBands+1) = f_center; % 将能量中心频率写入特征向量的后续分量中,提供频谱位置特征  
features(k,numFreqBands+2) = f_var; % 将能量分布方差作为带宽特征写入特征矩阵,用于反映频谱扩散程度  
features(k,numFreqBands+3) = specEntropy; % 将谱熵作为频谱复杂度特征加入特征矩阵,增加对频谱形态的描述维度  
end
标签构造与数据集划分示例
segLabels = label(1:numSegments); % 从原始标签向量中截取与样本窗口数量一致的部分,作为每个样本的分类标签
uniqueClasses = unique(segLabels); % 获取故障类别标签中的唯一值集合,便于了解总共有多少种不同故障类型
numClasses = numel(uniqueClasses); % 统计故障类别总数,用于后续设置分类模型相关参数或进行结果分析
cv = cvpartition(segLabels,'HoldOut',0.3); % 使用留出法对样本标签进行划分,将30%数据作为测试集其余作为训练集
trainIdx = training(cv); % 从划分结构中获取训练集中样本的逻辑索引,标记哪些样本属于训练集合
testIdx = test(cv); % 从划分结构中获取测试集中样本的逻辑索引,标记哪些样本用于性能评估
X_train = features(trainIdx,:); % 从特征矩阵中提取训练集样本对应的特征向量,作为模型训练数据输入
Y_train = segLabels(trainIdx); % 提取训练集样本的标签向量,作为随机森林训练的目标输出
X_test = features(testIdx,:); % 从特征矩阵中提取测试集样本对应的特征向量,作为模型测试数据输入
Y_test = segLabels(testIdx); % 提取测试集样本的标签向量,用于对训练好的模型进行泛化性能评估
Y_train = categorical(Y_train); % 将训练标签转换为分类类型,符合fitcensemble等分类函数的输入要求
Y_test = categorical(Y_test); % 将测试标签转换为分类类型,便于后续绘制混淆矩阵和计算分类指标
随机森林模型训练示例
t = templateTree('MaxNumSplits',50,'MinLeafSize',5); % 创建决策树模板,限制最多分裂节点数为50,最小叶节点样本数为5以控制树复杂度
numTrees = 100; % 设置随机森林中包含的决策树数量为100棵,在精度和训练时间之间取得平衡
rfModel = fitcensemble(X_train,Y_train, ... % 调用集成学习函数基于训练特征和分类标签构建随机森林模型
'Method','Bag', ... % 指定采用Bagging方法进行集成,通过自助采样构造每棵树降低方差
'NumLearningCycles',numTrees, ... % 指定集成中学习器的数量为预先定义的树数,确保森林规模
'Learners',t); % 使用前面定义的决策树模板作为基学习器结构,统一控制每棵树的深度与叶节点大小
Y_pred = predict(rfModel,X_test); % 使用训练好的随机森林模型对测试集特征进行预测,得到预测的故障类别标签
accuracy = mean(Y_pred == Y_test); % 计算预测标签与真实标签完全匹配的比例,得到测试集上的整体分类准确率
disp(['Test Accuracy = ',num2str(accuracy*100,'%.2f'),' %']); % 在命令窗口输出测试准确率的百分比形式,便于直观查看模型表现
oobPredictions = oobPredict(rfModel); % 调用袋外预测函数,利用未被某棵树看到的样本进行内部误差估计
oobErrorRate = mean(oobPredictions ~= Y_train); % 计算袋外预测标签与真实训练标签不相等的比例得到袋外误差率
disp(['OOB Error Rate = ',num2str(oobErrorRate*100,'%.2f'),' %']); % 输出袋外误差信息,用于评估随机森林泛化能力和是否过拟合
结果评估与可视化示例
figure; % 打开一个新的图形窗口,用于绘制分类结果的混淆矩阵
cm = confusionchart(Y_test,Y_pred); % 使用测试集真实标签与预测标签创建混淆矩阵对象,展示各类之间的分类情况
cm.Title = 'Confusion Matrix of STFT-RF Fault Diagnosis'; % 设置混淆矩阵的标题,说明该图对应短时傅里叶与随机森林诊断结果
cm.RowSummary = 'row-normalized'; % 将行归一化显示,使每一行的预测比例总和为1便于比较不同类别的召回率
cm.ColumnSummary = 'column-normalized'; % 将列归一化显示,使每一列总和为1便于比较不同类别的精确率
colormap(gcf,turbo); % 为当前图窗应用turbo颜色映射,使混淆矩阵各单元格颜色差异更加明显便于分析
figure; % 打开另一图形窗口,用于展示特征重要性排序结果
imp = oobPermutedPredictorImportance(rfModel); % 计算通过袋外数据进行特征置乱后造成误差增加的程度作为特征重要性值
bar(imp); % 使用柱状图将各特征的重要性数值可视化展示,便于比较不同特征贡献大小
xlabel('Feature Index'); % 设置横轴标签为特征索引,指示每根柱子对应的特征编号
ylabel('Importance'); % 设置纵轴标签为重要性数值,表示对模型性能影响的相对程度
title('Feature Importance of STFT-based Features in RF Model'); % 设置图形标题说明柱状图展示的是短时傅里叶特征在随机森林中的重要性排序
简易交互与可视化示例界面
fig = figure('Name','STFT-RF Fault Diagnosis Demo','NumberTitle','off'); % 创建一个图形窗口并设置名称,作为简易演示界面承载图像和控件
fig.Position(3:4) = [900 500]; % 调整图窗尺寸,将宽度设置为900像素高度设置为500像素以便同时显示图像与控件
ax1 = subplot(1,2,1,'Parent',fig); % 在图窗中创建左侧子图轴对象,用于显示原始信号波形或样本窗口信号
ax2 = subplot(1,2,2,'Parent',fig); % 在图窗中创建右侧子图轴对象,用于显示对应样本的短时傅里叶时频图
uicontrol('Parent',fig,'Style','text', ... % 在图窗中添加静态文本控件,用于标示样本编号输入框用途
'String','Segment index:', ... % 设置文本控件显示字符为样本索引提示,指引用户输入感兴趣的样本编号
'Units','normalized', ... % 将控件的位置单位设置为归一化,便于在不同屏幕分辨率下保持相对布局
'Position',[0.1 0.92 0.15 0.05]); % 指定文本控件在图窗中的归一化位置和大小,使其位于上方左侧区域
editIdx = uicontrol('Parent',fig,'Style','edit', ... % 创建可编辑文本框控件,用于接收用户输入的样本索引数值
'String','10', ... % 默认在编辑框中显示数字10,对应初始展示样本编号
'Units','normalized', ... % 使用归一化单位定义位置和尺寸,增强界面在不同显示设备上的自适应性
'Position',[0.25 0.92 0.1 0.05]); % 设置编辑框在图窗顶部中的位置和大小,方便用户操作
btnShow = uicontrol('Parent',fig,'Style','pushbutton', ... % 创建按钮控件,当被点击时触发绘图更新回调函数
'String','Show Segment', ... % 设置按钮显示的文字为“Show Segment”,提示按钮的功能为展示样本
'Units','normalized', ... % 使用归一化单位进行布局以适应不同分辨率
'Position',[0.38 0.92 0.15 0.05]); % 将按钮放置在编辑框右侧,方便用户先输入索引再点击按钮执行
btnShow.Callback = @(src,evt) showSegmentCallback(); % 将按钮的回调函数句柄设置为自定义的内部函数,实现信号与时频图的动态更新
function showSegmentCallback() % 定义回调函数,当用户点击按钮时根据当前编辑框内容更新图像展示
idxStr = get(editIdx,'String'); % 从编辑框控件读取当前输入的字符内容,用于解析成样本索引
idxVal = str2double(idxStr); % 将字符串形式的索引转换为数值,以便在数据矩阵中进行索引操作
if isnan(idxVal) || idxVal<1 || idxVal>numSegments % 判断输入是否为有效数字并在合法样本范围内,防止越界或非法索引
return; % 若输入不合法则直接返回不进行任何操作,保证界面响应安全性
end
idxVal = round(idxVal); % 将索引值四舍五入到最近整数,避免非整数输入造成索引错误
x_show = segData(:,idxVal); % 从样本矩阵中取出对应索引的时间信号片段,用于绘制时域波形和时频图
axes(ax1); % 将当前绘图目标切换到左侧时域轴对象,保证后续绘图在正确子图中显示
cla(ax1); % 清空左侧坐标轴中的已有图形内容,为绘制新波形释放空间
plot((0:segmentLength-1)/fs,x_show,'b'); % 绘制该样本窗口的时域波形,横轴为时间纵轴为归一化振动幅值
xlabel('Time / s'); % 设置左侧图的横轴标签为时间,以秒为单位说明波形时间范围
ylabel('Amplitude'); % 设置左侧图纵轴标签为幅值,表示振动信号的相对强度
title(ax1,['Time waveform of Segment ',num2str(idxVal)]); % 设置左侧图标题指明当前展示的样本窗口索引号
[S_show,f_show,t_show] = stft(x_show,fs,'Window',stftWin,'OverlapLength',stftOverlap,'FFTLength',nfft); % 对当前选中的样本信号执行短时傅里叶变换得到时频谱  
P_show = abs(S_show).^2; % 计算当前样本的功率谱矩阵,用于时频图能量展示  
axes(ax2); % 将绘图目标切换到右侧时频轴对象,确保时频图出现在对应子图中  
cla(ax2); % 清空右侧坐标轴上的旧图像,为新的时频图绘制腾出空间  
imagesc(ax2,t_show,f_show,10*log10(P_show)); % 在右侧轴上绘制当前样本窗口的时频功率谱对数图,以彩色图形式展示能量分布  
axis(ax2,'xy'); % 将右侧轴的坐标方向设置为常规形式,使频率轴自下而上增加时间轴自左而右增加  
xlabel('Time / s'); % 设置右侧图横轴标签,标明时频图的时间尺度  
ylabel('Frequency / Hz'); % 设置右侧图纵轴标签,标明时频图的频率范围  
title(ax2,['STFT of Segment ',num2str(idxVal)]); % 设置右侧图标题,显示当前时频图对应的样本索引号  
colormap(fig,turbo); % 为整个图窗设置turbo颜色映射,使时频能量分布在色彩上更易分辨  
end

更多详细内容请访问

http://故障诊断MATLAB实现基于STFT-RF短时傅里叶变换(STFT)结合随机森林(RF)进行故障诊断分类预测的详细项目实例(含完整的程序,GUI设计和代码详解)_Python实现BO-KNN分类预测资源-CSDN下载 https://download.csdn.net/download/xiaoxingkongyuxi/90241697

 https://download.csdn.net/download/xiaoxingkongyuxi/90241697

http:// https://download.csdn.net/download/xiaoxingkongyuxi/90241697

Logo

AtomGit 是由开放原子开源基金会联合 CSDN 等生态伙伴共同推出的新一代开源与人工智能协作平台。平台坚持“开放、中立、公益”的理念,把代码托管、模型共享、数据集托管、智能体开发体验和算力服务整合在一起,为开发者提供从开发、训练到部署的一站式体验。

更多推荐