Self-tuning Multivariate Variational Mode Decomposition
-
摘要: 多元变分模态分解(MVMD)作为变分模态分解(VMD)的多元扩展，在继承VMD优点的同时，也存在其分解性能很大程度上依赖于两个预置参数——模态数量K 和惩罚系数α的问题。为此，该文提出一种自整定MVMD(SMVMD)算法。SMVMD采取了匹配追踪法的思想，通过频域的能量占比和模态正交性分别自适应地更新K和α。对仿真信号与真实案例的分析结果表明，所提SMVMD方法不仅有效解决了原MVMD的参数整定问题，而且表现出以下优势，(1) 与MVMD相比，SMVMD抗模态混叠的能力更强，且对噪声和α值的变化都具有更好的鲁棒性。(2) 与多元经验模态分解、快速多元经验模态分解和多元变分模态分解这些经典算法相比，SMVMD算法的分解误差最小，分解效果最好。Abstract: The Multivariate Variational Mode Decomposition (MVMD), being an extension of the Variational Mode Decomposition (VMD), inherits the merits of VMD. However, it encounters an issue wherein its decomposition performance relies heavily on two predefined parameters, the number of modes (K) and the penalty factor (α). To address this issue, a Self-tuning MVMD (SMVMD) algorithm is proposed. SMVMD employs the notion of matching pursuit to adaptively update K and α based on energy occupation and mode orthogonality in the frequency domain, respectively. The experimented results of both simulated signals and real cases demonstrate that the proposed SMVMD not only effectively addresses the parameter rectification problem of the original MVMD, but also exhibits the following advantages: (1) SMVMD displays superior resilience to mode-mixing compared to MVMD, along with enhanced robustness to both noise and variations in α-value. (2) In comparison to the classical algorithms of multivariate empirical mode decomposition, fast multivariate empirical mode decomposition, and multivariate variational mode decomposition, SMVMD showcases the lowest decomposition error and the best decomposition effect.
-
1. 引言
多元/多变量/多通道信号是指具有两个或者两个以上通道的信号，表示为X(t)={x1(t),x2(t),⋯,xn(t),t=1,2,⋯,N}，其中n和N分别表示通道数目以及信号长度。随着计算机和传感器的进步，工程应用中出现了大量的非平稳多变量信号[1] ，例如多通道脑电信号[2] 、多通道心电信号[3]、控制系统厂级振荡信号[4]等。信号分解作为处理非平稳信号的重要方法，可以将一个复杂的信号分解成几个有规律的简单模态，进而获取数据中蕴含的有用信息。非平稳多变量数据的信号分解与时频分析应注意以下两点：(1)要具有模态齐整性，即不同通道中具有相同或者相似频率的模态应出现在同一尺度；(2)要能提取出多通道之间的相关信息[5] 。
若直接用单变量信号分解方法逐一分解多变量信号的各个通道，不仅会忽略掉通道之间的相关信息，也不满足模态齐整性[5]。多变量信号分解方法的发展始于复数经验模态分解(Complex Empirical Mode Decomposition, CEMD)。2007 年，Tanaka和Mandic [6]提出了CEMD，但CEMD只能处理复数信号，不具有普适性。接着Rilling 等人[7] 提出了双变量 EMD (Bivariate EMD, BEMD)，并且明确地给出了3维空间中二元信号的极值、均值和包络的有关定义。随后，Rehman和Mandic[8,9]受BEMD的启发，将二变量EMD拓展到三变量EMD (Trivariate EMD, TEMD)以及多变量EMD (Multivariate EMD, MEMD)。MEMD算法被提出后，受到了广泛的关注，目前已经被广泛应用于生物医学[10]、过程控制[11]、地球物理[12]等领域。
由于MEMD的分解耗时会随着信号通道数目的增多而显著增加，Lang等人[13]于2018年提出了快速MEMD算法 (Fast MEMD, FMEMD)。通过证明多变量信号投影与分解的可交换性，将MEMD涉及的多元样条插值问题转化为一元问题，极大地提高了算法的计算效率。以上多元信号分解方法均继承了EMD的优点，但依然有着对采样频率敏感、对噪声鲁棒性差、缺乏严格的数学证明等不足[14]。
与上述基于经验的多变量信号分解方法不同，基于VMD的多变量信号分解方法不仅具有完备的数学框架，还继承了标准VMD的许多理想特性。基于变分框架的多变量信号分解算法包括复数VMD (Complex Variational Mode Decomposition, CVMD)，2维VMD (Two-Dimensional VMD, 2D-VMD)和多变量VMD (Multivariate VMD, MVMD)[5,15,16]。然而CVMD和2D-VMD都不是通用的多变量信号分解算法。因此，Rehman等人[5]于2019年提出了MVMD，用于处理包含两个通道及以上的多变量信号。MVMD定义了多变量调制振荡的概念，并以完全重构输入信号为条件，通过最小化所有模态的带宽之和，构建了一个带约束的变分优化问题。与MEMD相比，MVMD的噪声鲁棒性和抗模态混叠能力更加优秀。由于其优良的特性，目前MVMD已经成功应用于脑电信号检测[17]、风机故障诊断[18]和信号消噪[19]等应用。
然而，MVMD作为一个完全非递归的分解模型，其有效性高度依赖于两个预置参数，即模态数K和惩罚系数α。若K值设置不当，就会导致模态的丢失或者混叠。此外，惩罚系数α的选取也会对MVMD的性能产生很大影响，而MVMD没有提供如何设置α的指导准则。这些都限制了MVMD的适应性。为解决上述问题，提升MVMD算法的应用有效性，本文提出一种自整定MVMD算法(Self-tuning MVMD, SMVMD)。
本文的结构如下：第2节详细介绍SMVMD算法的目标函数和算法流程；第3节对SMVMD的特性进行研究，包括滤波器组结构、对噪声和α值变化的鲁棒性，然后通过一个仿真案例验证该算法的有效性及优越性；第4节以两则真实案例为对象，进一步证明算法的实用性；第5节对全文进行总结与展望。
2. SMVMD
MVMD继承了标准VMD的许多理想性质，如具有良好的数学框架，在对多变量调制振荡的求解过程中通过镜像延拓方法避免了端点效应等。然而，MVMD仍存在如下不足：(1) 需要事先指定模态K的数目。实际应用中，复杂信号固有的模态数量难以预知，若K设置不合适，可能造成过分解或欠分解等问题。(2) 惩罚因子α是固定的。由MVMD的模态更新公式可知，α的大小决定了模态的带宽。在优化求解的初期，需要比较宽的滤波器带宽以捕获正确的中心模态，而在优化求解末期得到的模态已经在期望值附近，这时需要更窄的带宽以滤除噪声等干扰因素[20,21]。因此，在优化过程中更新α的值才是更为合理的。
为解决上述参数整定问题，在最近的相关工作的启发下[20,21]，本文提出SMVMD算法。该算法改进了原MVMD的目标函数，将联合优化框架改为递归框架。具体来说，SMVMD会逐一提取多变量调制振荡，而不是同时提取，模态数K由剩余信号的能量占比确定，无需事先指定。此外，SMVMD算法提出了基于模态正交性的带宽更新法则，使模态的带宽会随着优化过程不断调整。本节将从目标函数及具体步骤两个方面详细介绍所提出的SMVMD算法。
2.1 目标函数
MVMD是令所有通道、所有模态的带宽之和最小，SMVMD则是令每一层/每一个尺度的模态的带宽之和最小，以此构造出一个变分优化问题。
minimize{uk,c},{ωk}{α∑c||∂t[u+k,c(t)e−jωkt]||22+∑c‖ (1) 其中，角标k, c分别用于模态与通道计数， {x_c}(t) 是每个通道对应的信号， \alpha 为惩罚因子。SMVMD采用了维纳滤波器的思想，式(1)中的第1项对模态 {u_{k,c}}(t) 加了平滑度约束，代表每一层所有通道的模态的带宽之和，第2项代表每个通道中的残留信号之和。根据交替方向乘子法，上述复杂的优化问题可以转化成两个子优化问题，即交替更新模态与中心频率。同时根据模态 {u_{k,c}}(t) 与残余信号在频域的正交性，可以推导出惩罚系数 \alpha 的更新公式。有关 {u_{k,c}}(t) ， {\omega _k} 及 \alpha 更新的推导过程如下。
2.1.1 模态 {u_{k,c}}(t) 的更新
\begin{split} u_{k,c}^{n + 1} =\,& \mathop {\arg \min }\limits_{{u_{k,c}}} \Bigr\{ \alpha \left\| {{\partial _t}[u_{k,c}^ + (t){{\text{e}}^{ - {\text{j}}{\omega _k}t}}]} \right\|_2^2 \\ & \left. + \left\| {{x_c}(t) - {u_{k,c}}(t)} \right\|_2^2\right\} \end{split} (2) 变换到频域可得
\begin{split} u_{k,c}^{n + 1}(\omega ) =\,& \arg \min \left\{ \alpha \left\| {\text{j}}\omega [(1 + {{\mathrm{sgn}}} (\omega + {\omega _k})\right.\right.\\ & \left.\cdot {u_{k,c}}(\omega + {\omega _k})] \right\|_2^2 \\ & \left.+ \left\| {{x_c}(\omega ) - {u_{k,c}}(\omega )} \right\|_2^2\right\} \end{split} (3) 去掉式(3)中的二范数，再对右式关于 {u_{k,c}}(\omega ) 求偏导，整理可得
u_{k,c}^{n + 1}(\omega ) = \frac{{{x_c}(\omega )}}{{1 + 2\alpha {{(\omega - {\omega _k})}^2}}} (4) 2.1.2 中心频率 {w_k} 的更新
式(1)中只有第1项与 {\omega _k} 有关，因此在优化过程中只需关注第1项，
\omega _k^{n + 1} = \mathop {\arg \min }\limits_{{\omega _k}} \left\{ \sum\limits_c {\left\| {{\partial _t}[u_{k,c}^ + (t){{\text{e}}^{ - {\text{j}}{\omega _k}t}}]} \right\|_2^2} \right\} (5) 将上式先变换到频域，再去掉2范数可得
\omega _k^{n + 1}(\omega ) = \arg \min \sum\limits_c {\int\limits_0^\infty {{{(\omega - {\omega _k})}^2}|{u_{k,c}}(\omega ){|^2}} } {\text{d}}\omega (6) 对式(6)中的等式右侧关于 {\omega _k} 求偏导，并令导数得0，整理可得
\omega _k^{n + 1} \leftarrow \frac{{\displaystyle\sum\limits_c {\displaystyle\int\limits_0^\infty {\omega |u_{^{k,c}}^{n + 1}(\omega ){|^2}{\text{d}}\omega } } }}{{\displaystyle\sum\limits_c {\displaystyle\int\limits_0^\infty {|u_{k,c}^{n + 1}(\omega ){|^2}{\text{d}}\omega } } }} (7) 2.1.3 惩罚系数 \alpha 的更新
由时域正交的两个信号各自变换到频域后依然正交，可知 u_{k,c}^*(\omega ){r_{k,c}}(\omega ) = 0 。同时结合 {x_c}(\omega ) = {u_{k,c}}(\omega ) + {r_{k,c}}(\omega ) 可得
\frac{{u_{k,c}^*(\omega ){u_{k,c}}(\omega )}}{{u_{k,c}^*(\omega ){x_c}(\omega )}} = 1 (8) 将 {u_{k,c}}(\omega ) 的更新公式，即式(4)代入式(8)中，整理可得
\begin{split} & \alpha = \\ & \frac{{x_c^*(\omega ){{\left(\dfrac{1}{{1/\alpha + 2{{(\omega - {\omega _k})}^2}}}\right)}^*}\left(\dfrac{1}{{1/\alpha + 2{{(\omega - {\omega _k})}^2}}}\right){x_c}(\omega )}}{{x_c^*(\omega ){{\left(\dfrac{1}{{1/\alpha + 2{{(\omega - {\omega _k})}^2}}}\right)}^*}{x_c}(\omega )}} \end{split} (9) 由于 \dfrac{1}{{1/{\alpha ^n} + 2{{(\omega - \omega _k^n)}^2}}} = {\alpha ^n}\dfrac{1}{{1 + 2{\alpha ^n}{{(\omega - \omega _k^n)}^2}}} ，因此式(9)可以表示为
\begin{split} & {\alpha ^{n + 1}} =\\ & \frac{{{\alpha ^n}{{\left(\dfrac{{{x_c}(\omega )}}{{1 + 2{\alpha ^n}{{(\omega - \omega _k^n)}^2}}}\right)}^*}\left(\dfrac{{{x_c}(\omega )}}{{1 + 2{\alpha ^n}{{(\omega - \omega _k^n)}^2}}}\right)}}{{{{\left(\dfrac{{{x_c}(\omega )}}{{1 + 2{\alpha ^n}{{(\omega - \omega _k^n)}^2}}}\right)}^*}{x_c}(\omega )}} \end{split} (10) 将式(4)代入式(10)，得到
\alpha _c^{n + 1} = \alpha _c^n\frac{{u_{k,c}^{n + 1}{{(\omega )}^*}u_{k,c}^{n + 1}(\omega )}}{{u_{k,c}^{n + 1}{{(\omega )}^*}{r_{k,c}}(\omega )}} (11) 随着模态 {u_{k,c}}(\omega ) 的不断更新，式(11)中分母的值会逐渐减小，使 \alpha 值逐渐增大，进而模态的带宽随更新过程不断减小，满足前文带宽更新的要求。
2.2 SMVMD算法
SMVMD算法共有两层循环，外循环通过剩余信号的能量占比控制模态数K，外循环的终止参数 \delta 可以根据实际噪声水平调整。内循环与MVMD相似，主要进行模态、中心频率以及 \alpha 值的更新。SMVMD的具体算法见算法1。
表 1 SMVMD算法初始化： \{ u_{k,c}^1(\omega )\} ,\{ \omega _k^1\} ,n \leftarrow 0 令 k = 1,{{\mathbf{r}}_1}(\omega ) = {\mathbf{x}}(\omega ) while \left\| {{{\mathbf{r}}_{{k}}}(\omega )} \right\|_2^2/\left\| {{\mathbf{x}}(\omega )} \right\|_2^2 > \delta do 令n=1，初始化中心频率 \{ \omega _k^1\} while \displaystyle\sum\limits_c {\dfrac{{||u_{k,c}^n(\omega ) - u_{k,c}^{n - 1}(\omega )||_2^2}}{{||u_{k,c}^{n - 1}(\omega )||_2^2}}} \gt \varepsilon do n \leftarrow n + 1 for c = 1:C do u_{k,c}^{n + 1}(\omega ) \leftarrow \dfrac{{{r_{k,c}}(\omega )}}{{1 + 2\alpha _c^n{{(\omega - \omega _k^n)}^2}}} end for \omega _k^{n + 1} \leftarrow \dfrac{{\displaystyle\sum\limits_c {\int\limits_0^\infty {\omega |u_{^{_{k,c}}}^{n + 1}(\omega ){|^2}{\text{d}}\omega } } }}{{\displaystyle\sum\limits_c {\int\limits_0^\infty {|u_{k,c}^{n + 1}(\omega ){|^2}{\text{d}}\omega } } }} for c = 1:C do \alpha _c^{n + 1} = \alpha _c^n\dfrac{{u_{k,c}^{n + 1}{{(\omega )}^*}u_{k,c}^{n + 1}(\omega )}}{{u_{k,c}^{n + 1}{{(\omega )}^*}{r_{k,c}}(\omega )}} end for end while 令 {u_{k,c}}(\omega ) = u_{k,c}^{n + 1}(\omega ),{\omega _k} = \omega _k^{n + 1} ，同时更新 {{\mathbf{r}}_{k + 1}}(\omega ) = {\mathbf{x}}(\omega ) - {{\mathbf{u}}_k}(\omega ) end while 3. SMVMD的性质与比较
3.1 滤波器组结构
滤波器组结构是信号分解方法的一个重要特性，许多多变量信号分解算法都对其滤波器组特性进行过研究，如MEMD、MVMD[5,22]。类比MEMD与MVMD，本节对SMVMD算法进行该特性分析，并与MVMD的结果进行对比。
采用SMVMD分解数据长度为N=1024的四通道高斯白噪声过程(服从均值为0，标准差为 \sigma 的正态分布)。为了消除随机因素，重复100次实验，并将MVMD也应用于同一数据集。 为了方便比较，将MVMD与SMVMD的模态数K统一设为8，二者在100次实验下的平均功率谱密度(Power Spectral Density, PSD)如图1所示。可以看出，两种算法的本质都是一组维纳滤波器。对比二者的滤波器组特性，MVMD的相邻滤波器之间有明显的重叠，而SMVMD的几乎没有重叠。根据文献[23]，模态频谱之间的重叠面积越大则发生模态混叠的概率越高，因此SMVMD抗模态混叠的能力更佳。
3.2 对噪声及 \alpha 的鲁棒性实验
为调查SMVMD算法关于噪声以及 \alpha 值变化时的鲁棒性，本节采用一组经典的多变量信号对SMVMD进行测试，式(12)所示的信号经常用于仿真实验分析[1,5,20]。该信号由复杂的多分量组成，可以很好地测试算法的模态对齐性及噪声鲁棒性。
\left. \begin{aligned} & {{\boldsymbol{X}}_1}(t) = [{x_1}(t),{x_2}(t),{x_3}(t)] \\ & {x_1}(t) = (1 + 0.5\cos (2\pi \times t))\cos (2\pi \times 5t) \\ & \qquad\quad + \cos (2\pi \times 15t) + w(t) \\ & {x_2}(t) = \cos (2\pi \times t)\cos (2\pi \times 5t) + \cos (2\pi \times 15t) \\ & \qquad\quad + w(t) \\ & {x_3}(t) = 2\cos (2\pi \times 15t) + w(t) \end{aligned} \right\} (12) 其中t从0到1变化，采样频率为1000 Hz，信号的数据长度为N=1000， w\left( t \right) 是高斯白噪声，服从正态分布 w\left( t \right) \sim \mathcal{N}(0,{\sigma ^2}) 。
为公平比较，将MVMD与SMVMD的 \alpha 初始值设置为1000(这是文献[5]中的默认值)，同时设置MVMD的模态数为K=2。令高斯白噪声的标准差 \sigma 从0变化到0.5，间隔为0.1，通过观察两个余弦分量频率的收敛性来判断两个算法对噪声的鲁棒性。在100次实验的基础上进行集总平均，结果如图2(a)与图2(b)所示。可以看出，当噪声水平逐渐增强时，SMVMD最终收敛得到的中心频率更贴近其真实值。而MVMD分解信号X1(t)所得到的第1个模态的中心频率 {\omega _1} 偏离了正确的值。
其次，探讨MVMD与SMVMD对 \alpha 值变化的鲁棒性。不失一般性，固定信号X1(t)中噪声的标准差为0.1，同时令 \alpha 从1变化到10000，变化间隔为1。不同 \alpha 值对MVMD和SMVMD分解效果的影响如图2(c)与图2(d)所示。可以看到，SMVMD最终收敛得到的中心频率与正确的中心频率几乎重合。这是因为SMVMD算法的 \alpha 具有自整定性，在更新过程中会逐渐调整到合适的值，从而保证分解的准确性。相反，MVMD的收敛曲线没有达到预期的结果，对参数 \alpha 的值较为敏感。
3.3 与典型算法的对比研究
为进一步验证SMVMD算法的分解有效性，将其与MVMD, MEMD, FMEMD等经典的多元信号分解方法进行对比。上述算法涉及的参数均采用默认值[5,9,13]。采用上述4种算法分解信号X1(t)(其中噪声的标准差设为0.1)。在此，为了进行定量比较，引入绝对误差之和(Sum of the Absolute Difference, SAD)来准确地衡量分解性能，定义为
{\text{SAD}} = \sum\limits_{c = 1}^C {\sum\limits_{k = 1}^K {\sum\limits_{n = 1}^N {|{x_{c,k}}(n) - {u_{c,k}}(n)|} } } (13) 其中， {x_{c,k}}(n) 和 {u_{c,k}}(n) 分别代表第k个真实模态和分解后的模态，K, C和N分别是模态数、通道数和数据长度。算法对应的SAD值越小，分解效果越好。4种算法分解信号X1(t)对应的SAD值如表1所示，与MVMD, MEMD, FMEMD相比，所提SMVMD算法的分解误差最小，分解准确率最高。
表 1 4种算法的SAD值算法 MVMD MEMD FMEMD SMVMD SAD 858.30 854.64 993.42 226.11 4. 真实案例分析
4.1 EEG信号\alpha 节律提取
本节通过一个真实案例说明所提出的SMVMD算法在信号分解和时频分析应用中的实用性。该真实案例为一组4通道的脑电(Electroencephalogram, EEG)信号，采集自一项受试者保持闭眼和放松状态下的实验[5]。已经证实，在闭眼且放松状态下，EEG信号中可检测出频率范围为8～12 Hz的 \alpha 节律。
图3显示了上述EEG信号的时域图，以及由SMVMD分解得到的模态时域图和相应功率谱图。参考文献[5]，此处仅选择模态u2, u3, u5进行展示。由图3(b)可以看出，EEG信号中的 \alpha 节律位于SMVMD分解得到的模态u2中。同时还可以观察到，与u2相对应的所有通道均包含了 \alpha 节律，这验证了SMVMD在处理实际信号时的模态齐整性。该结果与MVMD的分解结果一致，进一步验证了所提算法的有效性。
值得注意的是，两种算法分解得到的模态u5是由交流电源引起的伪影，其频率约为50 Hz[1,5]。然而，MVMD检测到的电源频率约为40 Hz，与事实不符[5]。对比之下，SMVMD正确地提取了该电源频率，如图3(b4)所示。综上，本文所提出的SMVMD在这则真实案例中表现出了更好的性能。
4.2 控制系统厂级振荡分析
振荡是过程控制系统中最常见的异常现象之一。虽然振荡通常在1个源回路中产生，但它经常通过相互连接的回路传播，从而引起厂级范围的振荡。这会造成生产资源浪费、产品质量波动甚至会损害系统的稳定性和安全性[24,25]。因此，检测和分析厂级振荡是十分必要的。
此处将SMVMD应用于一个选矿厂浮选回路中观察到的厂级振荡数据集上[26]。该浮选回路由2个并联的库组成，每个库各有7个串联的浮选槽(库1：FT001至FT007层；库2：FT008至FT014层)。数据的采样间隔为10 s，为了便于显示和分析结果，这里选择了9个单元(即FT002至FT010)，分别标记为变量x1到x9。图4的第1行显示了9个变量的测量结果，其余两行是SMVMD分解的模态(u1, u2)。
参考文献[27]，本文在此选用归一化相关系数和稀疏指数(Sparsity Index, SI)两个指标对模态进行振荡检测。归一化相关系数(14)可以剔除由于模态齐整性产生的虚假模式[24]。
{\theta _{k,c}} = \frac{{{\rho _{k,c}}}}{{\max \{ {\rho _{1,c}},{\rho _{2,c}}, \cdots ,{\rho _{K,c}}\} }} (14) 其中， {\rho _{k,c}} 是模态 {u_{k,c}} 与信号 {x_c} 之间的相关系数。仅保留 {\theta _{k,c}} \gt {T_\theta } 的模态用于分析振荡特性，其中阈值一般设置为 {T_\theta } = 0.5 [27]。
稀疏指数的计算公式为
{\text{SI}} = \frac{{\sqrt N - \left(\displaystyle\sum\limits_{n = 1}^N {|\hat u _{k,c}^n|} /\sqrt {{{\displaystyle\sum\limits_{n = 1}^N {|\hat u _{k,c}^n|} }^2}} \right)}}{{\sqrt N - 1}} (15) 其中， u_{k,c}^{ \wedge} 是待检测信号 {u_{k,c}} 的傅里叶变换，N表示信号的长度。振荡的模态在频域中会表现出明显的峰值，对应的SI值接近于1。根据文献[28]，SI值大于0.58的模态可以被认定为是振荡的。
模态u1, u2对应的归一化相关系数与稀疏指数的取值如表2所示。我们对SMVMD分解后的模态检测到的振荡结果进行了可视化，结果如图5所示，其中黑色方块对应于振荡模态，白色方块对应于非振荡模态。可以看出，x1～x6通道的模态u2存在一个相同频率成分的振荡，而x7～x9通道的模态u1存在一个相同频率成分的振荡，这与先前工作的分析结果相一致[26]。SMVMD在该厂级振荡案例中的成功应用再次验证了所提算法的实用性和有效性。
表 2 归一化相关系数与稀疏指数模态 x1 x2 x3 x4 x5 x6 x7 x8 x9 {\theta _{k,c}} u1 0.446 0.497 0.398 0.441 0.371 0.497 1 1 1 u2 1 1 1 1 1 1 0.491 0.419 0.493 SI u1 – – – – – – 0.8936 0.8853 0.8769 u2 0.9117 0.8818 0.8814 0.8866 0.8361 0.8306 - - - 5. 结束语
本文提出了一种自整定多元变分模态分解算法。该方法采用了匹配追踪思想，用剩余信号的能量占比和模态在频域的正交性分别自适应地更新模态数K和惩罚因子 \alpha ，解决了原MVMD算法需要预置参数K和 \alpha 的问题。通过代表性的仿真与真实案例，不仅验证了所提算法分解得到的子信号具有模态齐整性，还证明了该算法的有效性及其相较其他经典多变量信号分解方法的优越性。本文提出的算法在多通道脑电信号和控制系统厂级振荡信号的分析中表现出显著优势，未来的研究工作将深入探索SMVMD算法在其他工程领域的应用。同时可参考FMEMD的思想，从将MVMD涉及的多元优化问题转化为一元问题入手，提高算法的分解速度，解决MVMD难以处理大通道数多元信号的问题。
-
1 SMVMD算法
初始化： \{ u_{k,c}^1(\omega )\} ,\{ \omega _k^1\} ,n \leftarrow 0 令 k = 1,{{\mathbf{r}}_1}(\omega ) = {\mathbf{x}}(\omega ) while \left\| {{{\mathbf{r}}_{{k}}}(\omega )} \right\|_2^2/\left\| {{\mathbf{x}}(\omega )} \right\|_2^2 > \delta do 令n=1，初始化中心频率 \{ \omega _k^1\} while \displaystyle\sum\limits_c {\dfrac{{||u_{k,c}^n(\omega ) - u_{k,c}^{n - 1}(\omega )||_2^2}}{{||u_{k,c}^{n - 1}(\omega )||_2^2}}} \gt \varepsilon do n \leftarrow n + 1 for c = 1:C do u_{k,c}^{n + 1}(\omega ) \leftarrow \dfrac{{{r_{k,c}}(\omega )}}{{1 + 2\alpha _c^n{{(\omega - \omega _k^n)}^2}}} end for \omega _k^{n + 1} \leftarrow \dfrac{{\displaystyle\sum\limits_c {\int\limits_0^\infty {\omega |u_{^{_{k,c}}}^{n + 1}(\omega ){|^2}{\text{d}}\omega } } }}{{\displaystyle\sum\limits_c {\int\limits_0^\infty {|u_{k,c}^{n + 1}(\omega ){|^2}{\text{d}}\omega } } }} for c = 1:C do \alpha _c^{n + 1} = \alpha _c^n\dfrac{{u_{k,c}^{n + 1}{{(\omega )}^*}u_{k,c}^{n + 1}(\omega )}}{{u_{k,c}^{n + 1}{{(\omega )}^*}{r_{k,c}}(\omega )}} end for end while 令 {u_{k,c}}(\omega ) = u_{k,c}^{n + 1}(\omega ),{\omega _k} = \omega _k^{n + 1} ，同时更新 {{\mathbf{r}}_{k + 1}}(\omega ) = {\mathbf{x}}(\omega ) - {{\mathbf{u}}_k}(\omega ) end while 
下载: 导出CSV
表 2 归一化相关系数与稀疏指数
模态 x1 x2 x3 x4 x5 x6 x7 x8 x9 {\theta _{k,c}} u1 0.446 0.497 0.398 0.441 0.371 0.497 1 1 1 u2 1 1 1 1 1 1 0.491 0.419 0.493 SI u1 – – – – – – 0.8936 0.8853 0.8769 u2 0.9117 0.8818 0.8814 0.8866 0.8361 0.8306 - - -
下载: 导出CSV
-
[1] CHEN Qiming, LANG Xun, XIE Lei, et al. Multivariate intrinsic chirp mode decomposition[J]. Signal Processing, 2021, 183: 108009. doi: 10.1016/j.sigpro.2021.108009. [2] ZAHRA A, KANWAL N, REHMAN N U, et al. Seizure detection from EEG signals using multivariate empirical mode decomposition[J]. Computers in Biology and Medicine, 2017, 88: 132–141. doi: 10.1016/j.compbiomed.2017.07.010. [3] HAN G, LIN B, and XU Z. Electrocardiogram signal denoising based on empirical mode decomposition technique: An overview[J]. Journal of Instrumentation, 2017, 12(3): P03010. doi: 10.1088/1748-0221/12/03/P03010. [4] 王宇红, 高志兴. 改进多维本质时间尺度分解的厂级振荡检测[J]. 控制工程, 2022, 29(10): 1835–1840. doi: 10.14107/j.cnki.kzgc.CAC2020-1557.WANG Yuhong and GAO Zhixing. Plant-wide oscillation detection based on improved multivariate intrinsic time-scale decomposition[J]. Control Engineering of China, 2022, 29(10): 1835–1840. doi: 10.14107/j.cnki.kzgc.CAC2020-1557. [5] REHMAN N U and AFTAB H. Multivariate variational mode decomposition[J]. IEEE Transactions on Signal Processing, 2019, 67(23): 6039–6052. doi: 10.1109/TSP.2019.2951223. [6] TANAKA T and MANDIC D P. Complex empirical mode decomposition[J]. IEEE Signal Processing Letters, 2007, 14(2): 101–104. doi: 10.1109/LSP.2006.882107. [7] RILLING G, FLANDRIN P, GONCALVES P, et al. Bivariate empirical mode decomposition[J]. IEEE Signal Processing Letters, 2007, 14(12): 936–939. doi: 10.1109/LSP.2007.904710. [8] REHMAN N U and MANDIC D P. Empirical mode decomposition for trivariate signals[J]. IEEE Transactions on Signal Processing, 2010, 58(3): 1059–1068. doi: 10.1109/TSP.2009.2033730. [9] REHMAN N and MANDIC D P. Multivariate empirical mode decomposition[J]. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2010, 466(2117): 1291–1302. doi: 10.1098/rspa.2009.0502. [10] ASGHAR M A, KHAN M J, RIZWAN M, et al. AI inspired EEG-based spatial feature selection method using multivariate empirical mode decomposition for emotion classification[J]. Multimedia Systems, 2022, 28(4): 1275–1288. doi: 10.1007/s00530-021-00782-w. [11] LANG Xun, ZHANG Yufeng, XIE Lei, et al. Detrending and denoising of industrial oscillation data[J]. IEEE Transactions on Industrial Informatics, 2023, 19(4): 5809–5820. doi: 10.1109/TII.2022.3188844. [12] HUANG Guoqing, PENG Liuliu, KAREEM A, et al. Data-driven simulation of multivariate nonstationary winds: A hybrid multivariate empirical mode decomposition and spectral representation method[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2020, 197: 104073. doi: 10.1016/j.jweia.2019.104073. [13] LANG Xun, ZHENG Qian, ZHANG Zhiming, et al. Fast multivariate empirical mode decomposition[J]. IEEE Access, 2018, 6: 65521–65538. doi: 10.1109/ACCESS.2018.2877150. [14] 蔡念, 黄威威, 谢伟, 等. 基于互补自适应噪声的集合经验模式分解算法[J]. 电子与信息学报, 2015, 37(10): 2383–2389. doi: 10.11999/JEIT141632.CAI Nian, HUANG Weiwei, XIE Wei, et al. Ensemble empirical mode decomposition base on complementary adaptive noises[J]. Journal of Electronics & Information Technology, 2015, 37(10): 2383–2389. doi: 10.11999/JEIT141632. [15] WANG Yanxue, LIU Fuyun, JIANG Zhansi, et al. Complex variational mode decomposition for signal processing applications[J]. Mechanical Systems and Signal Processing, 2017, 86: 75–85. doi: 10.1016/j.ymssp.2016.09.032. [16] ZOSSO D, DRAGOMIRETSKIY K, BERTOZZI A L, et al. Two-dimensional compact variational mode decomposition[J]. Journal of Mathematical Imaging and Vision, 2017, 58(2): 294–320. doi: 10.1007/s10851-017-0710-z. [17] 孟明, 闫冉, 高云园, 等. 基于多元变分模态分解的脑电多域特征提取方法[J]. 传感技术学报, 2020, 33(6): 853–860. doi: 10.3969/j.issn.1004-1699.2020.06.011.MENG Ming, YAN Ran, GAO Yunyuan, et al. Multi-domain feature extraction of EEG based on multivariate variational mode decomposition[J]. Chinese Journal of Sensors and Actuators, 2020, 33(6): 853–860. doi: 10.3969/j.issn.1004-1699.2020.06.011. [18] YAN Xiaoan, LIU Ying, XU Yadong, et al. Multichannel fault diagnosis of wind turbine driving system using multivariate singular spectrum decomposition and improved Kolmogorov complexity[J]. Renewable Energy, 2021, 170: 724–748. doi: 10.1016/j.renene.2021.02.011. [19] ZHANG Yijie, ZHANG Haoran, YANG Yang, et al. Seismic random noise separation and attenuation based on MVMD and MSSA[J]. IEEE Transactions on Geoscience and Remote Sensing, 2022, 60: 5908916. doi: 10.1109/TGRS.2021.3131655. [20] CHEN Qiming, CHEN Junghui, LANG Xun, et al. Self-tuning variational mode decomposition[J]. Journal of the Franklin Institute, 2021, 358(15): 7825–7862. doi: 10.1016/j.jfranklin.2021.07.021. [21] CHEN Shiqian, YANG Yang, PENG Zhike, et al. Adaptive chirp mode pursuit: Algorithm and applications[J]. Mechanical Systems and Signal Processing, 2019, 116: 566–584. doi: 10.1016/j.ymssp.2018.06.052. [22] REHMAN N U and MANDIC D P. Filter bank property of multivariate empirical mode decomposition[J]. IEEE Transactions on Signal Processing, 2011, 59(5): 2421–2426. doi: 10.1109/TSP.2011.2106779. [23] LANG Xun, ZHANG Yufeng, XIE Lei, et al. Use of fast multivariate empirical mode decomposition for oscillation monitoring in noisy process plant[J]. Industrial & Engineering Chemistry Research, 2020, 59(25): 11537–11551. doi: 10.1021/acs.iecr.9b06351. [24] CHEN Qiming, WEN Qingsong, WU Xialai, et al. Detection and time–frequency analysis of multiple plant-wide oscillations using adaptive multivariate intrinsic chirp component decomposition[J]. Control Engineering Practice, 2023, 141: 105715. doi: 10.1016/j.conengprac.2023.105715. [25] YU Wanke, ZHAO Chunhui, and HUANG Biao. MoniNet with concurrent analytics of temporal and spatial information for fault detection in industrial processes[J]. IEEE Transactions on Cybernetics, 2022, 52(8): 8340–8351. doi: 10.1109/TCYB.2021.3050398. [26] LINDNER B, AURET L, and BAUER M. A systematic workflow for oscillation diagnosis using transfer entropy[J]. IEEE Transactions on Control Systems Technology, 2020, 28(3): 908–919. doi: 10.1109/TCST.2019.2896223. [27] CHEN Qiming, LANG Xun, LU Shan, et al. Detection and root cause analysis of multiple plant-wide oscillations using multivariate nonlinear chirp mode decomposition and multivariate granger causality[J]. Computers & Chemical Engineering, 2021, 147: 107231. doi: 10.1016/j.compchemeng. 2021.107231. [28] AFTAB M F, HOVD M, and SIVALINGAM S. Detecting non-linearity induced oscillations via the dyadic filter bank property of multivariate empirical mode decomposition[J]. Journal of Process Control, 2017, 60: 68–81. doi: 10.1016/j.jprocont.2017.08.005. -
计量
- 文章访问数: 424
- HTML全文浏览量: 168
- PDF下载量: 47
- 被引次数: 0
下载:
