【熵与特征提取】基于“信息熵”的特征指标及其MATLAB代码实现（功率谱熵、奇异谱熵、能量熵）

《三体》中对监听员的日常工作有这样一段描述：

...1379 号监听站已经存在了上千年，像这样的监听站，在三体世界中有⼏千个，它们全神贯注地聆听着宇宙间可能存在的智慧⽂明的信息...
... 1379 号监听员最不愿意看的，就是显⽰器上缓缓移动的那条曲线，那是监听系统接收到的宇宙电波的波形，⽆意义的噪声...
...但今天，当监听员扫了一眼示波器后，发现有些异样。即使是专业人员，也很难仅凭肉眼看出波形是否携带信息，但监听员对宇宙噪声的波形太熟悉了，眼前移动的波形，似乎多了某种说不出来的东西，这条起伏的细线像是有了灵魂。他敢肯定，眼前的电波是被智能调剂的...

这个可怜的监听员盯了上千年的示波器，但是如果他学过关于熵的知识，可以很轻松地编写一段信息熵特征提取的程序，实现外星文明的自动监听。看来三体星的高等教育还是不够普及。

话说回来，三体监听员能用肉眼识别信号是否是白噪声，也实在是非同一般。下边这两张图中，其中一张是线性分布的白噪声信号，一张是包含“智能调剂”信息的信号，你能做出分辨吗？

两段信号都采集了200个数据点，你能看出哪段信号中包含了“智能调剂”过的信息吗？

人类的眼睛还不行，但是“信息熵”可以。

我们开始吧。

1.关于信息熵

信息熵借用了熵的概念。信息熵的来历和热力学的熵有所不同，之所以信息熵也叫“熵”据说是当年香农提出此概念时参考了热力学的熵，并且两者表达式上有所相似。

信息熵可以很好对“信息量”进行定量描述，其基本思想是：发生概率低的事件，包含的信息量更高。

1.1 举个例子

显然，发生概率最低的“被三体人抓走了”，包含的信息量是最多的，就像是“此处省略了一万字”。

香农先生把信息量进行了量化，事件的自信息[1]被定义为：

如果对数log的底数为e，那么 的单位就是奈特(nats)；如果以2为底数，单位就是比特(bit)。

信息被量化之后，就可以对信息熵下定义了：信息熵是信息量的期望。对于离散量其公式为：

具体到上边小明没上学的例子中，可以具体算出：

此时的信息熵总和就是

1.2 思考一下

确定事件的信息熵：如果小明没来上课的原因100%是因为睡过了，那信息熵是多少呢？代入公式后就能发现，是0。所以确定的事件信息熵值为零。

时间不确定时的信息熵：如果小明没来上课的三个原因概率相等，都等于33.33%，带入公式得到信息熵等于1.0986，这是这个例子中的信息熵极值。

所以我们可以得到一个结论：+

信息熵越大，代表不确定性越大，信号中包含的信息量越少。

换句话说，信号越无秩序/越接近于白噪声，则信息熵越大。

2. 信息熵的三种求法

回到三体监听员的例子。根据上述结论，当监听员检测到信号信息熵突然变小时，就说明此时的信号中包含了特别的信息。

但是，直接对时域信号求信息熵并不是一个好点子，要在求信息熵之前对原始信号进行一些处理，常用的预先处理方法包括①求功率谱、②奇异值分解以及③求信号分解分量的能量等。

2.1 功率谱熵

有些信号的特征比较容易从频域空间反映出来，比如包含了周期性波动的信号。通过频域计算方法可以获得信号在频域空间内的能量分布情况，将这种能量分布与信息熵结合即可以实现信号在频域空间的定量描述。

功率谱熵是用来表示信号能量在功率谱划分下的不确定性，当信号中频率组成简单时，功率谱集中在部分频率成分，对应的频率谱线也会比较少，与之对应的组分概率就会越少，从而导致得到的功率谱熵的值就会变小。与之相反，如果信号较为复杂时，信号对应的功率谱越分散，对应的功率谱线会增多，得到的功率谱熵值就会变大，因此，功率谱熵是对信号在频域上能量分布的复杂程度的定量描述。 [2]
信号越接近白噪声，功率谱越分散，功率谱熵越大。

下图就是某个宇宙信号的功率谱图，可以明显看到在200Hz的位置存在尖峰，这就是一个明显的“确定性”信息。

某个“宇宙信号”的功率谱

功率谱熵的求法分为三步：

a. 求信号的功率谱。

b. 求各个频段的能量占总能量的比重大小（也可以将每个频率点视作一个频段），得到概率密度函数 。

c. 将概率密度函数带入公式（2）即可。

验证一下：

两段“宇宙信号”中有一段是白噪声信号。

如果监听员试着不断把截取长度为200的信号放到功率谱熵公式里计算，会发现其功率谱熵在4.37-4.45的范围内跳动（探测器采集到的是白噪声）。此时监听员可以设置一个报警阈值，比如4.2，然后就可以去喝咖啡了。

当另外一段宇宙信号输入进来时，功率谱熵是3.75，低于报警阈值，检测器报警。

可靠又简便。

2.2 奇异谱熵

奇异值分解的概念在专栏里还没有介绍过，不过相关文章不少，推荐一篇建议先读一下：

漫漫成长：奇异值分解（SVD）

简单来说，奇异谱熵方法是一种通过奇异谱和信息熵相结合的分析方法。这个方法中将奇异值分解结果作为信息熵求解对象。其适用于采样点数较少、采集到的信号中含有噪声的应用场合。

信号的奇异谱分析是以动力学分析为基础的现代谱分析技术，其基本思想是通过对系统的时域信号序列进行相空间重构和奇异值分解，获取其内在复杂性特征。在信号的奇异谱分析基础上，计算奇异谱的信息熵（简称奇异谱熵），则可定量描述时间序列的复杂状态特征。 [3]

奇异谱熵的求法分为四步：

a. 将原始信号  处理成矩阵A：

x为信号数据值，N代表数据长度

b. 对矩阵A进行奇异值分解，得到奇异值谱 。

c. 求各个奇异值占全部奇异值之和的比重大小，得到概率密度函数 。

d. 将概率密度函数带入公式（2）即可。

2.3 能量熵

能量熵有多种求法，因为能量熵通常与信号分解算法相结合，比如EMD、CEEMD、VMD、小波分解等等。

假设一组信号经过分解得到M个分量［S1,S2, …,Sm］，计算每个分量的能量值。以每个分量的能量占总能量的比例作为概率密度函数，计算信息熵。

使用能量熵的方法可以将各种分解方法的特色融合进来，比如结合EMD分解，由于IMF分量包含了不同频段的信息，原始信号处于不同状态时，信号在某些分量上的能量分布就会发生改变；如果结合小波分析，可以将信号时域信息和频域信息相结合。

能量熵的求法分为三步：

a. 使用某种方法将原始信号进行分解

b.分别计算每个分解出来的分量的能量值（即平方和）。

c.计算每个分量的能量占总能量的比例，得到概率密度函数 。

d.将概率密度函数带入公式（2）即可。

2.4 来自地球的信号

揭晓一下谜底，第一张图中左边的图为加入了“智能调剂”，即来自地球的信号，你猜对了吗？

左边的信号中叠加了一个正弦波

使用上述功率谱熵、奇异谱熵、能量熵得到的值分别如下：

奇异值谱和能量谱在分辨上述案例的时候效果相对弱一些。这是跟处理的信号特征相关的，不代表奇异值谱法和能量谱在所有场合都不好用。

3.MATLAB代码实现

网上没有靠谱的功率谱熵、奇异谱熵、能量熵程序，笔者按照理论编写了相关程序。

为了特征提取代码的易用性，笔者对一系列熵特征提取进行了封装，包括上边添加注释的代码都集中到一起。由于搞科研写论文时，对特征提取的需要往往是集中性的、多种类的、需求各异的，所以我把信息熵特征值和后边将会讲到的几种熵特征进行了打包：

熵特征值共7个—— 功率谱熵、奇异谱熵、能量熵、近似熵、样本熵、排列熵、模糊熵

以上7种全都集中到一个封装函数里，实现一行代码完成特征提取。

比如提取数据“功率谱熵”就可以像这样写：

fea = genFeatureEn(data,{'eE'}) %对data求近似熵

如果提取数据“功率谱熵、奇异谱熵、能量熵、近似熵、样本熵、排列熵、模糊熵”这全部7种特征，就可以这样写：

fea =genFeatureEn(data,{'psdE','svdpE','eE', 'ApEn', 'SpEn','PeEn','FuzzyEn'});  
%调用genFeature函数，完成特征提取，算出的特征值会保存在fea变量里

也就是说需要提取哪个特征，在函数中直接指定就可以了。输出的fea变量里就会得到相应的这些特征值，顺序也是与输入的排序保持一致的。

这个函数的介绍如下：

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为排列熵的时间延迟
% featureNamesCell：拟进行特征提取的特征名称，该变量为cell类型，其中包含的特征名称为字符串，特征名称需要在下边列表中：
% 目前支持的特征（2022.5.23，共7种）：
%      psdE：功率谱熵
%      svdpE：奇异谱熵
%      eE：能量熵
%      ApEn：近似熵
%      SampleEn：样本熵
%      FuzzyEn：模糊熵
%      PerEn：排列熵
% 
% 输出：
% fea：数据data的特征值数组，其特征值顺序与featureNamesCell一一对应
% 

需要上边这个函数文件以及测试代码的同学，可以在下述链接中获取：

特征提取（熵特征指标） | 工具箱文档

这个函数会持续更新，有哪些想要加进去的特征指标，同学们可以在评论区留言，笔者会考虑纳入到这个“特征提取指标全家桶”中。

参考

1. ^信息论、交叉熵 - 韬小虾 - 博客园
2. ^矿用带式输送机托辊远程故障诊断系统
3. ^基于信息熵的实时信号测量技术及其应用