CN105099462A - 一种基于压缩感知的信号处理方法 - Google Patents

Abstract

本申请公开了一种基于压缩感知的信号处理方法，包括：采集原始信号；根据原始信号构造测量矩阵；将测量矩阵与原始信号相乘，得到降维的观测信号，并对观测信号进行后续信号处理。其中，在构造测量矩阵时，首先根据原始信号x的长度N，构造N1×N1(N1≥N)维的标准哈达玛(Hadamard)矩阵H_N1，再对哈达玛矩阵进行截取和采样构成亚采样矩阵，并对亚采样矩阵进行基于SVD分解的奇异值修正，从而得到测量矩阵。应用本申请，能够在解除对原始信号长度限制的基础上提高信号的重构精度。

Description

一种基于压缩感知的信号处理方法

技术领域

本申请涉及信号处理技术，特别涉及一种基于压缩感知的信号处理方法。

背景技术

传统Nyquist采样定理提出，要想从采样信号中精确重构出原始信号，采样速率必须大于信号的两倍带宽。但随着信息爆炸时代的到来以及各种大数据的冲击，信号带宽越来越大，对传统采样方法的速度要求也越来越高，这给硬件实现带来了巨大挑战。压缩感知(CompressedSensing)技术的出现解决了上述问题。依照压缩感知理论，只要原始信号本身或者在某个变换域上是稀疏的，就可以利用一个测量矩阵与原始信号相乘，将原始信号从高维投影到低维，从而使得采样速率摆脱了信号带宽的束缚，并仍能准确地重构出原始信号。压缩感知由于突破了信号处理的瓶颈，因而被引入到了各种研究应用领域：如雷达成像、图像处理、核磁共振、信道估计和模式识别。

压缩感知一个显著的特点是摒弃了传统先采样再压缩的方法，将信号的采样、压缩放在同一步进行，从而减轻了计算处理和硬件存储的压力。下面给出利用压缩感知对信号进行处理的大体过程包括：设x为采样得到的原始信号，长度为N，且为K稀疏(即信号只有K个元素非零)，通过压缩感知则可直接得到观测信号y，长度为M(M<N)，然后再对观测信号y进行存储、传输等后续信号处理；当需要使用原始信号时，例如接收端对信号处理时，可以根据观测信号y准确重构出原始信号x。

在上述处理过程中，原始信号与观测信号的关系(即压缩感知过程)可被表述为

y＝Φx

其中Φ称为感知矩阵或者测量矩阵，大小为M×N。

若信号x本身不稀疏，但在某个变换域上是稀疏的，即x＝Ψs，这里Ψ为稀疏基，s为K稀疏信号，则也可将压缩感知过程重新写为：

y＝Φx＝ΦΨs＝Θs

其中Θ＝ΦΨ。在接收端重构原始信号时，可以先利用某种重构算法重构出稀疏信号s，即

$\begin{array}{c}\min {\left|\right|s\left|\right|}_1\\ s.t.y=\mathrm{\&Theta;s}\end{array}---\left(1\right)$

然后再利用稀疏基Ψ的逆矩阵重建原始信号其中||·||₁为向量的l₁范数。

上述过程中，最重要的实现环节就是通过硬件，用一个与稀疏基不相关的M×N(M＜N)测量矩阵Φ对信号x进行线性投影，得到观测信号y。其目的是为了将原始数据均匀地分散在少量的随机向量上，以期将来能够从这些少量的随机测量值中高概率地恢复出原始数据，所以测量矩阵的性能直接决定了重构信号时重构误差的大小。

由于测量矩阵实现的是对原始信号的降维处理，在重构时所求解的优化方程为欠定方程，也就是方程个数远小于未知数的个数，方程无唯一解。所以，要从y精确找到x的唯一解，Θ需要满足有限等距性质(RestrictedIsometryProperty，RIP)，即对于任意K稀疏信号s和常数δ_k∈(0,1)，有：

$1+{\mathrm{\&delta;}}_k\&le;\frac{{\left|\right|\mathrm{\&Theta;s}\left|\right|}_2^2}{{\left|\right|s\left|\right|}_2^2}\&le;1+{\mathrm{\&delta;}}_k---\left(2\right)$

RIP性质保证信号变化后收敛，是方程组(1)有无确定解的充分条件。但是通过公式(2)证明测量矩阵是否满足RIP性质并不容易。针对这个问题，相关文献给出了等价于RIP性质的三个条件，被称为CS1～CS3。其中，CS1要求测量矩阵的列向量要满足一定的线性独立性；CS2要求测量矩阵的列向量体现了某种类似噪声的独立随机性；CS3则要求达到稀疏度的解时满足l₁范数最小的向量。

因此在考虑实际硬件条件时，优质的测量矩阵应该：1)需要满足RIP等价的CS1～CS3条件，且支持程度越高性能越好；2)矩阵元素以+1、0、-1为主，便于构造和硬件存储，减少成本开销；3)具有普适性，可以应用于多场景。

现有的测量矩阵主要包括亚高斯家族随机矩阵、Toeplitz矩阵和哈达玛(Hadamard)矩阵，其中Hadamard矩阵的整体性能最优。但Hadamard对原始信号有限制，要求其长度必须满足2ⁿ。对于长度不满足2ⁿ的原始信号，如果简单地抽取Hadamard矩阵中的若干列，形成部分Hadamard矩阵，将破坏原来矩阵元素的之间的线性无关性(即CS2条件)，使得测量矩阵无法完全符合RIP特性，重构误差增大。

发明内容

本申请提供一种基于压缩感知的信号处理方法，特别是对于长度不满足2n的原始信号，能够改善信号重构性能。

为实现上述目的，本申请采用如下的技术方案：

一种基于压缩感知的信号处理方法，包括：

采集原始信号，根据所述原始信号构造测量矩阵；

将测量矩阵与原始信号相乘，得到观测信号，并对所述观测信号进行信号处理；其特征在于，根据所述原始信号构造测量矩阵的方式包括：

a、根据原始信号x的长度N，构造N1×N1维的标准哈达玛矩阵H_N1，将其作为当前哈达玛矩阵；其中，N1＝min{n1|n1≥Nandn1＝2ⁿ,n∈Z⁺}；

b、提取当前哈达玛矩阵的前N×N部分，得到部分哈达玛矩阵，随机提取所述部分哈达玛矩阵N行中的M行构成M×N维的亚采样矩阵；其中，M为根据传感器硬件确定的正整数；

c、对所述亚采样矩阵进行SVD分解，将分解得到的矩阵 $\left[\begin{array}{cc}\mathrm{\&Sigma;}& 0\\ 0& 0\end{array}\right]$ 进行奇异值平均化，并对矩阵 $U\left[\begin{array}{cc}\mathrm{\&Sigma;}& 0\\ 0& 0\end{array}\right]V^T$ 进行归一化处理，得到所述测量矩阵。

较佳地，在步骤a和步骤b之间，该方法进一步包括：构造对角线元素满足伯努利分布的N1×N1维的对角矩阵R_N1；将所述标准哈达玛矩阵H_N1与所述对角矩阵R_N1相乘得到亚高斯化的哈达玛矩阵H_N1R_N1，将所述亚高斯化的哈达玛矩阵H_N1R_N1作为当前哈达玛矩阵。

较佳地，所述构造N1×N1维的标准哈达玛矩阵包括：选取矩阵 $U=\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]$ 作为种子矩阵，利用哈达玛种子迭代法构造所述标准哈达玛矩阵。

由上述技术方案可见，本申请在进行信号处理时，采集原始信号，根据原始信号构造测量矩阵；将测量矩阵与原始信号相乘，得到观测信号；最后可对观测信号进行进一步的信号处理。其中，在构造测量矩阵时，根据原始信号x的长度N，构造N1×N1维的标准哈达玛矩阵H_N1，再对哈达玛矩阵进行截取和采样构成亚采样矩阵，并对亚采样矩阵进行基于SVD分解的奇异值修正，从而得到测量矩阵。通过本申请的方式，在利用哈达玛矩阵构造测量矩阵时，不需要限定原始信号的长度；并通过基于SVD分解的奇异值修正，降低矩阵的线性相关性，并能够在重构时去除信号在传输过程中遇到的高斯白噪声，从而提高信号重构的性能。

附图说明

图1为本申请中基于压缩感知的信号处理方法流程图；

图2为本申请实例中风速信号的模拟曲线示意图；

图3为本申请的测量矩阵与现有测量矩阵的内部相关性比较示意图；

图4为利用本申请的SVD-SPHR矩阵与现有其他矩阵对原始信号处理后的重构误差对比曲线示意图；

图5为利用本申请的SVD-SPHR矩阵与现有其他矩阵对原始信号处理后的重构信噪比对比曲线示意图。

具体实施方式

为了使本申请的目的、技术手段和优点更加清楚明白，以下结合附图对本申请做进一步详细说明。

本申请中，在Hardamard矩阵的基础上构造测量矩阵，通过预处理和亚采样摆脱信号长度2的幂次方限制，同时利用奇异值分解降低矩阵行列之间的线性相关性，并能作用于重构算法进行间接去噪，从而提高信号重构的性能。进一步优选地，还可以引进伯努利R矩阵对哈达玛矩阵进行亚高斯化，从而提升对RIP等价的CS1和CS2条件支持程度，进一步改善信号重构的性能。

具体地，图1为本申请中基于压缩感知的信号处理方法的基本流程示意图。如图1所示，该方法包括：

步骤101，采集原始信号。

本步骤的处理与现有技术相同，利用各种传感器采集原始信号。对于不同的信号处理领域和系统，可以采集不同类别的信号。例如，对于图像处理，可以采集原始图像信号；对于信道估计，可以采集参考符号的信号等。

步骤102，根据采集的原始信号的长度N确定标准Hadamard矩阵的行数和列数N1。

本申请中，N1＝min{n1|n1≥Nandn1＝2ⁿ,n∈Z⁺}。

步骤103，构造N1×N1维标准Hadamard矩阵H_N1。

标准哈达玛矩阵需要利用一个种子矩阵，采用Hadamard种子迭代法构造完成。所选的“种子”矩阵U必须满足哈达玛正交条件UU^T＝E。本申请中，选 $U=\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right],$ 然后用循环直积的方法构造出N1×N1维标准Hadamard矩阵H_N1，即H₁＝U，i＝0,1,...,N1。H_i的行、列数是H_i-1的2倍。例如： $H_2=H_1\&CircleTimes;U=\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]\&CircleTimes;\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]=\left[\begin{array}{cccc}1& 1& 1& 1\\ 1& -1& 1& -1\\ 1& 1& -1& -1\\ 1& -1& -1& 1\end{array}\right].$ 这里，选取 $U=\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]$ 作为种子矩阵，有利于节省存储空间。在硬件实现时，只需在最初的硬件设备中存下2×2的、且矩阵元素都是+1或者-1的种子矩阵，则每次进行压缩感知处理时，就可以利用该种子矩阵进行循环直积得到N1×N1的标准Hadamard方阵。

步骤104，提取标准Hadamard矩阵的前N×N部分，构成部分Hadamard矩阵；再随机提取该部分Hadamard矩阵的N行中的M行，构造M×N维的亚采样矩阵S_M×N。

本申请中，M的取值由传感器的硬件决定，同时满足与Klog(N/K)在同一数量级的要求，即M＝ο(Klog(N/K))。

在目前的压缩感知技术中，要利用标准Hadamard矩阵，原始信号的维数必须满足2的幂次方要求。但实际信号并非都能满足此要求，于是本申请利用步骤102先由信号长度N计算出大于信号长度且满足2ⁿ的最小数N1，然后在本步骤进行截取，这样就摆脱了2的幂次方限制。具体本步骤中进行截取时，首先将生成的N1×N1维的标准Hadamard矩阵提取前N行N列，形成部分Hadamard矩阵，并经过一个亚采样过程，形成亚采样矩阵，并满足M×N的要求，使得信号经过压缩感知后其长度能够被压缩为M。

步骤105，对亚采样矩阵进行基于SVD分解的奇异值修正，得到测量矩阵。

在上述步骤104的处理中，对Hadamard矩阵的截取和亚采样破坏了矩阵行列之间的不相关性，为降低矩阵的线性相关性，在本步骤中，对亚采样矩阵进行基于SVD分解的奇异值修正。

具体地，将步骤104得到的亚采样矩阵S_M×N进行SVD分解，即 $S_{M\&times;N}=U\left[\begin{array}{cc}\mathrm{\&Sigma;}& 0\\ 0& 0\end{array}\right]V^T,$ 其中U和V都是正交矩阵，Σ矩阵中的对角线元素为σ_i，即Σ_ii＝σ_i，σ_i为矩阵的奇异值。进行奇异值修正时，保持U和V不变，对 $\left[\begin{array}{cc}\mathrm{\&Sigma;}& 0\\ 0& 0\end{array}\right]$ 进行奇异值修正为 $\left[\begin{array}{cc}{\text{\&Sigma;}}^{\&prime;}& 0\\ 0& 0\end{array}\right],$ 即将奇异值平均化。最后归一化得到测量矩阵 ${S^{\&prime;}}_{M\&times;N}=U\left[\begin{array}{cc}{\mathrm{\&Sigma;}}^{\&prime;}& 0\\ 0& 0\end{array}\right]V^T.$

经过上述特征值分解和修正，降低了矩阵行列间的线性相关性，同时，在信号重构时还可以去除信号在传输过程中遇到的高斯白噪声。

步骤106，将步骤105得到的测量矩阵与原始信号相乘，得到观测信号，并可对该观测信号进行进一步的信号处理。

通过测量矩阵与原始信号的相乘处理，实现了对原始信号的降维，从而方便后续的信号处理。具体本步骤的处理与现有压缩感知的信号处理过程相同，例如存储、传输等，这里就不再赘述。

至此，本申请中基于压缩感知的信号处理方法的基本流程结束。通过上述处理，一方面利用步骤102和步骤104的处理使原始信号的长度摆脱了2的幂次方的限制，同时利用步骤105的处理降低了由于Hadamard矩阵的截取和亚采样而引入的线形相关性，并去除了信号在传输过程中遇到的高斯白噪声，从而对于长度不符合2的幂次方的原始信号重构，明显提高了信号的重构性能。

上述图1所示的流程仅为本申请中的基本处理流程。在该流程基础上，优选地，还可以在对标准哈达玛矩阵进行截取和亚采样前，进一步对标准哈达玛矩阵进行亚高斯化处理，以更好地满足CS2条件，进一步提高信号重构的性能。具体地，可以在上述图1所示流程的步骤102和104之间进一步包括如下处理：

步骤102a，根据步骤102中得到的N1构造一个对角矩阵R_N1×N1，R_N1×N1的对角线元素满足伯努利分布。

对角矩阵R_N1×N1的对角线元素满足伯努利分布，即对角线元素取+1和-1的概率均为1/2，即p(R_ii＝±1)＝1/2，非对角线元素为0。本步骤的处理可以在步骤103之前、或与步骤103同时，或在步骤103之后执行。

步骤103a，将步骤103得到的标准哈达玛矩阵H_N1与步骤102a得到的对角矩阵R_N1×N1相乘，得到亚高斯化的哈达玛矩阵H_N1R_N1。

得到对角矩阵R_N1×N1后，将标准Hadamard矩阵与该对角矩阵R_N1×N1相乘，得到亚高斯化的Hadamard矩阵H_N1R_N1。在接下来的步骤104中，对亚高斯化的Hadamard矩阵进行截取和采样，得到亚采样矩阵，为与前述根据标准哈达玛矩阵进行截取和采样得到的亚采样矩阵相区别，这里将亚采样后的矩阵称为SPHR矩阵。后续步骤105中对SPHR矩阵进行基于SVD分解的奇异值修正，得到测量矩阵，称为SVD-SPHR矩阵。

由于亚高斯族的测量矩阵往往有很好的随机特性，比Hadamard矩阵更加满足CS2条件。因此，优选地，在上述处理中，通过步骤102a构建一个满足伯努利分布的对角矩阵R，伯努利矩阵R本身是亚高斯随机矩阵，因此经过R优化的标准Hadamard矩阵，其亚高斯噪声特性可以得到进一步的提升，使本申请能够很好地满足CS2条件，从而进一步提高信号重构性能。上述本申请的信号处理方法，不仅适用于长度不符合2的幂次方的原始信号，对于长度符合2的幂次方的原始信号也同样适用。

下面通过在智能电网中风速信号的处理为例，说明本申请的信号处理方法，其中，以引入亚高斯化的处理为例进行说明。

选择一维组合风速模型作为原始信号，该组合风速信号包含如下所示的渐变风、阵风、随机风、平均风。

(1)平均风：

V＝K，K为常数。

(2)阵风：

$V=\left\{\begin{array}{c}\frac{G_{\max }}{2}[1-\cos 2\mathrm{\&pi;}\left(\frac{t-T_{1g}}{T_g}\right)],T_{1g}<t<T_{1g}+T_g\\ 0,\mathrm{others}\end{array}\right.$

其中，T_g为阵风周期，T_1g为阵风开始时间，t为时间。

(3)阶跃风：

$V=\left\{\begin{array}{cc}0& t<T_{1t}\\ R_{\max }\frac{t-T_{1t}}{T_{2t}-T_{1t}}& T_{1t}<t<T_{2t}\\ R_{\max }& T_{2t}<t<T_{2t}+T_t\\ 0& t>T_{2t}+T_t\end{array}\right.$

(4)随机风：

其中，ω_i＝(i-0.5)Δω， $S_v\left({\mathrm{\&omega;}}_i\right)=\frac{2K_N{F}^2\left|{\mathrm{\&omega;}}_i\right|}{{\mathrm{\&pi;}}^2{[1+{(F{\mathrm{\&omega;}}_i/\mathrm{\&mu;\&pi;})}^2]}^{4/3}}.$

仿真50s内风速信号的情况。如图2所示，0～20s是平均风，K＝7，5～20s出现了阵风，20～40s是阶跃风，最后10s是随机风。每0.05s风速传感器收集一个数据，一共167个数据点，也就是采集长度为167的原始信号。所做仿真实验对本申请信号处理方法中构造的测量矩阵的内部相关性、原始信号的重构误差和重构信噪比三个方面进行验证。

由于原始信号长度为167，不符合2ⁿ，在进行压缩感知时无法利用标准的Hadamard矩阵作为测量矩阵，所以仿真实验比较了本申请所构造的SVD-SPHR矩阵和高斯随机矩阵、部分Toeplitz矩阵的内部相关性。图3所示的实验结果表明，SVD-SPHR矩阵的内部相关性比高斯随机矩阵和部分Toeplitz矩阵的都小。可见，SVD-SPHR矩阵仍然保持了部分Hadamard矩阵接近最优的CS1特性，而它又弥补了Hadamard矩阵对于原始信号限制的缺陷，不再要求2ⁿ的长度，因而普适性更好。

在其余的两个实验中，利用本申请及现有技术的几种测量矩阵对原始信号进行压缩感知处理，并采用离散小波基(DWT)为稀疏基对原始信号进行重构，重构前人工加入一个高斯白噪声，用最经典的OMP(OrthogonalMatchingPursuit)重构出原始信号，比较各种测量矩阵对原始信号进行处理后的重构误差和重构信噪比。由于构造矩阵时随机性较强，为保证公平性和可靠性，所有实验结果取1000次测量的平均值。图4反映出重构误差上SVD-SPHR矩阵比强制2ⁿ的部分Hadamard矩阵、高斯随机矩阵以及部分Toeplitz随机矩阵的都小。而图5证明了SVD-SPHR的重构信噪比比其他矩阵的都大。可见，本申请中所构造的SVD-SPHR矩阵能够起到间接去噪的作用，提高了压缩感知的重构准确率。

应用上述本申请的信号处理方法，能够去除对原始信号长度的限制，提升RIP特性，并能够结合重构算法间接去噪，同时，硬件存储+1/-1容易实现，大大改善了压缩感知处理中的信号重构性能。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明保护的范围之内。

Claims (3)

1.一种基于压缩感知的信号处理方法，包括：

采集原始信号，根据所述原始信号构造测量矩阵；

将测量矩阵与原始信号相乘，得到观测信号，并对所述观测信号进行信号处理；其特征在于，根据所述原始信号构造测量矩阵的方式包括：

a、根据原始信号x的长度N，构造N1×N1维的标准哈达玛矩阵H_N1，将其作为当前哈达玛矩阵；其中，N1＝min{n1|n1≥Nandn1＝2ⁿ,n∈Z⁺}；

b、提取当前哈达玛矩阵的前N×N部分，得到部分哈达玛矩阵，随机提取所述部分哈达玛矩阵N行中的M行构成M×N维的亚采样矩阵；其中，M为根据传感器硬件确定的正整数；

c、对所述亚采样矩阵进行SVD分解，将分解得到的矩阵 $\left[\begin{array}{cc}\mathrm{\&Sigma;}& 0\\ 0& 0\end{array}\right]$ 进行奇异值平均化，并对矩阵 $U\left[\begin{array}{cc}\mathrm{\&Sigma;}& 0\\ 0& 0\end{array}\right]V^T$ 进行归一化处理，得到所述测量矩阵。

2.根据权利要求1所述的方法，其特征在于，在步骤a和步骤b之间，该方法进一步包括：构造对角线元素满足伯努利分布的N1×N1维的对角矩阵R_N1；将所述标准哈达玛矩阵H_N1与所述对角矩阵R_N1相乘得到亚高斯化的哈达玛矩阵H_N1R_N1，将所述亚高斯化的哈达玛矩阵H_N1R_N1作为当前哈达玛矩阵。

3.根据权利要求1或2所述的方法，其特征在于，所述构造N1×N1维的标准哈达玛矩阵包括：选取矩阵 $U=\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]$ 作为种子矩阵，利用哈达玛种子迭代法构造所述标准哈达玛矩阵。