CN113300714B - 一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法 - Google Patents

一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法 Download PDF

Info

Publication number
CN113300714B
CN113300714B CN202110463103.1A CN202110463103A CN113300714B CN 113300714 B CN113300714 B CN 113300714B CN 202110463103 A CN202110463103 A CN 202110463103A CN 113300714 B CN113300714 B CN 113300714B
Authority
CN
China
Prior art keywords
matrix
signal
reconstruction
iteration
dimension reduction
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110463103.1A
Other languages
English (en)
Other versions
CN113300714A (zh
Inventor
刘素娟
江书阳
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing University of Technology
Original Assignee
Beijing University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing University of Technology filed Critical Beijing University of Technology
Priority to CN202110463103.1A priority Critical patent/CN113300714B/zh
Publication of CN113300714A publication Critical patent/CN113300714A/zh
Application granted granted Critical
Publication of CN113300714B publication Critical patent/CN113300714B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • HELECTRICITY
    • H03ELECTRONIC CIRCUITRY
    • H03MCODING; DECODING; CODE CONVERSION IN GENERAL
    • H03M7/00Conversion of a code where information is represented by a given sequence or number of digits to a code where the same, similar or subset of information is represented by a different sequence or number of digits
    • H03M7/30Compression; Expansion; Suppression of unnecessary data, e.g. redundancy reduction
    • H03M7/55Compression Theory, e.g. compression of random number, repeated compression
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02DCLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
    • Y02D30/00Reducing energy consumption in communication networks
    • Y02D30/70Reducing energy consumption in communication networks in wireless communication networks

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Complex Calculations (AREA)

Abstract

一种基于压缩感知理论的联合稀疏信号降维重构改进方法属于模拟信息转换领域。所述方法是将原多测量向量转化为低维度的多测量向量,然后再利用梯度追踪算法从低维多测量向量中恢复出稀疏解,实现对信号的降维重构。根据基于压缩感知理论的相关定理和参数设置要求,通过理论推导得知对于稀疏度已知的联合稀疏信号,信号宽度达到一定值(称为临界值L1),就可以恢复出唯一解,因此对高维信号进行重构时,考虑将信号宽度降为L1后再进行重构。重构时使用梯度追踪算法,利用无约束最优化方法中的梯度思想来替代逆矩阵或者广义逆矩阵的计算,不必使用QR分解。本发明降低计算复杂度和提高重构成功率,且信号宽度越大,优势更加明显。

Description

一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构 算法
技术领域
本发明涉及一种基于压缩感知(Compressed Sensing,CS)理论的联合稀疏信号(又称为多测量向量,即Multiple Measurement Vectors,MMV)降维梯度追踪重构算法(Dimension Reduction Gradient Pursuit Reconstruction Algorithm,DRGP),属于模拟信息转换、数字信号处理、图像处理等技术领域。
背景技术
基于压缩感知(CS)理论的模拟信息转换器(AIC)极大地缓解了基于传统采样定理的ADC的压力,使得采样方法不再受Shannon-Nyquist采样定理与ADC输入带宽的限制,压缩感知理论利用满足约束等距条件的测量矩阵将信号从高维空间映射到低维空间,使得信号的采样与压缩同时进行。压缩感知理论在图像处理、生物传感、无线通讯和模式识别等领域有着广泛的应用前景,并且很多应用涉及到多重关联信号的离散采集,因此多测量向量(Multiple Measurement Vectors,MMV)问题被广泛研究,MMV问题实质上是单测量向量(Single Measurement Vector,SMV)问题的推广,也称为联合稀疏恢复(Joint SparseRecovery)问题,主要是通过多个测量从相同的感知矩阵中恢复出未知的每一列稀疏向量。压缩感知主要分为信号的稀疏表示、信号的线性测量和信号的重构恢复算法三大模块。基于压缩感知理论的联合稀疏信号重建过程可以理解为对方程Y=AX的求解过程,因为这是一个欠定方程,有无数多个解,所以需要通过特殊的限定条件来确定唯一的解,这种确定唯一解的过程可以统称为压缩感知信号重建算法。
压缩感知的核心目标是提出一种优化算法从相对较少的线性测量中恢复出稀疏信号,不仅要考虑算法重构精度、重构速度,还要考虑重构算法的计算复杂度、硬件资源耗费等因素。现阶段,常用的CS重构算法有贪婪迭代算法和凸优化算法:贪婪迭代算法主要是将信号与原子字典之间的联系作为测量原子(系数)更加有效的一种方式,基本原则就是通过迭代的方式寻找稀疏向量的支撑集,并且使用受限支撑最小二乘估计来重构信号,贪婪迭代算法主要包括正交匹配追踪(Orthogonal Matching Pursuit,OMP)算法及其变体;凸优化算法通过将非凸问题转化为凸问题求解找到信号的逼近,最常用的方法为基追踪(Basic Pursuit,BP)算法,用l1范数替代l0范数来解决最优化问题,以便使用线性编程方法来执行,稀疏梯度投影(Gradient Projection for sparse Reconstruction,GPSR)算法也是一种比较常见的凸优化算法,GPSR算法使用梯度降的方法求解有界约束最优化问题,算法要求投影在可行域中以确保迭代过程的可行性。贪婪迭代类算法迭代过程简单、快速,应用比较广泛,但是对于信号宽度较大的联合稀疏信号,多测量向量直接导致平方复杂度的计算量,计算复杂度会明显增大,不便于硬件实现。凸优化类算法当中的梯度追踪算法和贪婪迭代类算法相比,利用无约束最优化方法中的梯度思想来替代逆矩阵或者广义逆矩阵的计算,进行硬件实现时不必使用QR分解,降低了计算复杂度,减小了存储空间,但梯度追踪算法不太适用于大规模数据的处理。
对于联合稀疏信号,又出现了新的优化算法—ReMBo(Reduce MMV and Boost),此算法将多列测量向量(MMV)问题简化为单列测量向量(SMV)问题再进行信号重构,虽然降低了计算复杂度提高了计算效率,但重构精度明显降低,重构成功率很低,尤其是对于宽度较大的联合稀疏信号,缺点更加突出。因此如何提高联合稀疏信号重构算法的重构成功率、降低算法的计算复杂度、减小功耗以及降低硬件实现难度是本发明的目的。
发明内容
本发明的目的在于提供一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法(Dimension Reduction Gradient Pursuit Reconstruction Algorithm,DRGP),该重构算法利用降维思想以及梯度思想,对ReMBo算法的主体进行降维改进,局部进行降维后的重构算法改进。直接将MMV问题转化为SMV问题使稀疏信号的宽度由原来的维度变为单维,会丢失信号恢复所需要的关键数据,本发明基于压缩感知理论的相关定理和参数要求,在保证能够求出唯一稀疏解的前提下尽可能地减小信号的宽度,即达到在保证高精度的同时尽可能地降低计算复杂度的效果。但无论将信号的维度降低到什么程度,其计算复杂度和直接转化为SMV问题相比还是比较高,因此从降维后所利用的重构算法入手,继续进行降低计算复杂度的改进,用梯度追踪算法的方向更新来替代贪婪迭代算法中正交投影的计算。本发明提出的算法在进行硬件实现时不必使用QR分解来计算逆矩阵或广义逆矩阵,降低了硬件实现的难度,并且和传统的重构算法相比,信号宽度较大时,本发明提出的算法在降低计算复杂度的前提下,重构成功率还能几乎达到100%。
本发明是采用以下技术方案实现的:
一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法,其思想是,由于利用传统贪婪迭代算法对联合稀疏信号直接进行信号重构的计算复杂度高,以及利用ReMBo(Reduce MMV and Boost)算法将多测量向量(MMV)问题转化为单测量向量(SMV)问题再进行信号重构的重构成功率低,本方法降低原联合稀疏信号的宽度之后,再用梯度追踪(Gradient Pursuit,GP)算法恢复出原信号,达到降低计算复杂度和提高精度的目标。在压缩感知理论的基础上,当感知矩阵A的最大线性无关列σ(A)、稀疏度K、观测矩阵的秩rank(Y)满足:σ(A)≥2K-(rank(Y)-1)时,那么就可以恢复出式(1)的唯一稀疏解,结合利用压缩感知理论时压缩比需要满足的要求,即式(2),可推导出在联合稀疏信号的宽度L、稀疏度K和信号长度N满足:L≥(2K+1)-cKlog(N/K)时,信号能够被成功恢复出来,因此考虑将高维信号降维到临界情况,使宽度L满足L=L1=[((2K+1)-cKlog(N/K))]取整+1,其中c是一个常数,近似处理为1。利用降维后的信号进行重构时不再使用传统的贪婪迭代算法,而是使用梯度追踪算法,利用梯度思想来替代最小二乘解的计算,不必计算逆矩阵或广义逆矩阵,也不必使用QR分解,便于硬件实现。
YM×L=AM×NXN×L (1)
M≥cKlog(N/K) (2)
L1=[((2K+1)-cKlog(N/K))]取整+1 (3)
基于CS理论的联合稀疏信号DRGPAlgorithm具体步骤如下:
步骤一、输入输出数据
1.1输入:Y,A,L,ε,K,MaxIters,c;
1.2输出:
Figure GDA0004272023930000041
Γn,flag;
输入数据中:Y是大小为M×L的观测矩阵;A是大小为M×N的感知矩阵;N是信号的长度,M为观测矩阵Y和感知矩阵A的行数,L为信号的宽度,K是联合稀疏信号的稀疏度,这就需要信号稀疏先验;MaxIter用来控制算法主体的最大迭代次数;ε是残差指标,即残差值需小于等于ε,残差是实际值与估计值之间的差,一般信号重构的残差满足小于e-15级数即可,e是自然对数的底数,其值是2.71828...,c是计算降维后的信号宽度时需要用到的一个常数,近似处理为1。
输出数据中:
Figure GDA0004272023930000042
是重构出的联合稀疏信号;Γn是最终支撑集,所有被选择索引的集合;flag是重构是否成功的标志,flag值是false表示没有重构成功,flag值是true表示重构成功。
步骤二、数据初始化
2.1 X0=0N×L
Figure GDA0004272023930000043
iter=1,flag=false;
稀疏解X0的初始值为N×L的零矩阵;索引集Γ0的初始状态为空集;降维迭代次数iter置为1;重构是否成功的标志flag初始状态为false。
2.2计算降维后的信号宽度L1,L1=[((2K+1)-cKlog(N/K))]取整+1;
将初始参数K、N的值带入式(3),式子(2K+1)-cKlog(N/K)的计算结果取整再加1,计算得出L1的值,作为降维维度。
步骤三、DRGP算法主体:While循环语句
主要用来控制每一次的降维迭代,稀疏度K是变化的量,不同的K对应不同的降维维度,本发明的DRGP算法主体是将稀疏度不同的信号降维到不同的维度,然后进行信号重构,若信号重构失败则重新进行降维再重构,直到信号重构成功或降维次数iter超过最大降维迭代次数MaxIter。DRGP算法主体运行结束有两种情况,一种情况是降维迭代次数iter不超过最大降维迭代次数MaxIter的情况下信号重构成功,另一种情况是,信号重构失败但降维迭代次数iter已经超过了最大降维迭代次数MaxIter。
3.1 While语句执行的条件是:迭代次数iter不超过最大迭代次数MaxIters并且重构是否成功的标志flag的状态是flag=false即还没有重构成功,两者不满足其一即停止While循环。
While语句执行的内容是:
3.2随机生成一个大小为L×L1的矩阵B;
3.3观测矩阵Y与随机矩阵B相乘得到降维后的观测矩阵
Figure GDA0004272023930000051
目的是将观测矩阵Y进行降维处理转化为宽度为L1的观测矩阵/>
Figure GDA0004272023930000052
3.4令残差初始值
Figure GDA0004272023930000053
然后利用梯度追踪更新索引集、残差和信号稀疏解,具体见步骤四、DRGP算法局部。
步骤四、DRGP算法局部:梯度追踪算法即For循环语句,在算法主体While循环之内
主要是利用梯度追踪算法进行信号重构,梯度追踪算法是一种迭代算法,迭代次数即For语句的循环次数取决于信号的稀疏度K,因此对于稀疏度K不同的信号,重构算法的迭代次数也不同,并且每一次迭代都要对索引集Γn、残差
Figure GDA0004272023930000054
以及信号的稀疏解/>
Figure GDA0004272023930000055
进行更新。
4.1 For语句控制梯度追踪迭代次数,n=1:K,循环执行K次;
For语句在While循环语句内,也就是While语句每循环一次,for语句被完整执行一次即循环K次,For语句执行过程中,带有n上标的符号都代表第n次迭代。
For语句执行的内容是:
4.2计算降维前的梯度矩阵,相当于正交匹配追踪OMP算法当中残差与感知矩阵的内积计算:gn=<Y,A>,gn是一个大小为N×L的矩阵,Y的每一列分别与感知矩阵A的每一列进行内积计算,得到L个观测列向量;
4.3计算降维后的梯度矩阵:
Figure GDA0004272023930000061
Figure GDA0004272023930000062
是一个大小为N×L1的矩阵,/>
Figure GDA0004272023930000063
的每一列分别与感知矩阵A的每一列进行内积计算,得到L1个观测列向量;
4.4对降维后的梯度矩阵每一行求二范数,行向量的二范数是每一行各个元素平方之和再开根号:
Figure GDA0004272023930000064
得到一个N×1的列向量In,针对这N个数据进行索引;
4.5索引选择
Figure GDA0004272023930000065
Figure GDA0004272023930000066
是列向量In的第i个元素,先对In所有元素取绝对值,然后取绝对值最大的元素所对应的位置索引记为in,i的大小范围是1--N;;
4.6索引支撑集扩充:Γn=Γn-1∪in,Γn-1是上一次迭代即第n-1次迭代的支撑集,Γn是本次迭代即第n次迭代的新的支撑集,将上一步4.5所得的索引in添加到支撑集当中进行支撑集更新;
使用最速下降梯度追踪方法更新方向:
4.7
Figure GDA0004272023930000067
dn初始状态为N×L的零矩阵,提取降维前的梯度矩阵gn中索引集Γn对应行的数据即/>
Figure GDA0004272023930000068
作为矩阵dn对应行的新数据,dn中索引集之外其它行保持为零,得到本次迭代中即第n次迭代更新稀疏解所需的方向矩阵/>
Figure GDA0004272023930000069
4.8
Figure GDA00042720239300000610
Figure GDA00042720239300000611
初始状态为N×L1的零矩阵,提取降维后的梯度矩阵/>
Figure GDA00042720239300000612
中索引集Γn对应行的数据,即/>
Figure GDA00042720239300000613
作为矩阵/>
Figure GDA00042720239300000614
对应行的新数据,/>
Figure GDA00042720239300000615
中索引集之外其它行保持为零,得到本次迭代中即第n次迭代更新残差所需的方向矩阵/>
Figure GDA00042720239300000616
使用最速下降梯度追踪方法更新步长:
4.9计算中间量
Figure GDA00042720239300000617
Figure GDA00042720239300000618
其中/>
Figure GDA00042720239300000619
为矩阵A保留索引集对应列的数据,其他列为零;
4.10计算步长向量
Figure GDA0004272023930000071
其中/>
Figure GDA0004272023930000072
为残差/>
Figure GDA0004272023930000073
与中间量/>
Figure GDA0004272023930000074
的内积计算,/>
Figure GDA0004272023930000075
的列向量与/>
Figure GDA0004272023930000076
的对应列向量进行内积计算,得到L1个内积值,组成一个1×L1的行向量an,/>
Figure GDA0004272023930000077
为上一次迭代即第n-1次迭代的残差值;
4.11计算步长向量an的模an_=|an|=sqrt(sum(abs(an).^2)),作为新的步长,求向量的模即计算向量二范数,向量的二范数是向量各元素的平方之和再开根号;
更新残差:
4.12
Figure GDA0004272023930000078
新的残差/>
Figure GDA0004272023930000079
与上一次迭代的残差/>
Figure GDA00042720239300000710
的差值为步长|an|与中间量/>
Figure GDA00042720239300000711
的乘积;
更新稀疏解:
4.13
Figure GDA00042720239300000712
每次更新稀疏解,就是在上一次迭代所得的稀疏解的基础上加上一项步长向量an的模an_与降维前的方向矩阵/>
Figure GDA00042720239300000713
的乘积;
4.14结束For语句。
步骤五、重构是否成功的标志flag的状态更新:判断残差是否符合要求,For语句之外,While语句之内
残差小于等于残差指标代表信号重构成功,重构是否成功的标志flag=true;残差大于残差指标代表信号重构失败,重构是否成功的标志flag=false;因此残差是否符合要求决定了重构是否成功的标志flag的状态。
5.1如果残差不符合要求,即
Figure GDA00042720239300000714
ε是残差指标,已知的一个常数;
5.2保留参数flag的状态,flag=false;
5.3降维重构迭代次数即While语句执行次数增1,iter=iter+1;
5.4如果残差符合要求,改变参数flag的状态,flag=true。
步骤六、结束While循环,输出结果
6.1 end while;
6.2
Figure GDA0004272023930000081
取最后一次迭代所求得的稀疏解/>
Figure GDA0004272023930000082
作为最终重构出的稀疏信号
Figure GDA0004272023930000083
6.3输出结果
Figure GDA0004272023930000084
索引集Γn以及重构是否成功的标志flag。
本发明的有益效果在于:本发明所述的方法能成功对压缩信号进行重构并能降低计算复杂度,硬件实现的时候不必使用QR分解。在联合稀疏信号的长度和稀疏度已知的情况下,利用推导公式L1=[((2K+1)-Klog(N/K))]取整+1计算出能够保证较高重构成功率的信号最低宽度,然后将信号转化为宽度为L1的低维联合稀疏信号,完成降维处理之后,使用梯度追踪算法对信号进行重构,替代贪婪算法中稀疏解和残差更新的步骤,采用梯度思想,用方向更新来替代最小二乘解的求解。从重构成功率的角度来分析,将信号降维到能够保证唯一解的确切低维度情况,与直接转化为单列信号相比,重构成功率有了明显的提高;从计算复杂度的角度来分析,采用贪婪迭代算法直接重构信号,导致平方复杂度的计算量,而本发明所采用的方法可以有效降低计算复杂度,不必计算逆矩阵或伪逆矩阵,硬件实现时不必使用QR分解,提高了重构速度,降低了硬件实现难度。联合稀疏信号的宽度越大,本发明在提高重构成功率和降低计算复杂度方面的优势更为显著。
附图说明
图1为压缩采样基本结构;
图2为压缩感知测量过程的多测量向量(MMV)模型;
图3中的(a)为ReMBo算法示意图;
图3中的(b)为本发明DRGP算法示意图;
图4为第一组信号的重构成功率对比图;
图5为第二组信号的重构成功率对比图。
具体实施方式
下面将结合附图与实例详细说明本发明的具体实施方式。
如图1所示为基于压缩感知理论的压缩采样基本结构,针对可稀疏表示的信号,该结构在对信号采样的同时对数据进行压缩,既保留了信号恢复所需要的关键数据,又缓解了采样系统的压力,将数据采集和压缩合二为一,采样得到的数据通过信号重构算法恢复出来。本发明属于基于压缩感知理论的模拟信息转换(AIC)技术领域,其基本原理是利用具有稀疏特性的模拟信号与大于奈奎斯特频率的跳变频率的随机序列相乘,然后经过积分器对调制后的信号进行压缩,即低通滤波器的滤波,再用低精度高速度的ADC进行采样,得到包含原始信号信息的少量测量数据,用于后续的信号恢复。
如图2所示为压缩感知测量过程的多测量向量(MMV)模型,即为本发明的研究对象,压缩感知过程实际上是通过测量矩阵将信号从高维空间映射到低维空间的过程,其数学表达式为YM×L=AM×NXN×L。其中矩阵Y为观测矩阵,即压缩测量所得的数据,大小为M×L,A为压缩感知理论中的感知矩阵,大小为M×N,X为联合稀疏信号,大小为N×L,N为信号的长度,L为信号的宽度,M为压缩后的信号长度。基于压缩感知理论的联合稀疏信号重建过程可以理解为对方程Y=AX的求解过程,因为这是一个欠定方程,有无数多个解,所以需要通过特殊的限定条件来确定唯一的解,即求解
Figure GDA0004272023930000091
这种确定唯一解的过程可以统称为信号重构算法。
如图3(a)所示为ReMBo算法示意图,主要分为两步,第一步将多测量向量转化为单测量向量,第二步利用贪婪类算法对信号进行重构;图3(b)所示为本发明DRGP算法示意图,和ReMBo算法相比,改进在两个方面,将多测量向量转化为固定宽度的低维多测量向量,然后再利用梯度追踪算法对信号进行重构。
本发明的DRGP算法伪代码如表1所示,表2是表1的参数说明。
表1
Figure GDA0004272023930000101
Figure GDA0004272023930000111
表2
Figure GDA0004272023930000112
Figure GDA0004272023930000121
Figure GDA0004272023930000131
以如下两组信号为例进行仿真:N=128,M=32,L=30,K=1~16;N=256,M=64,L=50,K=5~35。利用本发明降维梯度追踪重构算法(DRGP算法)对这两组信号进行重构的具体步骤如下:
第一步,输入数据。
1.1ε=5e-25,c=1;
残差指标ε为5×10-25,常数c为1。
第一组信号:
N=128,M=32,L=30,K=1~16,信号长度为128,信号宽度为30,感知矩阵的行数为32,稀疏度K的变化范围是1~16;
感知矩阵A是一个大小为32×128的正态分布随机矩阵;
信号矩阵X是一个大小为128×30的矩阵,其中随机选取K行,将数据变为正态分布的随机数,其他行的数据均为零;
Y32×30=A32×128*X128×30,观测矩阵Y是感知矩阵A和信号矩阵X相乘所得;
MaxIters=18,设置最大迭代次数为18。
第二组信号:
N=256,M=64,L=50,K=5~35,信号长度为256,信号宽度为50,感知矩阵的行数为64,稀疏度K的变化范围是5~35;
感知矩阵A是一个大小为64×256的正态分布随机矩阵;
信号矩阵X是一个大小为256×50的矩阵,其中随机选取K行,将数据变为正态分布的随机数,其他行的数据均为零;
Y64×50=A64×256*X256×50,观测矩阵Y是感知矩阵A和信号矩阵X相乘所得;
MaxIters=30,设置最大迭代次数为30。
第二步,数据初始化。
2.1 X0=0N×L
Figure GDA0004272023930000141
iter=1,flag=false;
稀疏解X0的初始值为N×L的零矩阵,索引集Γ0的初始状态为空集,降维迭代次数iter置为1,重构是否成功的标志flag初始状态为false。
第一组信号:X0=0128×30,稀疏解X0的初始值为128×30的零矩阵;
第二组信号:X0=0256×50,稀疏解X0的初始值为256×50的零矩阵。
2.2计算降维后的信号宽度L1,L1=[((2K+1)-cKlog(N/K))]取整+1;
将初始参数K、N的值带入式(3),式子(2K+1)-cKlog(N/K)的计算结果取整再加1,计算得出L1的值,作为降维维度。
第一组信号:依次将K=1~16,N=128,c=1带入求得每次迭代的L1的值;
第二组信号:依次将K=5~35,N=256,c=1带入求得每次迭代的L1的值。
第三步,利用While语句控制降维迭代次数。
3.1 while(iter≤MaxIters)and(flag is false)do
While语句执行的条件是:迭代次数iter不超过最大迭代次数MaxIters并且重构是否成功的标志flag的状态是false,两者不满足其一即停止While循环。
第一组信号:while(iter≤18)and(flag is false)do,第一组信号的最大降维迭代次数为18;
第二组信号:while(iter≤30)and(flag is false)do,第二组信号的最大降维迭代次数为30。
While语句执行的内容是:
3.2随机生成一个大小为L×L1的矩阵B;
第一组信号L=30,第二组信号L=50。
3.3观测矩阵Y与随机矩阵B相乘得到降维后的观测矩阵
Figure GDA0004272023930000151
目的是将观测矩阵Y进行降维处理转化为宽度为L1的观测矩阵/>
Figure GDA0004272023930000152
3.4令残差初始值
Figure GDA0004272023930000153
然后利用梯度追踪更新索引集、残差和信号稀疏解。
第四步,利用梯度追踪算法对信号进行重构。
4.1For语句控制梯度追踪迭代次数,n=1:K,循环执行K次;
For语句在While循环语句内,也就是While语句每循环一次,for语句被完整执行一次即循环K次。
For语句执行的内容是:
4.2计算降维前的梯度矩阵:gn=<Y,A>,gn是一个大小为N×L的矩阵,Y的每一列分别与感知矩阵A的每一列进行内积计算,得到L个观测列向量;
第一组信号:gn是一个大小为128×30的矩阵,得到30个观测向量;
第二组信号:gn是一个大小为256×50的矩阵,得到50个观测向量。
4.3计算降维后的梯度矩阵:
Figure GDA0004272023930000154
Figure GDA0004272023930000155
是一个大小为N×L1的矩阵,/>
Figure GDA0004272023930000156
的每一列分别与感知矩阵A的每一列进行内积计算,得到L1个观测列向量;
第一组信号:
Figure GDA0004272023930000157
是一个大小为128×L1的矩阵;
第二组信号:
Figure GDA0004272023930000158
是一个大小为256×L1的矩阵。
4.4对降维后的梯度矩阵每一行求二范数,每一行各个元素平方之和再开根号:
Figure GDA0004272023930000159
得到一个N×1的列向量,针对这N个数据进行索引;
第一组信号:得到一个128×1的向量,针对这128个数据进行索引;
第二组信号:得到一个256×1的向量,针对这256个数据进行索引。
4.5索引选择
Figure GDA00042720239300001510
取In绝对值的最大元素所对应的位置索引;
4.6索引支撑集扩充:Γn=Γn-1∪in,Γn-1是上一次迭代即第n-1次迭代的支撑集,Γn是本次迭代即第n次迭代的新的支撑集,将上一步4.5所得的索引in添加到支撑集当中进行支撑集更新;
使用最速下降梯度追踪方法更新方向:
4.7
Figure GDA0004272023930000161
dn初始状态为N×L的零矩阵,提取降维前的梯度矩阵gn中索引集Γn对应行的数据即/>
Figure GDA0004272023930000162
作为矩阵dn对应行的新数据,dn中索引集之外其它行保持为零,得到本次迭代中即第n次迭代更新稀疏解所需的方向矩阵/>
Figure GDA0004272023930000163
4.8
Figure GDA0004272023930000164
Figure GDA0004272023930000165
初始状态为N×L1的零矩阵,提取降维后的梯度矩阵/>
Figure GDA0004272023930000166
中索引集Γn对应行的数据,即/>
Figure GDA0004272023930000167
作为矩阵/>
Figure GDA0004272023930000168
对应行的新数据,/>
Figure GDA0004272023930000169
中索引集之外其它行保持为零,得到本次迭代中即第n次迭代更新残差所需的方向矩阵/>
Figure GDA00042720239300001610
使用最速下降梯度追踪方法更新步长:
4.9计算中间量
Figure GDA00042720239300001611
Figure GDA00042720239300001612
其中/>
Figure GDA00042720239300001613
为矩阵A保留索引集对应列的数据,其他列为零;
4.10计算步长向量
Figure GDA00042720239300001614
其中/>
Figure GDA00042720239300001615
为残差/>
Figure GDA00042720239300001616
与中间量/>
Figure GDA00042720239300001617
的内积计算,/>
Figure GDA00042720239300001618
的列向量与/>
Figure GDA00042720239300001619
的对应列向量进行内积计算,得到L1个内积值,组成一个1×L1的行向量an,其中/>
Figure GDA00042720239300001620
为上一次迭代即第n-1次迭代的残差值;
4.11计算步长向量an的模an_=|an|=sqrt(sum(abs(an).^2)),作为新的步长,求向量的模即计算向量二范数,向量的二范数是向量各元素的平方之和再开根号;
更新残差:
4.12
Figure GDA00042720239300001621
新的残差/>
Figure GDA00042720239300001622
与上一次迭代的残差/>
Figure GDA00042720239300001623
的差值为步长|an|与中间量/>
Figure GDA00042720239300001624
的乘积;
更新稀疏解:
4.13
Figure GDA00042720239300001625
每次更新稀疏解,就是在上一次迭代所得的稀疏解的基础上加上一项步长向量an的模an_与降维前的方向矩阵/>
Figure GDA0004272023930000171
的乘积;
4.14结束For语句。
第五步,判断残差是否符合要求,For语句之外,While语句之内。
5.1如果残差不符合要求,即
Figure GDA0004272023930000172
其中残差指标ε是已知的一个常数,大小为5e-25;
5.2保留参数flag的状态,flag=false;
5.3降维重构迭代次数即While语句执行次数增1,iter=iter+1;
5.4如果残差符合要求,改变参数flag的状态,flag=true。
第六步,结束While循环,输出结果。
6.1 end while;
6.2
Figure GDA0004272023930000173
取最后一次迭代所求得的稀疏解/>
Figure GDA0004272023930000174
作为最终重构出的稀疏信号
Figure GDA0004272023930000175
/>
6.3输出结果
Figure GDA0004272023930000176
索引集Γn以及重构是否成功的标志flag。
如图4所示为分别用贪婪算法(M-CoSaMP、M-SP、M-OMP、M-POMP)、ReMBo算法(ReMBo-OMP、ReMBo-CoSaMP)以及本发明DRGP算法对第一组信号进行重构的成功率对比图;如图5所示为分别用贪婪算法(M-CoSaMP、M-SP、M-OMP、M-POMP)、ReMBo算法(ReMBo-OMP、ReMBo-CoSaMP)以及本发明DRGP算法对第二组信号进行重构的成功率对比图。从图4、图5可以看出本发明提出的DRGP算法对信号的重构成功率达到理想状态。

Claims (1)

1.一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法,其特征在于,降低原联合稀疏信号的宽度之后,再用梯度追踪算法恢复出原信号;在压缩感知理论的基础上,当感知矩阵A的最大线性无关列σ(A)、稀疏度K、观测矩阵的秩rank(Y)满足:σ(A)≥2K-(rank(Y)-1)时,那么就能够恢复出式(1)的唯一稀疏解,其中YM×L表示矩阵维度为M×L的观测矩阵,AM×N表示矩阵维度为M×N的感知矩阵,XN×L表示矩阵维度为N×L的被测信号;结合利用压缩感知理论时压缩比需要满足的要求,即式(2),推导出在联合稀疏信号的宽度L、稀疏度K和信号长度N满足:L≥(2K+1)-cKlog(N/K)时,信号能够被成功恢复出来,因此考虑将高维信号降维到临界情况,使宽度L满足L=L1=[((2K+1)-cKlog(N/K))]取整+1,其中c是一个常数,近似为1;利用降维后的信号进行重构时使用梯度追踪算法;
YM×L=AM×NXN×L (1)
M≥cKlog(N/K) (2)
L1=[((2K+1)-cKlog(N/K))]取整+1 (3);
具体步骤如下:
步骤一、输入输出数据
1.1输入:Y,A,L,ε,K,MaxIters,c;
1.2输出:
Figure FDA0004272023920000011
Γn,flag;
输入数据中:Y是大小为M×L的观测矩阵;A是大小为M×N的感知矩阵;N是信号的长度,M为观测矩阵Y和感知矩阵A的行数,L为信号的宽度,K是联合稀疏信号的稀疏度;MaxIter用来控制算法主体的最大迭代次数;ε是残差指标,即残差值需小于等于ε,残差是实际值与估计值之间的差,当信号重构的残差满足小于e-15级数时视为成功重构,e是自然对数的底数,其值是2.71828...,c为1;
输出数据中:
Figure FDA0004272023920000012
是重构出的联合稀疏信号;Γn是最终支撑集,所有被选择索引的集合;flag是重构是否成功的标志,flag值是false表示没有重构成功,flag值是true表示重构成功;
步骤二、数据初始化
2.1X0=0N×L
Figure FDA0004272023920000021
iter=1,flag=false;
稀疏解X0的初始值为N×L的零矩阵;索引集Γ0的初始状态为空集;降维迭代次数iter置为1;重构是否成功的标志flag初始状态为false;
2.2计算降维后的信号宽度L1,L1=[((2K+1)-cKlog(N/K))]取整+1;
将初始参数K、N的值带入式(3),式子(2K+1)-cKlog(N/K)的计算结果取整再加1,计算得出L1的值,作为降维维度;
步骤三、DRGP算法主体:While循环语句
DRGP算法主体是将稀疏度不同的信号降维到不同的维度,然后进行信号重构,若信号重构失败则重新进行降维再重构,直到信号重构成功或降维次数iter超过最大降维迭代次数MaxIter;DRGP算法主体运行结束有两种情况,一种情况是降维迭代次数iter不超过最大降维迭代次数MaxIter的情况下信号重构成功,另一种情况是,信号重构失败但降维迭代次数iter已经超过了最大降维迭代次数MaxIter;
3.1While语句执行的条件是:迭代次数iter不超过最大迭代次数MaxIters并且重构是否成功的标志flag的状态是flag=false即还没有重构成功,两者不满足其一即停止While循环;
While语句执行的内容是:
3.2随机生成一个大小为L×L1的矩阵B;
3.3观测矩阵Y与随机矩阵B相乘得到降维后的观测矩阵
Figure FDA0004272023920000022
目的是将观测矩阵Y进行降维处理转化为宽度为L1的观测矩阵/>
Figure FDA0004272023920000023
3.4令残差初始值
Figure FDA0004272023920000024
步骤四、DRGP算法局部:梯度追踪算法,即For循环语句,在算法主体While循环之内;
4.1For语句控制梯度追踪迭代次数,n=1:K,循环执行K次;
For语句在While循环语句内,也就是While语句每循环一次,for语句被完整执行一次即循环K次,For语句执行过程中,带有n上标的符号都代表第n次迭代;
For语句执行的内容是:
4.2计算降维前的梯度矩阵,相当于正交匹配追踪OMP算法当中残差与感知矩阵的内积计算:gn=<Y,A>,gn是一个大小为N×L的矩阵,Y的每一列分别与感知矩阵A的每一列进行内积计算,得到L个观测列向量;
4.3计算降维后的梯度矩阵:
Figure FDA0004272023920000031
Figure FDA0004272023920000032
是一个大小为N×L1的矩阵,/>
Figure FDA0004272023920000033
的每一列分别与感知矩阵A的每一列进行内积计算,得到L1个观测列向量;
4.4对降维后的梯度矩阵每一行求二范数,行向量的二范数是每一行各个元素平方之和再开根号:
Figure FDA0004272023920000034
得到一个N×1的列向量In,针对这N个数据进行索引;
4.5索引选择
Figure FDA0004272023920000035
Figure FDA0004272023920000036
是列向量In的第i个元素,先对In所有元素取绝对值,然后取绝对值最大的元素所对应的位置索引记为in,i的大小范围是1-N;
4.6索引支撑集扩充:Γn=Γn-1∪in,Γn-1是上一次迭代即第n-1次迭代的支撑集,Γn是本次迭代即第n次迭代的新的支撑集,将上一步4.5所得的索引in添加到支撑集当中进行支撑集更新;
使用最速下降梯度追踪方法更新方向:
4.7
Figure FDA0004272023920000037
dn初始状态为N×L的零矩阵,提取降维前的梯度矩阵gn中索引集Γn对应行的数据即/>
Figure FDA00042720239200000312
作为矩阵dn对应行的新数据,dn中索引集之外其它行保持为零,得到本次迭代中即第n次迭代更新稀疏解所需的方向矩阵/>
Figure FDA0004272023920000038
4.8
Figure FDA0004272023920000039
Figure FDA00042720239200000310
初始状态为N×L1的零矩阵,提取降维后的梯度矩阵/>
Figure FDA00042720239200000311
中索引集Γn对应行的数据,即/>
Figure FDA0004272023920000041
作为矩阵/>
Figure FDA0004272023920000042
对应行的新数据,/>
Figure FDA0004272023920000043
中索引集之外其它行保持为零,得到本次迭代中即第n次迭代更新残差所需的方向矩阵/>
Figure FDA0004272023920000044
使用最速下降梯度追踪方法更新步长:
4.9计算中间量
Figure FDA0004272023920000045
Figure FDA0004272023920000046
其中/>
Figure FDA0004272023920000047
为矩阵A保留索引集对应列的数据,其他列为零;
4.10计算步长向量
Figure FDA0004272023920000048
其中/>
Figure FDA0004272023920000049
为残差/>
Figure FDA00042720239200000410
与中间量/>
Figure FDA00042720239200000411
的内积计算,/>
Figure FDA00042720239200000412
的列向量与/>
Figure FDA00042720239200000413
的对应列向量进行内积计算,得到L1个内积值,组成一个1×L1的行向量an,其中/>
Figure FDA00042720239200000414
为上一次迭代即第n-1次迭代的残差值;
4.11计算步长向量an的模an_=|an|=sqrt(sum(abs(an).^2)),作为新的步长,求向量的模即计算向量二范数,向量的二范数是向量各元素的平方之和再开根号;
更新残差:
4.12
Figure FDA00042720239200000415
新的残差/>
Figure FDA00042720239200000416
与上一次迭代的残差/>
Figure FDA00042720239200000417
的差值为步长|an|与中间量/>
Figure FDA00042720239200000418
的乘积;
更新稀疏解:
4.13
Figure FDA00042720239200000419
每次更新稀疏解,就是在上一次迭代所得的稀疏解的基础上加上一项步长向量an的模an_与降维前的方向矩阵/>
Figure FDA00042720239200000420
的乘积;
4.14结束For语句;
步骤五、重构是否成功的标志flag的状态更新:判断残差是否符合要求即是否在For语句之外,While语句之内;
残差小于等于残差指标代表信号重构成功,重构是否成功的标志flag=true;残差大于残差指标代表信号重构失败,重构是否成功的标志flag=false;因此残差是否符合要求决定了重构是否成功的标志flag的状态;
5.1如果残差不符合要求,即
Figure FDA00042720239200000421
ε是残差指标,已知的一个常数;
5.2保留参数flag的状态,flag=false;
5.3降维重构迭代次数即While语句执行次数增1,iter=iter+1;
5.4如果残差符合要求,改变参数flag的状态,flag=true;
步骤六、结束While循环,输出结果
6.1end while;
6.2
Figure FDA0004272023920000051
取最后一次迭代所求得的稀疏解/>
Figure FDA0004272023920000052
作为最终重构出的稀疏信号/>
Figure FDA0004272023920000053
6.3输出结果
Figure FDA0004272023920000054
索引集Γn以及重构是否成功的标志flag。
CN202110463103.1A 2021-04-23 2021-04-23 一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法 Active CN113300714B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110463103.1A CN113300714B (zh) 2021-04-23 2021-04-23 一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110463103.1A CN113300714B (zh) 2021-04-23 2021-04-23 一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法

Publications (2)

Publication Number Publication Date
CN113300714A CN113300714A (zh) 2021-08-24
CN113300714B true CN113300714B (zh) 2023-07-14

Family

ID=77320383

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110463103.1A Active CN113300714B (zh) 2021-04-23 2021-04-23 一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法

Country Status (1)

Country Link
CN (1) CN113300714B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114034755B (zh) * 2021-10-13 2024-01-12 郑州航空工业管理学院 一种基于发动机气路静电信号的异常颗粒物检测方法
CN114050832A (zh) * 2021-11-17 2022-02-15 重庆邮电大学 一种基于两步深度展开策略的稀疏信号重构方法
CN114111997B (zh) * 2021-11-22 2022-08-12 大连理工大学 一种基于叶端定时欠采样信号特性的叶片同步共振频率恢复方法
CN114375004A (zh) * 2021-12-30 2022-04-19 安徽大学 基于组梯度追踪的低复杂度多用户检测方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107330946A (zh) * 2017-06-05 2017-11-07 中国农业大学 一种基于压缩感知的图像处理方法及装置
CN108322409A (zh) * 2018-01-25 2018-07-24 杭州电子科技大学 基于广义正交匹配追踪算法的稀疏ofdm信道估计方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8861588B2 (en) * 2011-04-04 2014-10-14 The United States Of America As Represented By The Secretary Of The Army Apparatus and method for sampling and reconstruction of wide bandwidth signals below Nyquist rate
WO2018027584A1 (zh) * 2016-08-09 2018-02-15 深圳大学 一种目标属性辅助的压缩感知图像恢复方法及其系统

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107330946A (zh) * 2017-06-05 2017-11-07 中国农业大学 一种基于压缩感知的图像处理方法及装置
CN108322409A (zh) * 2018-01-25 2018-07-24 杭州电子科技大学 基于广义正交匹配追踪算法的稀疏ofdm信道估计方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于稀疏度自适应压缩感知的电容层析成像图像重建算法;吴新杰等;《电子与信息学报》;第40卷(第5期);1250-1256 *

Also Published As

Publication number Publication date
CN113300714A (zh) 2021-08-24

Similar Documents

Publication Publication Date Title
CN113300714B (zh) 一种基于压缩感知理论的联合稀疏信号降维梯度追踪重构算法
CN107689795B (zh) 一种基于实时压缩感知的多地区电力控制方法
CN107527371B (zh) 一种在压缩感知中逼近光滑l0范数的图像重建算法的设计构建方法
CN110084862B (zh) 基于多尺度小波变换与深度学习的图像压缩感知算法
CN111274525B (zh) 一种基于多线性增广拉格朗日乘子法的张量数据恢复方法
CN109447921A (zh) 一种基于重构误差的图像测量矩阵优化方法
CN106961104B (zh) 基于数据分析和组合基函数神经网络的风电功率预测方法
Treister et al. A multilevel iterated-shrinkage approach to $ l_ {1} $ penalized least-squares minimization
CN109756740A (zh) 基于最优测量矩阵的半张量图像压缩方法和图像恢复方法
CN113708771A (zh) 一种基于斯特拉森算法的半张量积压缩感知方法
Chen et al. Sparse linear regression with beta process priors
Mathew et al. Automated regularization parameter selection using continuation based proximal method for compressed sensing MRI
CN116859140A (zh) 基于云边协同的非侵入式负荷监测数据在线压缩感知方法
CN106899305B (zh) 一种基于第二代小波的原始信号重构方法
CN111475768B (zh) 一种基于低相干单位范数紧框架的观测矩阵构造方法
CN110830044B (zh) 基于稀疏最小二乘优化的数据压缩方法
CN107231155A (zh) 一种基于改进StOMP的压缩感知重构算法
CN114545066A (zh) 一种非侵入式负荷监测模型聚合方法和系统
CN112163611A (zh) 一种基于特征张量的高维地震数据插值方法
CN116451747A (zh) 一种基于多尺度特征网络的压缩感知ToF信号重建方法
CN104660259A (zh) 一种自适应电荷再分布模数转换器、转换方法和校准方法
CN114693823B (zh) 一种基于空频双域并行重建的磁共振图像重建方法
Liu et al. MMV Batch Look Ahead Orthogonal Matching Pursuit (M-BLAOMP) Algorithm for Joint Sparse Recovery
Seghouane Robust Structured Dictionary Learning For Block Sparse Representations Using α-Divergence
Malek-Mohammadi et al. Successive concave sparsity approximation: Near-oracle performance in a wide range of sparsity levels

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant