CN107679499A - 基于变分模态分解的单通道无线电混合信号分离方法 - Google Patents

Abstract

本发明公开了一种基于变分模态分解的单通道无线电混合信号分离方法。该方法通过改进的量子粒子群优化算法对变分模态分解算法的分量个数和惩罚因子两个参数同时进行优化，确定变分模态分解算法中的分量个数和惩罚因子的值。之后首次利用参数优化的变分模态分解算法对单通道无线电混合信号进行分离，由于待处理的单通道无线电混合信号的信噪比较低，分离后所得各个无线电信号信噪比有所提高，但仍有大量噪声的存在，因此将分离后的各个无线电信号依次使用参数优化的变分模态分解算法进行二次分离，二次分离后所得各个无线电信号比首次分离后所得相应无线电信号可信度高。

Description

基于变分模态分解的单通道无线电混合信号分离方法

技术领域

本发明涉及单通道无线电混合信号分离技术领域，具体是一种基于变分模态分解的单通道无线电混合信号分离方法。

背景技术

在实际应用中，由于受到设备造价、安装条件等问题的制约，往往存在仅仅是由单个传感器接收到多路有用信号和噪声混合成一路的单通道混合信号的问题。单通道混合信号在生物医学、水声信号、语音信号、通信信号等领域也广泛存在。由于无线电信号在我们日常生活中起着很重要的作用，因此面对复杂的电磁环境，研究单通道无线电混合信号的分离具有广泛的应用价值和重要意义。

目前对于单通道混合信号分离已有一些方法，如朱航等(朱航,张淑宁,赵惠昌.基于改进自适应分解法的单通道雷达引信混合信号分离[J].物理学报,2014,63(5)058401:1-14.)利用自适应算法对语音信号和伪随机信号进行分离，并利用相关增强提高分离精度；朱行涛等(朱行涛,刘郁林,何为,等.基于变换域滤波的直扩通信单通道混合信号分离抗干扰方法[J].系统工程与电子技术,2016,38(10):2405-2412.)利用匹配追踪对稳态调制成分和冲击调制成分进行分离，有效诊断出齿轮箱故障；张淑宁等(张淑宁,赵惠昌,熊刚,等.基于粒子滤波的单通道正弦调频混合信号分离与参数估计[J].物理学报,2014,63(15)158401:1-9.)利用粒子滤波对信噪比为20dB的正弦调频信号进行分离，完全分离出混合信号中的正弦调频信号；刘威等(刘威,胡爱群.基于相关增强的单信道信号分离[J].东南大学学报(自然科学版),2010,40(3):464-470.)利用变换域滤波对混有干扰和噪声的直扩信号进行分离，提升了通信信号的抗干扰能力。然而，在单通道混合信号个数大于2以及混合后信号信噪比较低的情况下，上述方法并不总能达到良好的分离效果。经验模态分解(EMD)作为一种自适应信号处理方法，可以较好地对多个混合信号以及较低信噪比情况下的信号进行分离。EMD是递归式模式分解方法，会根据混合信号不同尺度的波动或趋势逐级分解，分解成频率从高频到低频分布的一系列本征模态函数分量，但EMD分解过程中容易产生模态混叠现象，针对这种现象，提出了总体平均经验模态分解法(EEMD)，虽在一定程度上消除了模态混叠的影响，但添加的白噪声不能被完全中和，不具有完备性。作为自适应信号处理的一种新方法——变分模态分解(VMD)，也可以将低信噪比条件下的信号分离开来，它是完全非递归式模式分解方法，会将混合信号分解成一系列分量信号，并且每一个模态的绝大部分都紧紧围绕在中心频率的周围，但VMD算法的信号处理结果同时受到分量个数和惩罚因子这两个参数的影响。

发明内容

针对现有技术的不足，本发明拟解决的技术问题是，提供一种基于变分模态分解的单通道无线电混合信号分离方法。该方法通过改进的量子粒子群优化算法对变分模态分解算法的分量个数和惩罚因子两个参数同时进行优化，确定变分模态分解算法中的分量个数和惩罚因子的值。之后首次利用参数优化的变分模态分解算法对单通道无线电混合信号进行分离，由于待处理的单通道无线电混合信号的信噪比较低，分离后所得各个无线电信号信噪比有所提高，但仍有大量噪声的存在，因此将分离后的各个无线电信号依次使用参数优化的变分模态分解算法进行二次分离，二次分离后所得各个无线电信号比首次分离后所得相应无线电信号可信度高。

本发明解决所述技术问题的技术方案是，提供一种基于变分模态分解的单通道无线电混合信号分离方法，其特征在于该分离方法包括以下步骤：

第一步，引入单通道无线电混合信号；

第二步，使用改进的量子粒子群优化算法确定此单通道无线电混合信号的变分模态分解算法的分量个数和惩罚因子两个参数；

2-1.改进的量子粒子群优化算法的初始化：设置改进的量子粒子群优化算法的最大迭代次数T_max、种群数量M以及粒子的位置向量的寻优范围，其中粒子的位置向量为两个维度，一个维度为分量个数，另一个维度为惩罚因子，并在其限制范围内随机初始化粒子的位置向量，其表达式如式(1)：

x_i＝(x_i1,x_i2)(i＝1,2,...,M) (1)

x_i为粒子的位置向量，x_i1为分量个数，x_i2为惩罚因子；

2-2.最优适应度值的初始化：选取信噪比f_snr作为适应度函数，根据信噪比的定义，其表达式如式(2)：

s(t)为混合之前的某个源信号，u(t)为分离后与之对应的信号；

由于源信号未知，可以将分离后的各个时域信号变换成频域信号，并利用信号和噪声的频率特性，再通过积分的方法将信号功率和噪声功率计算出来，即可得到信噪比，其表达式如式(3)：

P_s为信号功率，P_n为噪声功率；

然后根据式(3)计算2-1中各个粒子对应的适应度值，将各个粒子自身的适应度值作为单个粒子最优适应度值，记为p_best(i)(i＝1,2,...,M)，将p_best(i)中的最大值作为全局最优适应度值，记为g_best；

2-3.更新最优适应度值；首先要更新粒子的位置，然后重新计算出各个粒子的适应度值，若得到当前代的单个粒子最优适应度值p_best(i)优于上一代粒子的单个粒子最优适应度值p_best(i)，则更新对应的单个粒子最优适应度值；若得到当前代的全局最优适应度值g_best优于上一代粒子的全局最优适应度值g_best，则更新对应的全局最优适应度值；各个粒子位置更新的表达式如式(4)：

x_i(t_a+1)＝P_i(t_a)±β|m_best(t_a)-x_i(t_a)|ln(1/u) (4)

t_a为当前迭代次数；P_i(t_a)为第t_a次迭代中第i个粒子的最优适应度值p_best(i)；u为[0,1]范围内变化的随机数；β为收缩-扩张系数；m_best(t_a)为平均最好位置，即第t_a次迭代时所有粒子个体最好位置的均值，其表达式如式(5)：

M表示种群数量；

量子粒子群优化算法的改进是对平均最好位置m_best和收缩-扩张系数β的改进；平均最好位置m_best的改进是在平均最好位置m_best中加入sobol随机序列(sobolrand())，则平均最好位置的表达式如式(7)：

m_best(t_a)＝m_best(t_a)+sobolrand() (7)

为了实时调节收缩-扩张系数，引入进化速度因子v和聚集度因子d；

进化速度因子的表达式如式(8)：

f_snr(x(t_a))和f_snr(x(t_a-1))分别是量子粒子群中的所有粒子当前代和前一代的全局最优适应度值，并且0＜v≤1，当v值较小时，说明全局最优适应度值变化较大，此时的进化速度较快，粒子距最优位置较远；当v值较大时，说明全局最优适应度值变化较慢，此时的进化速度较慢，粒子距最优位置较近；当v＝1时，此时的进化停止，说明算法搜索到了最优位置；

聚集度因子的表达式如式(9)：

0＜d≤0，当d值较小时，说明粒子的聚集程度较低，此时量子粒子群中所有粒子的群体多样性较高；当d值较大时，说明粒子的聚集程度较高，此时量子粒子群中所有粒子的群体多样性较低；当d＝1时，说明粒子的聚集程度非常高，此时的所有粒子具有同一性，f_snra是当前代所有粒子中个体最优适应度值的均值，表达式如式(10)：

为了使改进后的量子粒子群优化算法可以自适应的调节粒子位置，因此对于距最优位置较远的粒子，即v较小时，可以增大收缩-扩张系数β达到快速寻优的效果；对于距最优位置较近的粒子，即v较大时，可以减小收缩-扩张系数β达到精确寻优的效果；对于聚集程度较低的粒子，即d较小时，粒子不容易陷入局部最优，可以减小收缩-扩张系数β达到减小粒子搜索空间的效果；对于聚集程度较高的粒子，即d较大时，粒子容易陷入局部最优，可以增大收缩-扩张系数β达到增大粒子搜索空间的效果，因此改进后的量子粒子群优化算法的收缩-扩张系数β的表达式如式(11)：

β＝β₀-δ_vv+δ_dd (11)

β₀为β的初始值，一般β₀的取值为1；δ_v和δ_d分别为调整进化速度因子的系数和调整聚集度因子的系数，一般取δ_v＝0.5，δ_d＝0.2；

2-4.迭代结束，得到全局最优适应度值，根据全局最优适应度值，得到相应的粒子的最优位置，即为所需分量个数和惩罚因子两个参数进行优化后的参数组合；

第三步，将第二步确定的参数组合应用到变分模态分解算法中，然后用参数优化后的变分模态分解算法将此单通道无线电混合信号进行分离；

变分模态分解算法的表达式如式(12)：

其中：u_k(t)是重新定义的调幅-调频信号，也就是重新定义后的IMF分量；A_k(t)为u_k(t)的瞬时幅值，为u_k(t)的瞬时相位；设ω_k(t)为u_k(t)的瞬时角频率，则

变分模态分解算法处理信号的方式是通过构造并求解约束变分问题将单通道无线电混合信号自适应的分解成指定个数并紧密围绕在中心频率ω_k附近的IMF分量信号；假设将此待处理信号分解成K个IMF分量，其中K为第二步中确定的分量个数，那么相应的约束变分问题构造模型如式(13)：

其中{u_k}＝{u₁,...,u_K}为分解得到的K个IMF分量，{ω_k}＝{ω₁,...,ω_K}为各分量的中心频率，为u_k(t)通过Hilbert变换得到其解析信号；

为求取上述约束变分问题的最优解，需将其转换成非约束变分问题，因此引入式(14)形式的增广Lagrange函数，即：

其中α为二次项的惩罚因子，其值为在第二步中确定的惩罚因子；λ(t)为Lagrange乘子；x(t)为待处理信号；

然后利用乘子交替方向算法求取式(14)的鞍点，即为式(13)约束变分问题的最优解，从而将信号x(t)分解成K个IMF分量；

参数优化后的变分模态分解算法对此单通道无线电混合信号进行分离的步骤如下：

3-1.初始化λ¹，n为0；其中为第k个初始化后的IMF分量；为第k个初始化后的中心频率；λ¹是初始化后的Lagrange乘子；n为循环的次数；

3-2.令n＝n+1，执行整个循环；

3-3.令k＝1，更新u_k；u_k的更新问题为的取值问题，则的取值问题可表述为式(15)：

令和利用Parseval/Plancherel、傅里叶等距变换，将式(15)变换到频域，表达式如式(16)：

令式(16)中第一项的ω＝ω-ω_k，则有式(17)

式(17)非负频率区间积分的形式为式(18)：

然后通过二次优化得到的解为式(19)：

3-4.令k＝k+1，并重复步骤3-3，直到k＝K时停止；

3-5.继续令k＝1，更新ω_k；ω_k的更新问题为的取值问题，则的取值问题可表述为式(20)：

同步骤3-3的计算原理，优化得：

求解得：

3-6.令k＝k+1，并重复步骤3-5，直到k＝K时停止；

3-7.根据表达式更新λ，其中τ表示更新参数；

3-8.重复步骤3-2～3-7，直至满足终止条件时，结束整个循环，输出的结果为K个IMF分量信号；

第四步，由于混合的无线电信号个数小于等于分量个数K，因此根据分离出的IMF分量信号的时域图和频域图以及IMF分量信号的信噪比筛选出各个无线电信号；

第五步，对筛选出的各个无线电信号分别按照第二步到第四步进行二次分离，得到各个无线电信号。

与现有技术相比，本发明有益效果在于：

(1)本发明方法在传统量子粒子群优化算法的基础上，在平均最好位置上进行sobol序列扰动，并引入进化速度因子和粒子聚集度因子，进而实时调节收缩-扩张系数，从而高效的搜索最优解，有效地避免了粒子陷入局部最优解。

(2)本发明方法采用改进的量子粒子群优化算法对变分模态分解中的分量个数和惩罚因子进行组合优化，避免了人为或者根据经验确定的参数组合，存在很大的偶然性和随机性，使分析效果难以保障。

(3)本发明方法选取信噪比作为改进量子粒子群优化算法的适应度函数，获得一种直观的确定变分模态分解算法参数的方法。

(4)本发明方法将参数优化后的变分模态分解应用于低信噪比的单通道无线电混合信号的分离，分离出的各个无线电信号在时域、频域以及信噪比等方面取得了较好的分离效果。

(5)本发明方法两次采用参数优化的变分模态分解算法对单通道无线电混合信号进行分离，二次分离后的各个无线电信号在时域、频域以及信噪比等方面均比首次分离后所得相应无线电信号有很大程度的提高。

附图说明

图1为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1单通道无线电混合信号的时域图。

图2为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1单通道无线电混合信号的频域图。

图3为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1首次分离后无线电信号1的时域图。

图4为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1首次分离后无线电信号1的频域图。

图5为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1首次分离后无线电信号2的时域图。

图6为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1首次分离后无线电信号2的频域图。

图7为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1首次分离后无线电信号3的时域图。

图8为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1首次分离后无线电信号3的频域图。

图9为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1二次分离后无线电信号1的时域图。

图10为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1二次分离后无线电信号1的频域图。

图11为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1二次分离后无线电信号2的时域图。

图12为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1二次分离后无线电信号2的频域图。

图13为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1二次分离后无线电信号3的时域图。

图14为本发明基于变分模态分解的单通道无线电混合信号分离方法实施例1二次分离后无线电信号3的频域图。

具体实施方式

下面结合附图和本实施例对本发明进一步说明，但并不以此作为对本申请权利要求保护范围的限定。

本发明提供了一种基于变分模态分解的单通道无线电混合信号分离方法(简称分离方法)，其特征在于该分离方法包括以下步骤：

第一步，引入单通道无线电混合信号；

第二步，使用改进的量子粒子群优化算法确定此单通道无线电混合信号的变分模态分解算法的分量个数和惩罚因子两个参数；

由于变分模态分解算法同时受到分量个数和惩罚因子的影响，且分量个数越多，惩罚因子越大，分量信号的带宽越窄；分量信号的带宽过宽或过窄，均会影响到信号处理的结果，而且人为或者根据经验确定的分量个数和惩罚因子，会使分析效果难以保障，并且存在很大的偶然性和随机性，因此选定合适的参数组合，是利用变分模态分解算法分离单通道混合信号的关键。

由于传统量子粒子群优化算法的收缩-扩张系数是线性递减的，这种收缩-扩张系数的确定方法独立于算法本身较为固定，但实际搜索过程往往是非线性的且高度复杂的，应能反映粒子运行过程中的实际状况而自适应变化，而且使用传统的量子粒子群优化算法也容易陷入局部最优的缺陷。为了解决上述问题，本发明提出改进的量子粒子群优化算法，对变分模态分解算法中的分量个数和惩罚因子两个参数同时进行寻优；改进的量子粒子群优化算法的寻优步骤如下：

2-1.改进的量子粒子群优化算法的初始化：设置改进的量子粒子群优化算法的最大迭代次数T_max、种群数量M以及粒子的位置向量的寻优范围，其中粒子的位置向量为两个维度，一个维度为分量个数，另一个维度为惩罚因子，并在其限制范围内随机初始化粒子的位置向量，其表达式如式(1)：

x_i＝(x_i1,x_i2)(i＝1,2,...,M) (1)

x_i为粒子的位置向量，x_i1为分量个数，x_i2为惩罚因子；

2-2.最优适应度值的初始化：最优适应度值包括单个粒子最优适应度值和全局最优适应度值，其中所有粒子的最优适应度值中的最大值也就是全局最优适应度值；由于信噪比f_snr可以反映对有用信号的筛选和提纯能力，所以选取信噪比作为适应度函数，根据信噪比的定义，其表达式如式(2)：

s(t)为混合之前的某个源信号，u(t)为分离后与之对应的信号；

在实际应用中，由于源信号未知，可以将分离后的各个时域信号变换成频域信号，并利用信号和噪声的频率特性，再通过积分的方法将信号功率和噪声功率计算出来，即可得到信噪比，其表达式如式(3)：

P_s为信号功率，P_n为噪声功率；

然后根据式(3)计算2-1中各个粒子对应的适应度值，将各个粒子自身的适应度值作为单个粒子最优适应度值，记为p_best(i)(i＝1,2,...,M)，将p_best(i)中的最大值作为全局最优适应度值，记为g_best；

2-3.更新最优适应度值；首先要更新粒子的位置，然后重新计算出各个粒子的适应度值，若得到当前代的单个粒子最优适应度值p_best(i)优于上一代粒子的单个粒子最优适应度值p_best(i)，则更新对应的单个粒子最优适应度值；若得到当前代的全局最优适应度值g_best优于上一代粒子的全局最优适应度值g_best，则更新对应的全局最优适应度值；各个粒子位置更新的表达式如式(4)：

x_i(t_a+1)＝P_i(t_a)±β|m_best(t_a)-x_i(t_a)|ln(1/u) (4)

t_a为当前迭代次数；P_i(t_a)为第t_a次迭代中第i个粒子的最优适应度值p_best(i)；u为[0,1]范围内变化的随机数；β为收缩-扩张系数；m_best(t_a)为平均最好位置，即第t_a次迭代时所有粒子个体最好位置的均值，其表达式如式(5)：

M表示种群数量；

传统的量子粒子群优化算法中，收缩-扩张系数β一般按照式(6)进行变化：

量子粒子群优化算法的改进是对平均最好位置m_best和收缩-扩张系数β的改进，改进后的量子粒子群优化算法提高了搜索最优解的效率，有效地避免了粒子陷入局部最优解；平均最好位置m_best的改进是在平均最好位置m_best中加入sobol随机序列(sobolrand())，在粒子群的不断进化过程中，对平均最好位置进行随机扰动，可使粒子群一直具有活动性，由于sobol随机序列是准随机序列中的一种，比均匀分布的随机序列均匀性好，耗时少，因此选择sobol随机序列对平均最好位置进行随机扰动，则平均最好位置的表达式如式(7)：

m_best(t_a)＝m_best(t_a)+sobolrand() (7)

由于传统的量子粒子群优化算法中控制收缩-扩张系数β的方法是基于简单的线性方法来实现的，随着迭代次数的增加，对量子粒子群中的所有粒子进行统一的参数调整，就不能根据粒子的实际情况而自适应的变化粒子的位置，有着易陷入局部最优的缺陷，因此为了实时调节收缩-扩张系数，引入进化速度因子v和聚集度因子d；

进化速度因子的表达式如式(8)：

f_snr(x(t_a))和f_snr(x(t_a-1))分别是量子粒子群中的所有粒子当前代和前一代的全局最优适应度值，并且0＜v≤1，当v值较小时，说明全局最优适应度值变化较大，此时的进化速度较快，粒子距最优位置较远；当v值较大时，说明全局最优适应度值变化较慢，此时的进化速度较慢，粒子距最优位置较近；当v＝1时，此时的进化停止，说明算法搜索到了最优位置；

聚集度因子的表达式如式(9)：

0＜d≤0，当d值较小时，说明粒子的聚集程度较低，此时量子粒子群中所有粒子的群体多样性较高；当d值较大时，说明粒子的聚集程度较高，此时量子粒子群中所有粒子的群体多样性较低；当d＝1时，说明粒子的聚集程度非常高，此时的所有粒子具有同一性，f_snra是当前代所有粒子中个体最优适应度值的均值，表达式如式(10)：

为了使改进后的量子粒子群优化算法可以自适应的调节粒子位置，因此对于距最优位置较远的粒子，即v较小时，可以增大收缩-扩张系数β达到快速寻优的效果；对于距最优位置较近的粒子，即v较大时，可以减小收缩-扩张系数β达到精确寻优的效果；对于聚集程度较低的粒子，即d较小时，粒子不容易陷入局部最优，可以减小收缩-扩张系数β达到减小粒子搜索空间的效果；对于聚集程度较高的粒子，即d较大时，粒子容易陷入局部最优，可以增大收缩-扩张系数β达到增大粒子搜索空间的效果，因此改进后的量子粒子群优化算法的收缩-扩张系数β的表达式如式(11)：

β＝β₀-δ_vv+δ_dd (11)

β₀为β的初始值，一般β₀的取值为1；δ_v和δ_d分别为调整进化速度因子的系数和调整聚集度因子的系数，一般取δ_v＝0.5，δ_d＝0.2；

2-4.迭代结束，得到全局最优适应度值，根据全局最优适应度值，得到相应的粒子的最优位置，即为所需分量个数和惩罚因子两个参数进行优化后的参数组合；

第三步，将第二步确定的参数组合应用到变分模态分解算法中，然后用参数优化后的变分模态分解算法将此单通道无线电混合信号进行分离；

变分模态分解算法是由Dragomiretskiy等人于2014年提出的，并将本征模态函数分量(IMF分量)重新定义为调幅-调频信号，其表达式如式(12)：

其中：u_k(t)是重新定义的调幅-调频信号，也就是重新定义后的IMF分量；A_k(t)为u_k(t)的瞬时幅值，为u_k(t)的瞬时相位。设ω_k(t)为u_k(t)的瞬时角频率，则

变分模态分解算法处理信号的方式是通过构造并求解约束变分问题将单通道无线电混合信号自适应的分解成指定个数并紧密围绕在中心频率ω_k附近的IMF分量信号。假设将此待处理信号分解成K个IMF分量，其中K为第二步中确定的分量个数，那么相应的约束变分问题构造模型如式(13)：

其中{u_k}＝{u₁,...,u_K}为分解得到的K个IMF分量，{ω_k}＝{ω₁,...,ω_K}为各分量的中心频率，为u_k(t)通过Hilbert变换得到其解析信号；

为求取上述约束变分问题的最优解，需将其转换成非约束变分问题，因此引入式(14)形式的增广Lagrange函数，即：

其中α为二次项的惩罚因子，其值为在第二步中确定的惩罚因子；λ(t)为Lagrange乘子；x(t)为待处理信号；

然后利用乘子交替方向算法求取式(14)的鞍点，即为式(13)约束变分问题的最优解，从而将信号x(t)分解成K个IMF分量；

参数优化后的变分模态分解算法对此单通道无线电混合信号进行分离的步骤如下：

3-1.初始化λ¹，n为0；其中为第k个初始化后的IMF分量；为第k个初始化后的中心频率；λ¹是初始化后的Lagrange乘子；n为循环的次数，当满足3-8中的终止条件，循环结束；

3-2.令n＝n+1，执行整个循环；

3-3.令k＝1，更新u_k；u_k的更新问题可以理解为的取值问题，则的取值问题可表述为式(15)：

为了简单起见，认为和然后利用Parseval/Plancherel、傅里叶等距变换，将式(15)变换到频域，表达式如式(16)：

式(16)中的字母是式(15)中时域变换到频域后的表达式；

令式(16)中第一项的ω＝ω-ω_k，则有式(17)

式(17)非负频率区间积分的形式为式(18)：

然后通过二次优化得到的解为式(19)：

3-4.令k＝k+1，并重复步骤3-3，直到k＝K时停止；

3-5.继续令k＝1，更新ω_k；ω_k的更新问题可以理解为的取值问题，则的取值问题可表述为式(20)：

同步骤3-3的计算原理，优化得：

求解得：

3-6.令k＝k+1，并重复步骤3-5，直到k＝K时停止；

3-7.根据表达式更新λ，其中τ表示更新参数；

3-8.重复步骤3-2～3-7，直至满足终止条件时，结束整个循环，输出的结果为K个IMF分量信号；

第四步，由于混合的无线电信号个数小于等于分量个数K，因此根据分离出的IMF分量信号的时域图和频域图以及IMF分量信号的信噪比筛选出各个无线电信号；

第五步，对筛选出的各个无线电信号分别按照第二步到第四步进行二次分离，得到在时域、频域以及信噪比方面比首次分离后所得相应无线电信号有很大程度提高的各个无线电信号。

实施例1

本实施例基于参数优化变分模态分解的单通道无线电混合信号分离方法，该方法的具体步骤如下：

第一步，引入单通道无线电混合信号；

待分离的单通道无线电信号由三个无线电信号和高斯白噪声混合而成，无线电信号1频率为30Hz；无线电信号2频率为50Hz；无线电信号3频率为80Hz；混合后的信噪比为-13.2192dB，图1和图2分别为单通道无线电混合信号的时域图和频域图。

第二步，使用改进的量子粒子群优化算法确定此单通道无线电混合信号的变分模态分解算法的分量个数和惩罚因子两个参数；

2-1.改进的量子粒子群优化算法的初始化，设置改进的量子粒子群优化算法的最大迭代次数T_max为30次、种群数量M为10个，粒子位置向量的维度为2，一个维度为分量个数，寻优范围是[3,20]，另一个维度为惩罚因子，寻优范围是[200,5000]；

2-2.最优适应度值的初始化；

2-3.更新最优适应度值；

2-4.迭代结束，得到全局最优适应度值为2.8156dB，然后根据全局最优适应度值，得到相应的粒子的最优位置[20,5000]，即所需分量个数为20个，和惩罚因子为5000；

第三步，将第二步确定的参数组合[20,5000]应用到变分模态分解算法中，然后用参数优化后的变分模态分解算法将此单通道无线电混合信号进行分离；

第四步，由于混合的无线电信号个数小于等于分量个数，因此根据分离出的IMF分量信号的时域图和频域图，以及信噪比筛选出各个无线电信号；分离出的IMF2为无线电信号1，其信噪比为3.8352dB；分离出的IMF3为无线电信号2，其信噪比为1.8722dB；分离出的IMF4为无线电信号3，其信噪比为2.7393dB；图3和图4分别为首次分离后筛选出无线电信号1时域图和频域图，图5和图6分别为首次分离后筛选出无线电信号2时域图和频域图，图7和图8分别为首次分离后筛选出无线电信号3时域图和频域图。

第五步，开始二次分离，对筛选出的IMF2，IMF3，IMF4分别使用改进的量子粒子群优化算法确定处理IMF2，IMF3，IMF4的变分模态分解算法的分量个数和惩罚因子两个参数；

改进的量子粒子群优化算法的寻优步骤如下：

5-1.改进的量子粒子群优化算法的初始化，设置改进的量子粒子群优化算法的最大迭代次数T_max为30次、种群数量M为10个，粒子位置向量的维度为2，一个维度为分量个数，寻优范围是[2,15]，另一个维度为惩罚因子，寻优范围是[200,5000]；

5-2.最优适应度值的初始化；

5-3.更新最优适应度值；

5-4.迭代结束，得到全局最优适应度值分别为8.6116dB，9.6724dB，9.1604dB，根据全局最优适应度值，得到相应的粒子的最优位置分别为[5,3463]，[8,3156]，[15,1004]，即IMF2所需分量个数为5个和惩罚因子为3463；IMF3所需分量个数为8个和惩罚因子为3156；IMF4所需分量个数为15个和惩罚因子为1004；

第六步，将第五步确定的参数组合分别应用到变分模态分解算法中，然后用参数优化后的变分模态分解算法分别将IMF2，IMF3，IMF4进行分离；

第七步，由于混合的无线电信号个数小于等于分量个数，因此根据分离出的IMF分量信号的时域图和频域图，以及信噪比筛选出各个无线电信号，分离出的无线电信号1的信噪比为8.6116dB，分离出的无线电信号2的信噪比为9.6724dB，分离出的无线电信号3的信噪比为9.1604dB，图9和图10分别为二次分离后筛选出无线电信号1时域图和频域图，图11和图12分别为二次分离后筛选出无线电信号2时域图和频域图，图13和图14分别为二次分离后筛选出无线电信号3时域图和频域图。

表1首次和二次分离后的3个无线电信号信噪比对比

根据表可以得出，二次分离后得到各个无线电信号在时域、频域以及信噪比方面均比首次分离后所得相应无线电信号有很大程度的提高。

上述步骤均采用MATLAB软件实现。

上述本实施例中所用的MATLAB软件是本技术领域的技术人员所熟知的。

本发明未述及之处适用于现有技术。

Claims (1)

1.一种基于变分模态分解的单通道无线电混合信号分离方法，其特征在于该分离方法包括以下步骤：

第一步，引入单通道无线电混合信号；

第二步，使用改进的量子粒子群优化算法确定此单通道无线电混合信号的变分模态分解算法的分量个数和惩罚因子两个参数；

2-1.改进的量子粒子群优化算法的初始化：设置改进的量子粒子群优化算法的最大迭代次数T_max、种群数量M以及粒子的位置向量的寻优范围，其中粒子的位置向量为两个维度，一个维度为分量个数，另一个维度为惩罚因子，并在其限制范围内随机初始化粒子的位置向量，其表达式如式(1)：

x_i＝(x_i1,x_i2)(i＝1,2,...,M) (1)

x_i为粒子的位置向量，x_i1为分量个数，x_i2为惩罚因子；

2-2.最优适应度值的初始化：选取信噪比f_snr作为适应度函数，根据信噪比的定义，其表达式如式(2)：

s(t)为混合之前的某个源信号，u(t)为分离后与之对应的信号；

由于源信号未知，可以将分离后的各个时域信号变换成频域信号，并利用信号和噪声的频率特性，再通过积分的方法将信号功率和噪声功率计算出来，即可得到信噪比，其表达式如式(3)：

P_s为信号功率，P_n为噪声功率；

然后根据式(3)计算2-1中各个粒子对应的适应度值，将各个粒子自身的适应度值作为单个粒子最优适应度值，记为p_best(i)(i＝1,2,...,M)，将p_best(i)中的最大值作为全局最优适应度值，记为g_best；

2-3.更新最优适应度值；首先要更新粒子的位置，然后重新计算出各个粒子的适应度值，若得到当前代的单个粒子最优适应度值p_best(i)优于上一代粒子的单个粒子最优适应度值p_best(i)，则更新对应的单个粒子最优适应度值；若得到当前代的全局最优适应度值g_best优于上一代粒子的全局最优适应度值g_best，则更新对应的全局最优适应度值；各个粒子位置更新的表达式如式(4)：

x_i(t_a+1)＝P_i(t_a)±β|m_best(t_a)-x_i(t_a)|ln(1/u) (4)

t_a为当前迭代次数；P_i(t_a)为第t_a次迭代中第i个粒子的最优适应度值p_best(i)；u为[0,1]范围内变化的随机数；β为收缩-扩张系数；m_best(t_a)为平均最好位置，即第t_a次迭代时所有粒子个体最好位置的均值，其表达式如式(5)：

M表示种群数量；

量子粒子群优化算法的改进是对平均最好位置m_best和收缩-扩张系数β的改进；平均最好位置m_best的改进是在平均最好位置m_best中加入sobol随机序列(sobolrand())，则平均最好位置的表达式如式(7)：

m_best(t_a)＝m_best(t_a)+sobolrand() (7)

为了实时调节收缩-扩张系数，引入进化速度因子v和聚集度因子d；

进化速度因子的表达式如式(8)：

f_snr(x(t_a))和f_snr(x(t_a-1))分别是量子粒子群中的所有粒子当前代和前一代的全局最优适应度值，并且0＜v≤1，当v值较小时，说明全局最优适应度值变化较大，此时的进化速度较快，粒子距最优位置较远；当v值较大时，说明全局最优适应度值变化较慢，此时的进化速度较慢，粒子距最优位置较近；当v＝1时，此时的进化停止，说明算法搜索到了最优位置；

聚集度因子的表达式如式(9)：

0＜d≤0，当d值较小时，说明粒子的聚集程度较低，此时量子粒子群中所有粒子的群体多样性较高；当d值较大时，说明粒子的聚集程度较高，此时量子粒子群中所有粒子的群体多样性较低；当d＝1时，说明粒子的聚集程度非常高，此时的所有粒子具有同一性，f_snra是当前代所有粒子中个体最优适应度值的均值，表达式如式(10)：

为了使改进后的量子粒子群优化算法可以自适应的调节粒子位置，因此对于距最优位置较远的粒子，即v较小时，可以增大收缩-扩张系数β达到快速寻优的效果；对于距最优位置较近的粒子，即v较大时，可以减小收缩-扩张系数β达到精确寻优的效果；对于聚集程度较低的粒子，即d较小时，粒子不容易陷入局部最优，可以减小收缩-扩张系数β达到减小粒子搜索空间的效果；对于聚集程度较高的粒子，即d较大时，粒子容易陷入局部最优，可以增大收缩-扩张系数β达到增大粒子搜索空间的效果，因此改进后的量子粒子群优化算法的收缩-扩张系数β的表达式如式(11)：

β＝β₀-δ_vv+δ_dd (11)

β₀为β的初始值，一般β₀的取值为1；δ_v和δ_d分别为调整进化速度因子的系数和调整聚集度因子的系数，一般取δ_v＝0.5，δ_d＝0.2；

2-4.迭代结束，得到全局最优适应度值，根据全局最优适应度值，得到相应的粒子的最优位置，即为所需分量个数和惩罚因子两个参数进行优化后的参数组合；

第三步，将第二步确定的参数组合应用到变分模态分解算法中，然后用参数优化后的变分模态分解算法将此单通道无线电混合信号进行分离；

变分模态分解算法的表达式如式(12)：

其中：u_k(t)是重新定义的调幅-调频信号，也就是重新定义后的IMF分量；A_k(t)为u_k(t)的瞬时幅值，为u_k(t)的瞬时相位；设ω_k(t)为u_k(t)的瞬时角频率，则

变分模态分解算法处理信号的方式是通过构造并求解约束变分问题将单通道无线电混合信号自适应的分解成指定个数并紧密围绕在中心频率ω_k附近的IMF分量信号；假设将此待处理信号分解成K个IMF分量，其中K为第二步中确定的分量个数，那么相应的约束变分问题构造模型如式(13)：

其中{u_k}＝{u₁,...,u_K}为分解得到的K个IMF分量，{ω_k}＝{ω₁,...,ω_K}为各分量的中心频率，为u_k(t)通过Hilbert变换得到其解析信号；

为求取上述约束变分问题的最优解，需将其转换成非约束变分问题，因此引入式(14)形式的增广Lagrange函数，即：

其中α为二次项的惩罚因子，其值为在第二步中确定的惩罚因子；λ(t)为Lagrange乘子；x(t)为待处理信号；

然后利用乘子交替方向算法求取式(14)的鞍点，即为式(13)约束变分问题的最优解，从而将信号x(t)分解成K个IMF分量；

参数优化后的变分模态分解算法对此单通道无线电混合信号进行分离的步骤如下：

3-1.初始化λ¹，n为0；其中为第k个初始化后的IMF分量；为第k个初始化后的中心频率；λ¹是初始化后的Lagrange乘子；n为循环的次数；

3-2.令n＝n+1，执行整个循环；

3-3.令k＝1，更新u_k；u_k的更新问题为的取值问题，则的取值问题可表述为式(15)：

令和利用Parseval/Plancherel、傅里叶等距变换，将式(15)变换到频域，表达式如式(16)：

令式(16)中第一项的ω＝ω-ω_k，则有式(17)

式(17)非负频率区间积分的形式为式(18)：

然后通过二次优化得到的解为式(19)：

3-4.令k＝k+1，并重复步骤3-3，直到k＝K时停止；

3-5.继续令k＝1，更新ω_k；ω_k的更新问题为的取值问题，则的取值问题可表述为式(20)：

同步骤3-3的计算原理，优化得：

求解得：

3-6.令k＝k+1，并重复步骤3-5，直到k＝K时停止；

3-7.根据表达式更新λ，其中τ表示更新参数；

3-8.重复步骤3-2～3-7，直至满足终止条件时，结束整个循环，输出的结果为K个IMF分量信号；

第四步，由于混合的无线电信号个数小于等于分量个数K，因此根据分离出的IMF分量信号的时域图和频域图以及IMF分量信号的信噪比筛选出各个无线电信号；

第五步，对筛选出的各个无线电信号分别按照第二步到第四步进行二次分离，得到各个无线电信号。