特征提取（熵特征指标）

最后更新于：2024-04-07 21:27:03

一、代码运行环境

MATLAB2018a及更新版本。

二、程序介绍

注：图标代表该m文件为脚本文件，可以直接运行；图标代表函数文件，在没有输入变量的情况下无法直接运行。更详细的解释可以看这里。

1.testGenFeaEn.m文件

特征提取函数genFeatureEn的测试文件，可以直接运行。程序运行完成后，在MATLAB的工作区，双击fea变量，可以查看求得的具体数值。

2.genFeatureEn.m文件

熵特征信号特征提取函数。7种全都集中到一个封装函数里，只需要输入待提取数据、采样频率和特征类型，实现一行代码完成特征提取。

function fea = genFeatureEn(data,featureNamesCell,options)

% 熵相关算法的信号特征提取函数
% 输入：
% data：待特征提取的时域信号，可以是二维数据，维度为m*n，其中m为数据组数，n为每组数据的长度。即每行数据为一组。行列方向不可出错
% options：其他设置，使用结构体的方式导入。目前可设置变量包括：
%   -svdpEn：即奇异值的窗口长度。
%   -Apdim：近似熵参数，Apdim为近似熵的模式维度
%   -Apr：近似熵参数，Apr为近似熵的阈值
%   -Spdim：样本熵参数，Spdim为样本熵的模式维度
%   -Spr：Spr为样本熵的阈值
%   -Fuzdim：模糊熵参数，Fuzdim为模糊熵模式维度
%   -Fuzr：模糊熵参数，Fuzr为模糊熵的阈值
%   -Fuzn：模糊熵参数，Fuzn为模糊熵权重
%   -Pedim：排列熵参数，Pedim为排列熵模式维度
%   -Pet：排列熵参数，Pet为排列熵的时间延迟
%   -fs：采样频率，采样频率即每秒钟采集的数据点数，按照实际情况设置，该参数目前在包络熵特征采集中用到
% featureNamesCell：拟进行特征提取的特征名称，该变量为cell类型，其中包含的特征名称为字符串，特征名称需要在下边列表中：
% 目前支持的特征（2022.7.10，共8种）：
%      psdE：功率谱熵
%      svdpE：奇异谱熵
%      eE：能量熵
%      ApEn：近似熵
%      SampleEn：样本熵
%      FuzzyEn：模糊熵
%      PerEn：排列熵
%      enveEn：包络熵
% 
% 输出：
% fea：数据data的特征值数组，其特征值顺序与featureNamesCell一一对应

该文件中还包含了功率谱熵、奇异谱熵、能量熵、近似熵、样本熵、排列熵、模糊熵、包络熵的子函数，函数介绍如下：

function ie = kInformationEntopy(sig,SegmentNum)

% 计算信号的信息熵
% 参考《矿用带式输送机托辊远程故障诊断系统》
% 输入：
% sig：输入信号
% SegmentNum：拟分组数，如果不输入，则自动使用斯特格斯（Sturges）经验公式计算。
% 输出：
% ie：信息熵求解结果

function svdpE = kSingularSpectrumEntropy(data,n)

% 求信号的奇异谱熵
% 参考《矿用带式输送机托辊远程故障诊断系统》
% 输入：
% data：待分析信号
% n：窗口长度，注意：2<=n<=lenght(data)-1
% 输出：
% svdpE：奇异谱熵值

function psdE = kPowerSpectrumEntropy(data)

% 求信号的功率谱熵
% 功率谱使用周期图法计算得到
% 参考《矿用带式输送机托辊远程故障诊断系统》
% 输入：
% data：待分析信号
% 输出：
% psdE：功率谱熵值

function eE = kEnergyEntropy(data)

% 求信号基于emd分解算法的能量熵
% 参考《面向高铁走行部故障诊断算法的研究与实现》
% 输入：
% data：待分析信号
% 输出：
% eE：能量熵值

function ApEn = kApproximateEntropy(data, dim, r)


%  计算近似熵(ApEn)，在网上流传的ApproximateEntropy函数基础上修改
%  输入：
%  data  - 待分析数据，需要是一维数据
%  dim - 模式维度
%  r  -  阈值大小，一般选择r=0.1~0.25，（网上流传的ApproximateEntropy需要再乘以data数据的标准差，本程序不需要，只用输入0.1-0.25间的数字即可）
%  输出：
%  ApEn：近似熵值大小
%  (网上流传的ApproximateEntropy还有一个降采样参数tau，此函数中不进行降采样，删去该参数)
%  论文：Pincus S M . Approximate entropy as a measure of system complexity[J]. Proceedings of the National Academy of Sciences ,1991,88(6):2297—2301.
% 理论讲解：https://zhuanlan.zhihu.com/p/494761890/

function kSampleEnValue = kSampleEn(data, dim, r)


% 计算data序列的样本熵
% 输入：
% data：一维数据
% dim：重构维数，一般选择1或2，优先选择2，一般不取m>2
% r：阈值大小，一般选择r=0.1~0.25，（网上流传的SampEn需要再乘以data数据的标准差，本程序不需要，只用输入0.1-0.25间的数字即可）
% 输出：
% kSampleEnValue：样本熵值大小
% 论文：Richman J S, Moorman J R. Physiological time-series analysis using approximate entropy and sample entropy[J]. American Journal of Physiology-Heart and Circulatory Physiology, 2000.
% 理论讲解：https://zhuanlan.zhihu.com/p/519042619/

function FuzEn = kFuzzyEntropy(data,dim,r,n)


% 计算模糊熵，在网传FuzzyEntropy.m基础上修改并注释
% 输入：
% data：一维数据
% dim：重构维数，一般选择1或2，优先选择2，一般不取m>2
% r：阈值大小，一般选择r=0.1~0.25，（网上流传的SampEn需要再乘以data数据的标准差，本程序不需要，只用输入0.1-0.25间的数字即可）
% n：权重，取较小的整数，如2或者3。推荐设置为2
% 输出：
% FuzEn：模糊熵值大小

% 论文：Chen W , Wang Z , Xie H , et al. Characterization of Surface EMG Signal Based on Fuzzy Entropy[J]. IEEE Transactions on Neural Systems and Rehabilitation Engineering, 2007, 15(2):266-272.
% 理论讲解：https://zhuanlan.zhihu.com/p/519633544/

function pe = kPermutationEntropy(data,m,t)


% 计算排列熵，在网传PermutationEntropy.m基础上修改并注释
% 输入：
% data：一维数据
% dim：重构维数，推荐值为3~7，如果能够接受程序运行时间久，可以将m设置得尽量大一些。
% t：时间延迟，取整数，设置为1的情况较多。
% 输出：
% pe：排列值大小

% 论文：Bandt C，Pompe B. Permutation entropy:a natural complexity measure for time series[J]. Physical Review Letters,2002,88(17):174102.
% 理论讲解：https://zhuanlan.zhihu.com/p/519809254/

function enveEn = kEnvelopeEntropy(data,fs)


% 计算包络熵
% 输入：
% data：一维数据
% fs：采样频率，采样频率即每秒钟采集的数据点数，按照实际情况设置
% 输出：
% enveEn：包络熵值
% 论文《自适应最大相关峭度解卷积方法及其在轴承早期故障诊断中的应用》

3.testGenFeaEn2D.m文件（完整版独有）

二维数据特征提取的测试文件，完整版代码中支持对二维数据（维度为m*n，其中m为数据组数，n为每组数据的长度。即每行数据为一组）进行特征提取计算。计算结果中每一行对应一组数据的特征结果。

三、快速开始

1.运行测试脚本

先在MATLAB里打开下载好的文件夹，然后运行testGenFeaEn.m程序，程序运行完毕后如果没报错，且工作区成功得到fea变量，则说明运行环境正常，程序正确。

2.导入数据

复制一个testGenFeaEn.m的文件副本，在这个副本里做如下修改：

根据你的文件类型的不同（excel,txt,csv等），将数据导入MATLAB的方法有所不同。同学们可以看博主针对常用文件的导入方法的这个教程，教程上未包含的数据类型，大家可以再参考这个文档。

如果你已经实现了数据导入，这时候应该拥有了一个一维数据变量，这时候就可以调用函数进行特征提取了。

3.实现特征提取

参照testGenFeaEn.m文件中第二步，genFeatureEn函数的调用方法，根据实际采样频率和特征提取类型进行修改设置，运行程序即可。

四、关于完整版与公开版代码

五、获取公开版程序（需使用电脑浏览器打开）

注：公开版代码需使用MATLAB2022a及以上版本。

六、获取完整版程序（使用电脑浏览器或者手机浏览器打开）

获取通道一（淘宝）：点击此处获取完整版程序

获取通道二（本页面）：点击下面“立即支付”按钮，付款后获取完整版代码下载链接和售后联系方式~本通道处于测试阶段，使用该通道可以额外优惠（仅需43元）。付款完成后刷新一下本页面即可看到下载链接。

（注意支付跳转失败的话，请使用浏览器打开本页面）

七、完整版代码重要更新

20230211 （重要更新）修正求能量熵时可能计算结果为0的bug。

20221203 修正了提取能量熵但未提取样本熵/包络熵时，也会同时计算样本熵/包络熵的bug。

20220830 （重要更新）修正求包络熵过程中的bug隐患

20220710 （重要更新）加入包络熵特征

20220225 完成初版代码

八、常见问题

1.程序运行了很长时间，一直跑不出结果。

原因：这通常是由于分析数据的数据量比较大，并且对奇异谱熵（svdpE）特征值进行提取导致的。由于程序运行奇异值分解时非常吃电脑资源，所以如果提取这个特征，运行时间可能会非常长。

建议：1.如果不是必须求奇异谱熵这个特征，那么建议在程序运行时间过长时不再提取这个特征，这样做可以立竿见影，也不会影响其他特征的提取。2.如果一定要计算奇异谱熵，那么建议对数据进行重采样或者数据截段，减小待分析数据的数据量。

2.求出来的能量熵（eE）的值为0

原因：可能是分析对象数据是经emd分解后的imf分量。需要注意，请不要对数据进行emd分解，然后把emd分解得到的某一imf分量求能量熵，因为在能量熵计算函数里会对数据进行emd分解，imf分量再次进行emd分解将只能得到它本身这一个分量，而对于一个常数求信息熵得到的结果就是0.

建议：请不要对经emd分解后的单一imf分量求能量熵。