基于去噪复数FastICA 和稀疏重构的相干信号欠定DOA 估计

侯进，李昀喆，李天宇

（1.西南交通大学信息科学与技术学院，四川 成都 611756；2.西南交通大学综合交通大数据应用技术国家工程实验室，四川 成都 611756；3.西南交通大学唐山研究生院，河北 唐山 063000）

1 引言

实际的通信环境中，多径效应导致存在大量相干信号，使入射信号的数量超过传感器的数量的现象普遍存在。在这样的复杂环境中，现有的到达角（DOA,direction of arrival）估计方法将失败。DOA估计还有一个重要的目标就是在大量噪声的情况下，定位间隔较小的信号源。许多经典的DOA 估计技术通过利用少量源的存在来实现超分辨率。例如，Krim 等[1]研究的多重信号分类（MUSIC,multiple signal classification）算法的关键部分就是低维子空间的假设。因此，需要针对以上情况寻找一种新的DOA 估计方法。

目前，DOA 估计方法可分为三类：基于子空间的方法、基于独立分量分析（ICA,independent component analysis）的方法、基于稀疏重构的方法。当传感器和时间样本足够时，除MUSIC 算法，Esfandiari 等[2]和Zhuang 等[3]研究了其他基于子空间的方法，该类方法利用信号和噪声子空间之间的正交性来估计信号源的DOA。但是，当入射信号的数量大于传感器的数量时，这些算法因为噪声子空间的消失都会失去最佳性能。

当入射信号的数量不大于传感器的数量时，基于ICA的方法可以实现信号的盲分离和DOA估计，其中分离再测向有助于多信源的 DOA 估计。Zamani 等[4]和Johnson[5]通过ICA 与波束合成结合，提高算法的收敛速度，但这些方法只能解决独立信号的DOA 估计问题。对于相干信号，基于ICA 的方法的性能将严重下降。基于稀疏重构的方法可以根据不同的稀疏重构算法分为三类：Trinh-Hoang等[6]探究的基于匹配追踪（MP,matching pursuit）的方法、Malioutov 等[7]和Soubies 等[8]探究的基于lp范数的方法、Hu 等[9]和Huang 等[10]探究的基于稀疏贝叶斯学习（SBL,sparse bayesian learning）的方法。与基于子空间的方法和基于ICA 的方法相比，基于稀疏重构的方法在低信噪比（SNR,signal-to-noise ratio）、小样本和信号相干的情况下可以获得更好的性能。但是需要入射信号的数量小于传感器的数量。对于欠定情况，仍然无法解决相干信号DOA 估计问题。

针对现有算法的不足，本文提出了一种基于ICA 和稀疏重构的DOA 估计方法。ICA 算法中的快速独立主成分分析（FastICA,fast independent component analysis）算法具有收敛速度快、易使用等优点，故本文使用FastICA 算法进行信号盲分离。首先，通过FastICA 算法估计阵列流形；然后，将阵列流形中的导向矢量视为阵列输出的单个时间样本，并通过稀疏重构来估计分离后的独立信源的DOA，独立信源可能包含相干信号和非相干信号。理论分析和仿真与实测实验结果表明，该方法可以解决欠定情况下相干信号的DOA 估计问题，并能进行多个非相干圆或非圆信号的DOA 估计。

2 系统模型

假设有M个阵元的均匀圆阵（UCA,uniform circular array），如图1 所示。该圆形阵列接收N个入射复数独立信源，包含NI个非相干复信号源和K组相干复信号源，则独立信源总数为N=NI+K，这些信号源位于不同方向的阵列远场中。其混合模型可以写成n

θ为第n个信号的波达方向角，λ为信号的波长，r为均匀圆阵的半径；n(t)为均值为0 的加性白高斯噪声。

当信号为相干信号时，假设第k个入射源信号sk(t)包含Lk个相干信号，这些相干信源从θlk方向入射到阵列上，其中lk=1,2…Lk。假设在θlk中有2 个相关信号sik(t)和sjk(t)，则两信源的相关系数表示为

由Schwartz 不等式可知|μlk|≤ 1，故信号间的相关性定义为

由以上定义可知，当信源相干时，相干信源之间只相差一复常数。阵列相干信号的数据模型应调整为

若信源中存在K组相干信号，相干信源总数为，则入射信源总数即L=NI+Nk。

在实际环境中，由于相干信号存在，很可能导致入射信号的数量大于阵列传感器的数量（即L>M）的现象。在这种情况下，现有方法无法成功估算DOA。基于子空间的方法和基于稀疏重构的方法都仅在LM时，一些不确定的DOA 估计方法[11]可以成功估计DOA，但是这些方法都是基于源信号稀疏的假设，这是一个限制性假设。针对这一问题本文考虑了更实用的假设，即源信号是独立的，且一组相干信源可以被看作一个独立信源。

通过分析可知，存在相干信号情况下的阵列数据模型只有阵列流形矩阵A不同，显然，当独立信源数小于阵元数，即N≤M时，可以使用ICA 方法估计A。但是由于相干信号的影响，DOA 的估计变得更加复杂，不能用传统的方法估计DOA。针对这一问题，需要寻求一种新的DOA 估计方法。

观察可以发现，式(2)类似于阵列数据模型。当Lk个相干信号同时入射到拥有M个阵元的阵列时，导向矢量an可以看作阵列输出的单个时间样本。由于阵列输出只有一个时间样本，基于子空间的方法和基于ICA 的方法都不能工作。近年来提出的基于稀疏重构的方法在单时间样本条件下仍能成功地估计到达角。

如果Lk

3 复数FastICA 算法

3.1 除噪复数FastICA 算法

一般复数信号ICA 模型为

ICA 的目标是找到一个解混矩阵W，当式(4)中y的所有分量都是独立的时，W达到其最优值。y为M个信源s的估计值，故WA=I,A=W−1，此时W−1为阵列流行矩矩阵A的估计。

真实的应用中往往存在噪声，含噪声复数信号ICA 模型为

一般情况下，获取的观测数据都具有相关性，需要对数据进行去相关处理。白化是一种常见的去相关方式，主成分分析（PCA,principal component analysis）可用于白化处理，其具有降维功能，可以去除信号噪声的干扰。在观测信号的维度数大于实际独立信源的个数，即M>N时，可以经PCA 降维使观测信号与实际信源数相等。经过白化后观测数据的协方差是单位矩阵，各维数据不相关。故经过PCA 处理可以降低观测数据的维度，且去除各观测数据间的相关性，从而降低在ICA 中问题的复杂度。

使用PCA 对数据白化，首先，对观测数据x各个分量进行中心化，即

使用中心化后的数据计算协方差矩阵，即R=E{xxH}∈CM×M，对R进行奇异值分解R=UΣVH，其中U、V分别为M×M的左、右奇异向量矩阵，且均为酉矩阵；Σ为M×M的奇异值对角矩阵。

ICA 的前提条件是独立信源的准确估计，当非相干信源与相干信源共存时，独立信源数N=NI+K。Bazzi等[12]提出了一种改进的MDL（MMDL,modified minimum description length）算法，在含噪声观测数据、快拍数较少的情况下，对信源估计的准确率优于传统MDL。使用MMDL 方法估计出独立信源数N，取U中前N列的向量组成UN，对数据进行PCA 降维处理，即

经过PCA 降维后，对降维数据进行白化处理，白化后的数据满足。降维后的每一维数据是独立的，为确保每一维数据具有单位方差，需对每一维数据除以标准差来进行缩放，则白化后的数据为

降维白化矩阵可表示为

数据经过白化处理后，在复数FastICA 算法[13]中，代价函数为

其中，G为非线性函数，常见的形式有G1(y)=y2，；y为自变量；b为任取的常数值，通常b=0.1。

Bingham 等[13-14]提出复数FastICA 算法的圆性信号学习规则为

Novey 等[15]提出复数FastICA 算法的非圆信号学习规则为

要进行N个信源的估计，需进行N次迭代计算，最终解混矩阵WH=[w1,w2,…,wN]∈CN×N，在完成所有迭代后，对W进行对称正交变换，去除W间的相关性，避免wn收敛到相同极值，即。

除噪复数FastICA算法具体步骤如算法1所示。

算法1除噪复数FastICA 算法

3.2 去噪复数FastICA 算法

当SNR 较大的情况下，噪声对信号的影响可以忽略，使用上文中的除噪算法，在计算时只选取信号子空间，使用PCA 进行降维直接去除噪声的影响。当信噪比SNR 较小时，需要考虑到噪声对信号的影响。

使用式(6)对观测数据进行中心化，并计算协方差矩阵R=E{xxH}∈CM×M。对协方差矩阵进行奇异值分解，即

奇异值对角矩阵Σ理论上应满足如下关系：。

奇异值前N个对应信号子空间，后M−N个对应噪声方差δ2。但在实际计算时，这些奇异值往往不是相同的值，故无法直接得到噪声方差δ2。

若已知独立信源数N，通过噪声子空间特性，可近似估计噪声方差

则N个信源组成的信号子空间奇异值矩阵为

依据白化原则，对观测数据进行降维和去相关线性变换，因为去噪白化操作后的数据不再满足协方差矩阵为单位矩阵，用“近白化”代替白化，则近白化矩阵为

其中，UN为矩阵U的前N列。

近白化后的观测矩阵为

经过去噪变换后，近白化后的数据协方差矩阵不再是正交矩阵，结合式(14)推导及，可得

Bingham等[13]和Novey等[15]在对式(12)和式(13)的推导中均使用了近似关系E{uuHf(u)}≈E{uuH}E{f(u)}，其中，u为白化处理后的数据E{uuH}=I，故

依据式(20)推导结果，在去噪情况下式(22)应为

根据以上推导，去噪复数FastICA 算法的学习规则为

信源s的估计表示为

去噪复数FastICA算法具体步骤如算法2所示。

算法2去噪复数FastICA 算法

上述2 种去噪算法中，算法1 直接使用PCA去除噪声部分，针对信号强度相当、信噪比高、噪声影响小的情况，利用白化处理把观测信号全部投影到信号子空间，在复数FastICA 中也没有考虑噪声的影响。算法2 适用于存在微弱信号或噪声过强等强干扰的情况，根据噪声高斯、独立性，对噪声进行了估计，在白化处理中，在信号子空间中减去噪声方差，滤除信号子空间中噪声的干扰，并在复数FastICA 算法的迭代中，考虑了噪声的影响对学习规则进行改进，使盲分离算法更适用于存在噪声的情况和微弱信号的提取。

4 基于稀疏重构的DOA 估计

从ICA 的角度来看，接收到的阵列数据是由不同方向接收到的独立信号混合而成的。混合矩阵为A，可以用复数FastICA 估计分离矩阵，则A可以用分离矩阵来估计。

对A中的每一列进行稀疏重构，稀疏重构是由观测向量重构稀疏信号的过程。只要满足信号稀疏或可压缩条件，就可用一个观测矩阵对矩阵进行稀疏表示，稀疏表示后的信号越稀疏则重建后的精度就越高，所包含的信源DOA 信息越全面。对于稀疏表示的信号，通过求解一个优化问题就能近似地重构信号。观测信号中信源具有稀疏性，且每次都对A中一列即单个独立信源进行重构，稀疏重构理论正是运用信号在观测空间上的稀疏特性来实现对入射信号角度估计的。

经过以上2 种算法分析计算，假设已经使用算法 2 估计去噪白化处理后的阵列流行矩阵，其中#表示矩阵的伪逆。虽然ICA不能识别出矩阵A的导向向量顺序，但在本文所提DOA 估计方法中并不会产生影响。假设an是的第n列即第n个信源的导向向量。算法的目标就是通过估计的导向向量an，求得其中非相干信源的到达角度nθ或相干组信号到达角度θ1,θ2,…,θlk。

为将导向向量an作为稀疏表示问题，引入一个包括所有可能入射信号DOA 的过完备字典D，假设Θ=[α1,…,αQ]是所有可能的DOA 组成的抽样网络。矩阵Θ中的列表示潜在入射信号DOA 对应的导向矢量，D=[d(α1),d(α2),…,d(αQ)]，

D是已知的，并且与信号源的实际位置无关。的空间过完备表示为

其中，为小噪声部分，v=[v1,v2,…,vQ]。如果α q=θlk，则v q=μlk；否则vq=0。将DOA 估计问题转化成v的稀疏谱问题，v中包含入射信号的真实DOA 处主峰。可以解决式(27)的逆问题，通过L1 范数的方法将其正则化，使其有利于信号稀疏。最小化v的空间谱，则代价函数为

由于目标函数的数据是复数值，因此，无论是线性规划还是二次规划都不适用于其优化。此时，可以采用二阶锥规划（SOCP，second order cone programming）[16]，把复数凸优化问题转化成SOCP问题，使在内点方法（IPM,interior point method）框架下实现高效的全局收敛算法。

把L1 范数形式转化为SOCP 形式。根据Lobo等[17]对SOCP 的推导与分析。把式(28)转换成SOCP形式，使用辅助向量e和f把目标函数进行线性表示，将非线性部分代入约束条件，则优化目标可转换成如下形式

故SOCP 最终优化问题表示为

式(30)中重要的部分是正则化参数β的选择，其平衡了实际问题中数据的拟合和稀疏性。Malioutov 等[6]所提的L1-SVD 算法中已经给出关于β的选择方法。β参数的选择是稀疏重构的关键，在文献[6]中β的选择与信噪比有关，需要根据噪声的大小分情况多次计算，特别是在信噪比较小的情况下，参数选择过程较复杂。但在本文中，阵列流行矩阵A的估计经过去噪处理，其中噪声已进行滤除，存在较小的噪声。此时，只需保证β尽可能小，且满足，β选择一个较小的固定值进行计算，即可进行信号稀疏。

本节算法具体步骤如算法3 所示。

算法3基于复数FastICA 与稀疏重构的DOA估计

5 性能分析

基于ICA 的算法最多能处理M个独立信源。在单个时间样本的情况下，稀疏重构算法最多能处理个入射源信号[18]，因此本文所提算法最多能够处理个入射信号。实际上，如果源信号的幅值满足随机分布，那么稀疏重构算法能够处理s的入射信号数可以达到M−1个。因此，本文所提算法最多能够处理M(M−1)个入射信号。例如，当阵列传感器数量为3 时，该方法可最多可以处理6 个入射信号，其中3 组相干信号，每组相干信号包含2 个相干信号。

在复数FastICA 算法之前，使用PCA 预处理方式。对于阵列流行矩阵A而言，若信源数与阵元数相等，A需要估计M2个元素，经过白化后A的待估计个数降为M(M−1)/2，且A的估计为，白化矩阵经过PCA降维处理为M×N，故W只需估计N2个元素，即A只需估计N2，大幅减少了FastICA 计算量，提高了算法的运算效率。

6 实验仿真分析

为验证ICA 算法的分离效果，仿真实验使用Amari 指数[19]对除噪FastICA 算法和去噪FastICA算法的分离信号性能进行评估。Amari 指数的定义为

其中，N代表信源数；P=WHCwhiteA=FA，即分离矩阵与混合矩阵A的积。IA∈[0,1]，10logΙA的值代表分离效果，值越小分离的效果越好，若10logIA<−10，则表示分离效果较好。

为验证本文所提2 种去噪算法在不同信噪比下的分离效果，使用非线性复数FastICA 算法[15]进行对比实验。进行仿真实验，经过500 次蒙特卡罗实验，在九阵元均匀圆阵，3 个非相干入射信号[30°,260°,120°]，中心频率为500 MHz，快拍数L=500，阵元半径为0.58 m，对3 个信号进行盲分离，分离效果如图2 所示。

从图2 可以看出，当信噪比为[−20 dB,20 dB]时，3 种FastICA 算法均能达到较好的分离效果。当信噪比为[−20 dB,−8 dB]时，去噪复数FastICA 算法优于除噪复数FastICA 算法，平均要好约−1 dB。当信噪比大于−5 dB 时，2 种去噪算法的性能趋于一致。非线性复数FastICA 算法的性能由于没有去噪降维处理，性能整体弱于2 种去噪算法，平均差约−4 dB。

对算法收敛性进行分析，SNR=−5 dB，经过500 次蒙特卡罗实验，在门限设置为10−12、最大迭代500 次情况下，收敛次数对比如表1 所示。非线性FastICA 算法迭代次数平均为68 次，远高于所提的2 种去噪算法。除噪复数FastICA 算法与去噪复数FastICA 算法均表现出良好的算法性能，在迭代15 次以内就取得良好的收敛性。

表1 3 种FastICA 算法收敛次数对比

通过非线性FastICA 算法与2 种去噪FastICA算法的分离效果和收敛性能分析实验，可以看出去噪确实提高了对盲信号的分离性能，特别是算法的收敛速度得到明显提高，分离效果也得以提高。去噪复数 FastICA 算法整体表现优于除噪复数FastICA 算法。

为验证 DOA 估计算法的有效性，使用L1-SVD[6]、MUSIC[1]算法和极大似然估计（DML,deterministic maximum likelihood）[1]与本文所提算法进行对比，分别进行实测数据实验、相干信号欠定仿真实验、非相干信号超定仿真实验、低信噪比和信源小间隔实验等，对算法的性能进行分析验证。

使用九阵元圆阵接收机，接收机阵元半径为0.58 m，中心频率为600 MHz，快拍数为8 192。采集实际数据，其中包含一个非相干信号103°及一组相干信号[256°,320°]，根据MMDL 判断独立信源数N为2。

表2 给出了真实采集数据的情况下，子空间经典算法、L1-SVD 与本文所提算法的DOA 估计结果。因为有相干信号与非相干信号同时存在，在没有进行解相干的情况，直接判断信源个数，故观测信号中只包含2 个独立信源数，MUSIC 算法和DML 算法都依赖于信源数的准确性，故只能测出2 个信源的DOA。L1-SVD 法和本文所提算法在只提供独立信源数的情况下，仍然可以得到3 个信源的DOA，在本次实验中，L1-SVD算法对其中一个信源320°的结果存在1°偏差，本文所提算法的测向结果较为准确。

表2 实测3 个信源数据在不同算法下DOA 估计结果

3 种算法在包含相干与非相干实测信号空间谱图如图3 所示。从图3 可以看出，3 种算法均在3 个目标信源DOA 处产生谱峰。但是在只有2 个独立信源的前提下，MUSIC 只能得到2 个信源DOA 结果。

表3 是实测6 个信源数据的DOA 结果，使用MMDL 估计信源数N为6。在实际应用中，DML 对多信号的测向效果较差，MUSIC 及L1-SVD 算法存在不同程度的错误，本文所提算法仍保持良好的性能，对于多信号测向表现出较强的分辨力和稳定性。

表3 实测6 个信源数据不同算法DOA 结果

为验证算法对相干信号欠定情况，本文对M(M−1)个相干信源的分辨性能进行了仿真实验，在阵元数为3 的均匀圆阵，信号的中心频率为500 MHz，阵元半径为0.58 m，快拍数L=200，信噪比为10 dB，生成3 组相干信号。通过去噪复数FastICA 算法估计出阵列流行矩阵。对的每一列数据采用L1 范数稀疏重构算法，得到三组相干信号空间谱。欠定情况下基于子空间的MUSIC 等方法失效，因此，实验对比了L1-SVD 与本文所提算法性能，如图4 所示。

L1-SVD 算法只生成了5 个谱峰，在[21°,220°]目标角度正确出现谱峰，由于[82°,92°]和[150°,153°]两两之间角度间隔较小，此算法只在两组信源中间位置产生峰值，在360°处产生一错误峰值，L1-SVD法在本实验中表现较差。本文所提算法在三组相干信源的真实DOA 处均出现峰值，每组相干信号的真实DOA 处分别产生2 个峰值。

下面的仿真实验验证了本文所提算法对最大数量非相干信源分辨性能。在六阵元均匀圆阵，信号的中心频率为500 MHz，阵元半径为0.58 m，快拍数L=200，信噪比为−5 dB，生成6 个非相干信源[30°,60°,120°,150°,180°,210°]。当信源数与阵元数一致时，基于子空间的MUSIC 算法将失效，因此，实验对比了L1-SVD 与本文所提算法的性能，如图5 所示。

从图5 可以看出，2 种算法均可测出6 个信源的DOA，但L1-SVD 算法在最后一个信源，即210°处有少许偏差，本文所提算法对6 个信号均能准确定位。

本文所提算法有一个重要的特性就是可以分辨间隔较小的信源，且对噪声有良好的稳健性。将这2个性能结合在一起进行分析。对于超分辨方法来说，分辨率取决于信噪比，如超分辨率方法MUSIC[1]在高信噪比时显示出良好的分辨率，但一旦噪声变得显著，它们的分辨率开始下降。DML[1]在初始化良好的情况下性能较好，但在低信噪比下存在类似的稳健性问题。由于去噪复数FastICA 和稀疏重构算法的结合使用，本文所提算法可以承受更高水平的噪声。

在九阵元均匀圆阵，3 个入射信号分别为15°、20°、25°，快拍数L=200，信噪比为−10 dB，中心频率为500 MHz，阵元半径为0.58 m 条件下，MUSIC、L1-SVD 和本文所提算法进行低信噪比、小间隔仿真实验，结果如 图6～图8 所示。

从图6 和图7 可以看出，在低信噪比、小间隔情况下，MUSIC 和L1-SVD 算法都只能得到2 个信源对应的角度，无法得到正确DOA。从图8 可以看出，本文所提算法分别得到3 个信源对应的正确DOA，对噪声稳健性强，且可定位小间隔的信源。

7 结束语

本文提出了一种基于去噪复数FastICA 和稀疏重构的DOA 估计方法，与传统DOA 方法相比，该方法可以解决入射信号总数欠定情况下，相干信号的DOA 估计问题。当阵列传感器数量为M时，本文所提算法最多可以处理M(M−1)个入射信号。该算法也同样适用于阵元数与信源数相等，即M=N情况下非相干信源的DOA 估计，以及非相干与相干信号共存的DOA 估计；而且，在低噪声、小信源间隔的情况下仍表现出良好的分辨率。经过理论分析、实际数据实验和仿真实验均可证明该算法的有效性。