CN103036573B - 基于滤波的压缩感知信号采集方法 - Google Patents

Abstract

一种基于滤波的压缩感知信号采集方法，步骤包括：首先，用感知设备在独立的采样周期内对目标信号x(t)进行采集，并用A/D方式对信号进行数字量化；然后，对量化后的信号x(i)进行降维；最后，对降维后的信号进行重构；其中t为采样时刻，i为量化后的信号排序；对量化后的信号进行降维，具体是对量化后的信号通过有限脉冲响应滤波器的差分方程i＝1,…,M，其中h(0),…,h(L-1)为滤波器系数，设计基于滤波的压缩感知信号采集框架，构造如下托普利兹测量矩阵：则观测i＝1,…,M，其中b₁,…,b_L看作滤波器系数；子矩阵Φ_FT的奇异值是格拉姆矩阵G(Φ_F,T)=Φ′_FTΦ_FT特征值的算术根，验证G(Φ_F,T)的所有特征值λ_i∈(1-δ_K,1+δ_K),i＝1,…,T，则Φ_F满足RIP，并通过求解如下l₁最优化问题来重构原信号：

Description

基于滤波的压缩感知信号采集方法

技术领域

本技术方案属于信号处理技术领域，具体是一种基于滤波的压缩感知信号采集方法。

背景技术

随着数字信号处理的发展，系统获取数据的能力不断提高，需要处理的数据量也不断增大。基于香农采样定理的传统信号处理领域存在以下两个关键难点：(1)对宽带信号而言奈奎斯特采样频率过高，导致采样数据量太大；(2)许多系统采用先采样再压缩的数据获取模式，既浪费传感元，又浪费时间、存储空间和带宽资源。这些在一定程度上限制了传统信号处理方法。近几年，D.Donoho、E.Candès及T.Tao等人提出了一种新兴采样技术－压缩感知(CompressedSensing,CS)理论，在模拟/信息采样模式下以远低于奈奎斯特频率的速度成功实现了对信号的同时采样与压缩^[1-2]。CS不再受香农采样定理的限制，无需复杂的信号编码就可通过重构算法从少量采样信号中精确恢复出原信号，大大降低了信号获取时间、采样速率、存储空间和观测设备数量，它已成为一个新的信息科学研究方向，并在信号检测、传感器网络、图像处理、雷达遥感、频谱检测等领域获得了越来越多的应用^[3-7]。

目前，CS研究热点主要涉及三个问题：信号稀疏表示，测量矩阵选择和重构算法设计，其中决定采样硬件实现的测量矩阵的好坏不仅对信号采样和压缩过程有着重要影响，而且其性质直接影响着信号重构速度和效果，因此测量矩阵需满足的性质和硬件构造方式一直是国内外学者研究的重点和难点。2006年，文献[8]首次提出了有限等距特性(Restricted Isometry Principle,RIP)理论，并成为最广泛的判别工具。根据RIP，一些具有普适性和实用性的随机和确定性测量矩阵相继被提出来。常用的随机测量矩阵包括高斯、贝努利和傅立叶随机矩阵等，它们已被验证满足RIP，但因自由元素过多，使得硬件实现存在较大困难，不利于广泛应用。例如文献[3]指出服从伯努利分布的随机信号发生器仍需在奈奎斯特频率下产生随机数据，才能以较低频率进行A/D转换而实现CS采样，这不可避免地增加了硬件成本。为此，许多应用领域都强调降低测量矩阵的随机性，构造确定性矩阵，以降低硬件成本。文献[9]构造多项式测量矩阵，其大小不任意，限制了压缩率，影响了应用性；文献[10]提出了结构化随机矩阵，但其与高斯随机测量矩阵在重构效果上存在差距；文献[11]通过一行元素构造了循环矩阵，将之应用于稀疏信道估计，但精确重构所需的观测个数并未明显减少；文献[12]利用模拟滤波器进行下采样，设计了伪托普利兹测量矩阵，但与随机矩阵相比，信号重构效果没有明显提高。

压缩感知CS基本理论

不同于传统的均匀采样，压缩感知CS理论的核心是利用特定矩阵把一个高维稀疏信号投影到一个低维空间上，然后利用先验的稀疏条件，通过某种线性或非线性模型重构出原信号。具体地，CS通过测量矩阵Φ∈R^M×N对稀疏信号x∈R^N进行非自适应的线性观测获得观测向量y∈R^M，利用l₀最优化问题来精确重构原信号

其中||x||₀＝K，||·||₀测量信号的稀疏性。因为问题(1)是一个NP难题，D.Donoho等学者提出将之转化为l₁最优化问题

$\begin{array}{ccc}\underset{x}{\min }{\left|\right|x\left|\right|}_1& s.t.& y=\mathrm{\Φx}---\left(1\right)\end{array}$

即通过线性规划方法来重构原信号，并称之为BP算法。信号重构的准确性和稳定性决定于Φ是否满足如下的RIP条件。

定义1(RIP)假设稀疏信号x∈R^N，||x||₀＝K，存在常数δ_K∈(0,1)使得测量矩阵Φ∈R^M×N满足

$\left({1-\mathrm{\δ}}_K\right){\left|\right|z\left|\right|}_2^2\≤{\left|\right|{\mathrm{\Φ}}_{T}z\left|\right|}_2^2\≤(1+{\mathrm{\δ}}_K){\left|\right|z\left|\right|}_2^2,\∀z\∈R^{\left|T\right|}---\left(3\right)$

则称Φ具有K阶RIP条件，其中Φ_T∈R^M×T是由集合T中元素所指向的Φ列向量构成的子矩阵。

因为 ${\left|\right|{\mathrm{\Φ}}_{T}z\left|\right|}_2^2=z^{\′}{\mathrm{\Φ}}_T^{\′}{\mathrm{\Φ}}_{T}z,$ RIP表明当Φ的所有子矩阵Φ_T的奇异值 ${\mathrm{\σ}}_i\∈(\sqrt{{1-\mathrm{\δ}}_K},\sqrt{{1+\mathrm{\δ}}_K}),i=1,\·\·\·,T$ 时,重构信号x^*的非零元素个数不大于K，并可求得问题(2)的最优稀疏解。尽管RIP定性和定量地给出了Φ需满足的特性，但它并未涉及硬件设计所需的要求，而Φ的硬件实现才是将CS推向实用的必备条件。例如莱斯大学研制了单像素相机^[13]，麻省理工学院研制了MRI RF脉冲设备、编码孔径相机^[14]，伊利诺伊州立大学研制了DNA微阵列传感器^[15]，中科院研制了压缩感知滤波器和混沌器^[16]等，这些硬件实现将CS向实用化推进了一大步。

发明内容

压缩感知中常选择随机矩阵作为测量矩阵来进行随机线性投影采样，但过多自由元素使得随机矩阵硬件实现、存储和计算困难，因此设计易于硬件实现的测量矩阵是将压缩感知推向实用化的关键。根据信号通过有限脉冲响应滤波器的差分方程，本发明提出一种新的基于滤波的压缩感知信号采集方法，实现了信号在托普利兹测量矩阵下有用信息的高效获取。仿真实验说明托普利兹测量矩阵比随机矩阵更易实现信号采样和重构，并具有硬件实现简单、存储量小、计算复杂度低的优点。

本发明具体技术方案如下：

一种基于滤波的压缩感知信号采集方法，步骤包括：首先，用感知设备在独立的采样周期内对目标信号x(t)进行采集，并用A/D方式对信号进行数字量化；然后，对量化后的信号x(i)进行降维；最后，对降维后的信号进行重构；其中t为采样时刻，i为量化后的信号排序；

对量化后的信号进行降维，具体是对量化后的信号通过有限脉冲响应滤波器的差分方程其中h(0),…,h(L-1)为滤波器系数(L是信号长度，i表示第i次采样)，设计基于滤波的压缩感知信号采集框架，构造如下托普利兹测量矩阵：

则观测其中b₁,…,b_L看作滤波器系数；

子矩阵Φ_FT的奇异值是格拉姆矩阵特征值的算术根，验证G(Φ_F,T)的所有特征值λ_i∈(1-δ_K,1+δ_K),i＝1,…,T，则Φ_F满足RIP，并通过求解如下l₁最优化问题来重构原信号：

$\begin{array}{ccc}\underset{x}{\min }{\left|\right|x\left|\right|}_1& s.t.& y=\mathrm{\Φx},\end{array}$

y为观测向量，x为稀疏信号，

即通过线性规划方法来重构原信号，亦即BP算法。

针对实际可压缩信号，如语音或图像信号的采集，则修改Φ_F为如下形式：

${\mathrm{\Φ}}_S=\left[\begin{array}{cccccccccc}{b}_1& ...& b_L& 0& ...& ...& ...& ...& ...& 0\\ 0& ...& 0& b_1& ...& b_L& 0& ...& ...& 0\\ ...& & ...& & ...& & ...& & ...& \\ 0& ...& 0& ...& 0& ...& 0& b_1& ...& b_L\end{array}\right];$

如果信号在变换基矩阵Ψ上具有稀疏性，则通过求解如下l₁最优化问题，精确重构出原信号：

$\begin{array}{ccc}\underset{\mathrm{\α}}{\min }{\left|\right|\mathrm{\α}\left|\right|}_1& s.t.& y=\mathrm{\Φx}=\mathrm{\Φ\Ξ\α}=\mathrm{\Ξ\α},\end{array}$ 其中Φ与Ψ不相关，Ξ称为CS矩阵。

Φ为测量矩阵，Ψ为基矩阵，y为观测向量，α为稀疏表示系数，x为信号，x＝Ψα表示信号x在Ψ上进行稀疏表示。

本发明基于现有技术，结合信号通过有限脉冲响应滤波器的差分方程，构造了托普利兹测量矩阵，此矩阵不仅以高概率满足RIP，保证信号重构的精确度，而且存储量小、计算复杂度低、易于硬件实现，克服了随机矩阵存在的问题，从而便于实现基于滤波的压缩感知信号采集

附图说明

图1：采用本基于滤波的压缩感知信号采集方法的信号采集原理框架图；

图2(a)～图2(c)为稀疏信号及其重构信号时域波形图，其中，

图2(a)为原始信号时域波形示意图；

图2(b)为基于Φ_F的信号重构波形示意图；

图2(c)为基于Φ_R的信号重构波形示意图；

图3(a)～图3(d)为不同K和R下托普利兹和随机测量矩阵的平均采样时间和平均相对重构误差图，其中，

图3(a)为K下托普利兹和随机测量矩阵的平均采样时间示意图；

图3(b)为K下托普利兹和随机测量矩阵的平均相对重构误差示意图；

图3(c)为R下托普利兹和随机测量矩阵的平均采样时间示意图；

图3(d)为R下托普利兹和随机测量矩阵的平均相对重构误差示意图；

图4(a)～图4(f)为女声和男声语音信号及其重构信号时域波形图，其中，

图4(a)为女声语音信号时域波形示意图；

图4(b)为男声语音信号时域波形示意图；

图4(c)为基于Φ_S的女声重构信号时域波形示意图；

图4(d)为基于Φ_S的男声重构信号时域波形示意图；

图4(e)为基于Φ_R的女声重构信号时域波形示意图；

图4(f)为基于Φ_R的男声重构信号时域波形示意图；

图5(a)～图5(d)为不同R和N下托普利兹和随机测量矩阵的平均采样时间和重构语音信号平均分段信噪比SegSNR图；

图5(a)为R下托普利兹和随机测量矩阵的平均采样时间示意图；

图5(b)为R下托普利兹和随机测量矩阵的重构语音信号平均分段信噪比SegSNR示意图；

图5(c)为N下托普利兹和随机测量矩阵的平均采样时间示意图；

图5(d)为N下托普利兹和随机测量矩阵的重构语音信号平均分段信噪比SegSNR示意图。

具体实施方式

下面结合附图与具体实施方式对本技术方案进一步说明如下：

基于滤波的压缩感知信号采集方法

CS为了获取“边采样边压缩”的观测，典型的物理实现方法有随机下采样^[17]、模拟信息转换器采样^[18]和随机滤波器采样^[19]等。文献[3]基于经典CS实现原理设计了双路A/D随机协同采样，但随机数寄存器存储量大及降维随机投影计算量大，影响了采样效率。

本发明考虑信号x∈R^N通过有限脉冲响应滤波器的差分方程

(其中h(0),…,h(L-1)为滤波器系数^[20])，

并根据文献[3]设计基于滤波的压缩感知信号采集框架，实现多路A/D协同信号采样，如图1。考虑到实际应用中硬件实现的复杂性，可以优先选取双路A/D协同采样。

由图1的信号采样过程构造如下托普利兹测量矩阵

则观测其中b₁,…,b_L可看作滤波器系数。由定义1，子矩阵Φ_FT的奇异值是格拉姆矩阵(Grammian Matrix)G(Φ_F,T)＝Φ′_FTΦ_FT特征值的算术根。如果验证G(Φ_F,T)的所有特征值λ_i∈(1-δ_K,1+δ_K),i＝1,…,T，则Φ_F满足RIP，并可通过l₁最优化以高概率重构出原信号。为此，由文献[11][21]的结论获得本发明中Φ_F满足RIP的定理。

引理1^[21](盖尔斯果林圆盘定理)设矩阵H∈R^N×N，则其所有特征值都落在平面的N个圆盘 $D_i\left(H\right)=\left\{z\right||z-h_{\mathrm{ii}}|\≤\underset{j=1,j\≠i}{\overset{N}{\mathrm{\Σ}}}\left|h_{\mathrm{ij}}\right|\},i=1,...,N$ 的并集中，其中h_ii为圆心， $r_i=\underset{j=1,j\≠i}{\overset{N}{\mathrm{\Σ}}}\left|h_{\mathrm{ij}}\right|$ 为半径。

假设存在整数K≥1以及正数δ_d、δ_o使得δ_d+δ_o＝δ_K∈(0,1)，并且G(Φ_F,T)的对角元素(diagonalelement)和非对角元素(off-diagonal element)分别满足|G_ii(Φ_F,T)-1|<δ_d和|G_ij(Φ_F,T)|<δ_o/K，则由引理1可知：当半径时，G(Φ_F,T)的特征值λ_i∈(1-δ_d-δ_o,1+δ_d+δ_o)＝(1-δ_K,1+δ_K),i＝1,…,T。

引理2^[11]设{u_i}是i.i.d的随机变量序列，并满足|u_i|≤a,E(u_i)＝0,则有

$P\left(\right|\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{u}_i^2-M{\mathrm{\&sigma;}}^2|\&GreaterEqual;t)\&le;2\exp (-\frac{{2t}^2}{M{a}^4})---\left(6\right)$

引理3^[11]设{u_i}和{v_i}是i.i.d的随机变量序列，并满足|u_i|≤a,|v_i|≤a,|u_iv_i|≤a²,E(u_i)＝E(v_i)＝0，则有

$P\left(\right|\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{u}_i{v}_i|\&GreaterEqual;t)\&le;2\exp (-\frac{t^2}{2{\mathrm{Ma}}^4})---\left(7\right)$

定理1设{b_l}是i.i.d的随机变量序列，并满足则当N≥3时使得对 $\&ForAll;M\&GreaterEqual;\frac{{32K}^2{c}^2}{{\mathrm{\&delta;}}_K^2-32c^2{c}_1}\log 3(N^2-N)$ 有

P(Φ_F满足 $\mathrm{RIP}\left({K,\mathrm{\&delta;}}_K\right))\&GreaterEqual;1-\exp (-\frac{c_{1}M}{K^2})---\left(8\right)$

证明由式(5)可知Gram矩阵G＝Φ′_FΦ_F的对角元素和非对角元素分别为

$G_{\mathrm{ij}}=\left\{\begin{array}{cc}{b}_i{b}_j+b_{i-1}{b}_{j-1}+...+b_{{i-t}_j+1}{b}_{{j-t}_j+1}& 1\&le;i\&le;L,1\&le;j\&le;N\\ b_L{b}_{j+L-i}+b_{L-1}{b}_{j+L-i-1}+...+b_{{L-t}_j+1}{b}_{j+{L-i-t}_j+1}& L<i\&le;N,i-L+1\&le;j\&le;i\\ 0& L<i\&le;N,1\&le;j<i-L+1\end{array}\right.$

$G_{\mathrm{ii}}=\left\{\begin{array}{cc}{b}_i^2+b_{i-1}^2+...+b_{{i-s}_i+1}^2& 1\&le;i\&le;L\\ b_L^2+b_{L-1}^2+...+b_{L-s_i+1}^2& L<i\&le;N\end{array}\right.---\left(9\right)$

其中G_ij＝G_ji，整数1≤s_i≤M,1≤t_j≤M,i,j＝1,…,N。

根据引理2和式(9)，G_ii(1≤i≤L)满足

$\begin{array}{c}P\left(\right|G_{\mathrm{ii}}-E\left(G_{\mathrm{ii}}\right)|\&GreaterEqual;{\mathrm{\&delta;}}_d)=P\left(\right|\underset{l=i}{\overset{{i-s}_i+1}{\mathrm{\&Sigma;}}}{b}_l^2-E\left(\underset{l=i}{\overset{{i-s}_i+1}{\mathrm{\&Sigma;}}}{b}_l^2\right)|\&GreaterEqual;{\mathrm{\&delta;}}_d)\\ \&le;P\left(\right|\underset{l=i}{\overset{{i-s}_l+1}{\mathrm{\&Sigma;}}}{b}_l^2-1|\&GreaterEqual;{\mathrm{\&delta;}}_d)\&le;2\exp (-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{s_i{c}^2})\end{array}---\left(10\right)$

$\begin{array}{cc}P\left(\right|G_{\mathrm{ii}}-E\left(G_{\mathrm{ii}}\right)|\&GreaterEqual;{\mathrm{\&delta;}}_d)=P\left(\right|\underset{l=L}{\overset{{L-s}_i+1}{\mathrm{\&Sigma;}}}{b}_l^2-E\left(\underset{l=L}{\overset{{L-s}_i+1}{\mathrm{\&Sigma;}}}{b}_l^2\right)|\&GreaterEqual;{\mathrm{\&delta;}}_d)& \&le;P\left(\right|\underset{l=L}{\overset{{L-s}_i+1}{\mathrm{\&Sigma;}}}{b}_l^2-1|\&GreaterEqual;{\mathrm{\&delta;}}_d)\&le;2\exp (-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{s_i{c}^2})\end{array}---\left(11\right)$

由式(10)和(11)可得

$\begin{array}{c}P\left(\underset{i=1}{\overset{N}{\&cup;}}\right|G_{\mathrm{ii}}-1|\&GreaterEqual;{\mathrm{\&delta;}}_d)\&le;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}2\exp (-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{s_i{c}^2})\\ \&le;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}2\exp (-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{{\mathrm{Mc}}^2})=2\mathrm{Nexp}(-\frac{{2\mathrm{M\&delta;}}_d^2}{c^2})\end{array}---\left(12\right)$

另一方面，由式(9)可见当1≤i≤L,i-j＝1时，中相邻两项(如b_ib_i-1和b_i-1b_i-2)不独立，无法直接利用引理3，于是考虑t_j分别为奇数和偶数的两种情况。当t_j为奇数时，将G_ij拆为

$G_{\mathrm{ij}}=G_{\mathrm{ij}}^1+G_{\mathrm{ij}}^2=(b_i{b}_{i-1}+b_{i-2}{b}_{i-3}+...+b_{{i-t}_j+1}{b}_{{i-t}_j})+(b_{i-1}{b}_{i-2}+b_{i-3}{b}_{i-4}+...+b_{i-t_j}{b}_{{i-t}_j-1}),$ 其中中相邻两项(如b_ib_i-1和b_i-2b_i-3或者b_i-1b_i-2和b_i-3b_i-4)都互相独立，和分别有和项，[·]为取整函数；当t_j为偶数时，将G_ij拆为

$G_{\mathrm{ij}}=G_{\mathrm{ij}}^1+G_{\mathrm{ij}}^2=(b_i{b}_{i-1}+b_{i-2}{b}_{i-3}+...+b_{{i-t}_j}{b}_{{i-t}_j-1})+(b_{i-1}{b}_{i-2}+b_{i-3}{b}_{i-4}+...+b_{i-t_j+1}{b}_{{i-t}_j}),$ 其中和都有 $q_1=q_2=\frac{t_j}{2}$ 项。利用引理3有

$\begin{array}{c}P\left(\right|G_{\mathrm{ij}}|\&GreaterEqual;\frac{{\mathrm{\&delta;}}_o}{K})\&le;P\left(\right\{\left|G_{\mathrm{ij}}^1\right|>\frac{{\mathrm{\&delta;}}_o}{2K}\left\}\mathrm{or}\right\{\left|G_{\mathrm{ij}}^2\right|>\frac{{\mathrm{\&delta;}}_o}{2K}\left\}\right)\\ \&le;2\max \left\{P\left(\right|G_{\mathrm{ij}}^1|>\frac{{\mathrm{\&delta;}}_o}{2K})\mathrm{orP}\left(\right|G_{\mathrm{ij}}^2|>\frac{{\mathrm{\&delta;}}_o}{2K})\right\}\end{array}$

$\begin{array}{c}\&le;2\max \{2\exp (-\frac{{({\mathrm{\&delta;}}_o/2K)}^2}{{2q}_1(c^2/M^2)}),2\exp (-\frac{{({\mathrm{\&delta;}}_o/2K)}^2}{{2q}_2(c^2/M^2)})\}\\ \&le;4\exp (-\frac{{({\mathrm{\&delta;}}_o/2K)}^2}{2M(c^2/M^2)})=4\exp (-\frac{{\mathrm{M\&delta;}}_o^2}{{8K}^2{c}^2})\end{array}---\left(13\right)$

当1≤i≤L,i-j>1时，G_ij中相邻两项都独立，可直接利用引理3，但为了便于计算，也将G_ij按上述方法进行拆分。同理，当L<i≤N时可得如式(13)的结果。因为G_ij＝G_ji，由式(13)可得

$P\left(\underset{i=1}{\overset{N}{\&cup;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&cup;}}\right\{\left|G_{\mathrm{ij}}\right|\&GreaterEqual;\frac{{\mathrm{\&delta;}}_o}{K}\left\}\right)\&le;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}\left(\underset{j=1,j\&NotEqual;i}{\overset{N}{\mathrm{\&Sigma;}}}4\exp (-\frac{{\mathrm{M\&delta;}}_o^2}{{8K}^2{c}^2})\right)\&le;2(N^2-N-M^2+M)\exp (-\frac{{\mathrm{M\&delta;}}_o^2}{{8K}^2{c}^2})---\left(14\right)$

令δ_d＝δ_o＝δ_K/2,N≥3，则有

P(Φ_F不满足RIP(K,δ_K))

$\begin{array}{c}=P\left(\underset{i=1}{\overset{N}{\&cup;}}\right|G_{\mathrm{ii}}-1|\&GreaterEqual;{\mathrm{\&delta;}}_d)+P\left(\underset{i=1}{\overset{N}{\&cup;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&cup;}}\right\{\left|G_{\mathrm{ij}}\right|\&GreaterEqual;\frac{{\mathrm{\&delta;}}_o}{K}\left\}\right)\\ \&le;2\mathrm{Nexp}(-\frac{{2\mathrm{M\&delta;}}_d^2}{c^2})+2(N^2-N-M^2+M)\exp (-\frac{{\mathrm{M\&delta;}}_o^2}{{8K}^2{c}^2})\\ \&le;3(N^2-N)\exp (-\frac{{\mathrm{M\&delta;}}_K^2}{{32K}^2{c}^2})\end{array}---\left(15\right)$

综上，对 $\&ForAll;M\&GreaterEqual;\frac{{32K}^2{c}^2}{{\mathrm{\&delta;}}_K^2-{32c}^2{c}_1}\log 3(N^2-N)$ 有

P(Φ_F满足RIP(K,δ_K))

$=P\left(\underset{i=1}{\overset{N}{\&cup;}}\right|G_{\mathrm{ii}}-1|<{\mathrm{\&delta;}}_d)+P\left(\underset{i=1}{\overset{N}{\&cup;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&cup;}}\right\{\left|G_{\mathrm{ij}}\right|<\frac{{\mathrm{\&delta;}}_o}{K}\left\}\right)\&GreaterEqual;1-\exp (-\frac{c_{1}M}{K^2})---\left(16\right)$

其中 $c_1<{\mathrm{\&delta;}}_K^2/{32c}^2.$

定理1说明Φ_F以高概率满足RIP，因此由y＝Φ_Fx可通过问题(2)以高概率重构出原信号，其中Φ_F只需要生成和存储L个随机数，Φ_Fx需要M×L个乘法操作。然而，被广泛应用的随机矩阵Φ_R需要生成和存储M×N个随机数，Φ_Rx需要M×N个乘法操作。可见，基于滤波的CS采样过程不仅易于硬件实现而且可大大减少系统计算开销，从而利于CS的实际应用。

优选的，CS针对实际语音、图像等可压缩信号提出：如果信号x在变换基矩阵Ψ上具有稀疏性，则可通过求解如下l₁最优化问题

$\begin{array}{ccc}\underset{\mathrm{\&alpha;}}{\min }{\left|\right|\mathrm{\&alpha;}\left|\right|}_1& s.t.& y=\mathrm{\&Phi;x}=\mathrm{\&Phi;\&Psi;\&alpha;}=\mathrm{\&Xi;\&alpha;}\end{array}---\left(17\right)$

精确重构出原信号，其中Φ与Ψ不相关，Ξ称为CS矩阵。文献[16]定义Φ与Ψ的相干度(coherence)为并指出当Φ和Ψ不相干时，即μ很小时，Ξ以高概率满足RIP且保证问题(17)存在稀疏解。为了使Φ_F与正交矩阵的相干度尽可能小，本发明修改Φ_F为如下形式：

${\mathrm{\&Phi;}}_S=\left[\begin{array}{cccccccccc}{b}_1& ...& b_L& 0& ...& ...& ...& ...& ...& 0\\ 0& ...& 0& b_1& ...& b_L& 0& ...& ...& 0\\ ...& & ...& & ...& & ...& & ...& \\ 0& ...& 0& ...& 0& ...& 0& b_1& ...& b_L\end{array}\right]---\left(18\right)$

稀疏矩阵Φ_S通过采集原信号的局部信息，就可精确重构出原信号。这使得CS可应用于一些实际场景，例如能量、资源受限的感知网络。每个感知设备在独立的采样周期内对目标进行信号采集，则所有采样信号为y＝Φ_Sx，汇聚中心通过问题(17)重构出原信号，从而可减少每个感知设备的工作时间和能耗。

下面仿照定理1验证Φ_S也以高概率满足RIP。

定理2设{b_l}是i.i.d的随机变量序列，并满足 $\left|b_l\right|\&le;\sqrt{c/M}(c\&GreaterEqual;1),E\left(b_l\right)=0,E\left(b_l^2\right)=1,$ 则 $\&Exists;c_2<{\mathrm{\&delta;}}_K^2/{8c}^2,$ 当N≥2K时使得对 $\&ForAll;M\&GreaterEqual;\frac{{8K}^2{c}^2}{{\mathrm{\&delta;}}_K^2-{8c}^2{c}_2}\log \frac{{2N}^2}{K}$ 有

P(Φ_S满足 $\mathrm{RIP}\left({K,\mathrm{\&delta;}}_K\right))\&GreaterEqual;1-\exp (-\frac{c_{2}M}{K^2})---\left(19\right)$

证明考虑Φ_S的Gram矩阵

其中 $B=[b_1,...,b_L],B^{\&prime;}B=\left[\begin{array}{cccc}{b}_1^2& b_1{b}_2& ...& b_1{b}_L\\ b_2{b}_1& b_2^2& ...& b_2{b}_L\\ ...& ...& ...& ...\\ b_L{b}_1& b_L{b}_2& ...& b_L^2\end{array}\right].$ 根据定理1可得

$P\left(\underset{i=1}{\overset{N}{\&cup;}}\right|G_{\mathrm{ii}}-1|\&GreaterEqual;{\mathrm{\&delta;}}_d)\&le;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}2\exp (-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{c^2})=2\mathrm{Nexp}(-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{c^2})---\left(21\right)$

$P\left(\underset{i=1}{\overset{N}{\&cup;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&cup;}}\right\{\left|G_{\mathrm{ij}}\right|\&GreaterEqual;\frac{{\mathrm{\&delta;}}_o}{K}\left\}\right)\&le;\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}\left(\underset{j=1,j\&NotEqual;i}{\overset{N}{\mathrm{\&Sigma;}}}2\exp (-\frac{M^2{\mathrm{\&delta;}}_o^2}{{2K}^2{c}^2})\right)\&le;(\frac{N^2}{M}-N)\exp (-\frac{M^2{\mathrm{\&delta;}}_o^2}{{2K}^2{c}^2})---\left(22\right)$

不妨令δ_d＝δ_o＝δ_K/2,N≥2K，则有

P(Φ_S不满足 $\mathrm{RIP}\left(K,{\mathrm{\&delta;}}_K)\right)=P\left(\underset{i=1}{\overset{N}{\&cup;}}\right|G_{\mathrm{ii}}-1|\&GreaterEqual;{\mathrm{\&delta;}}_d)+P\left(\underset{i=1}{\overset{N}{\&cup;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&cup;}}\right\{\left|G_{\mathrm{ij}}\right|\&GreaterEqual;\frac{{\mathrm{\&delta;}}_o}{K}\left\}\right)$

$\&le;2\mathrm{Nexp}(-\frac{{2M}^2{\mathrm{\&delta;}}_d^2}{c^2})+(\frac{N^2}{M}-N)\exp (-\frac{M^2{\mathrm{\&delta;}}_o^2}{{2K}^2{c}^2})\&le;\frac{{2N}^2}{K}\exp (-\frac{{\mathrm{M\&delta;}}_K^2}{{8K}^2{c}^2})---\left(23\right)$

综上，对 $\&ForAll;M\&GreaterEqual;\frac{{8K}^2{c}^2}{{\mathrm{\&delta;}}_K^2-{8c}^2{c}_2}\log \frac{{2N}^2}{K}$ 有

P(Φ_S满足 $\mathrm{RIP}\left({K,\mathrm{\&delta;}}_K\right))=P\left(\underset{i=1}{\overset{N}{\&cup;}}\right|G_{\mathrm{ii}}-1|<{\mathrm{\&delta;}}_d)+P\left(\underset{i=1}{\overset{N}{\&cup;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&cup;}}\right\{\left|G_{\mathrm{ij}}\right|<\frac{{\mathrm{\&delta;}}_o}{K}\left\}\right)\&GreaterEqual;1-\exp (-\frac{c_{2}M}{K^2})---\left(24\right)$

其中 $c_2<{\mathrm{\&delta;}}_K^2/{8c}^2.$

如果Ψ为单位矩阵，比较定理1和定理2易见：Φ_S比Φ_F更稀疏，所以Φ_S下所需观测个数的下界小于Φ_F，从而通过Φ_S进行CS采样可提高信号采集效率。

仿真实验结果及分析

下面在Intel Core 2处理器上进行仿真实验，对比托普利兹和随机测量矩阵下CS采样的计算开销和重构效果，以验证托普利兹测量矩阵的可行性和有效性。

1、稀疏信号的重构

假设某稀疏信号x₀∈R^N(N＝160)有K＝9个非零元素，随机产生非零元素的数值和位置。CS中观测个数M与非零元素个数K有密切关系，通常取M＝cK(c＝3～4)，并定义压缩比为R＝M/N。据此产生观测y＝Φ_Fx₀，其中M＝4K，R＝1/4，Φ_F中随机元素b_l～N(0,1/M)。根据式(4)，图2中基于Φ_F和Φ_R的降维滤波采样和降维随机采样的计算时间分别为2.4013e-004秒和3.2328e-004秒，相对重构误差ε＝||x₀-x^*||₂/||x₀||₂分别为1.6005e-011和6.7728e-014，均具有很好的重构结果。按照上述方法，对某个固定K随机产生50个稀疏信号进行实验。针对不同K和R，图3给出了基于Φ_F和Φ_R的平均采样时间和平均相对重构误差由图3可见随着K和R的增大，所需观测个数M也相应增加，因此利用Φ_F和Φ_R进行采样的和也相应增加，并且的增幅远大于显然，此实验结果与前述关于Φ_Fx和Φ_Rx计算开销的理论分析结果相一致。再由图3(b)可知基于Φ_F和Φ_R的观测均可精确重构出原信号，这表明利用硬件实现简单、存储量小的托普利兹矩阵进行低计算开销的滤波采样，可获得与随机采样相同的重构性能。针对不同R，图3(d)显示当K＝10,M≥4K时，基于Φ_F和Φ_R的和均很小，而当M<4K时，和均很大。这说明M小即采样过少时，无法获得原信号的大部分信息，因此无法精确重构出原信号。图3表明了基于滤波的CS采样具有可行性、实用性。

2、可压缩信号的重构

为了测试基于Φ_S进行CS采样的可行性，采用语音信号作为实验信号，其中随机元素b_l服从随机均匀分布。实验环境为安静环境，实验对象为中国科学院自动化所录制的4位说话人语音，男性两位，女性两位，采样率为16kHz。实验针对测试语音分别为女声和男声的“涡轮风扇发动机”，采用矩形窗分帧，帧间无重叠，帧长为160样点/帧，定义重构语音信号平均分段信噪比SegSNR(dB)为：

$\mathrm{SegSNR}=\frac{1}{\mathrm{Nframe}}\underset{i=1}{\overset{\mathrm{Nframe}}{\mathrm{\&Sigma;}}}10\&times;\lg \left(\frac{x_i^T{x}_i}{{(x_i-x_i^{*})}^T(x_i-x_i^{*})}\right)---\left(25\right)$

其中Nframe为原信号的总帧数。将重构信号按信噪比的大小进行划分，其中信噪比不大于0dB、0～10dB、10～20dB、20～50dB、大于50dB的信号分别划分为重构质量很差、差、中等、好、很好。语音信号在DCT域具有近似稀疏性，因此选DCT基为Ψ^[16]。当压缩比R＝1/4时，图4中基于Φ_S和Φ_R的SegSNR分别为女声：29.0811dB、14.9225dB，男声：46.0058dB、30.5197dB，这说明基于Φ_S的语音信号重构效果均优于Φ_R，其原因为女声和男声下Φ_S与DCT基的相干度2.4918、2.5281小于Φ_R与DCT基的相干度3.8996、3.8272。相干度越小，越易利用式(17)搜索到稀疏解α^*而获得重构信号x^*＝Ψα^*，以弥补语音信号的DCT分解系数α为近似稀疏的不足。另一方面，女声的DCT分解系数高频成份较多，即不如男声系数集中于低频，使得其近似稀疏性弱于男声，因此重构效果不如男声，但基于Φ_S的重构效果仍明显优于Φ_R。

图5着重研究了不同R和N下基于Φ_S和Φ_R的平均采样时间和SegSNR，测试语音为由上述4位说话人分别随机选取2句短语，共8句短语组成的一长段语音，采用矩形窗分帧，帧间无重叠。由图5易见随着R和N的增大，和也相应增加，并且的增幅远大于选取帧长为160样点/帧，图5(a)(b)研究R对信号重构性能的影响。当R越大，即观测个数M越多，则SegSNR越高，Φ_R下信号重构效果越好；而Φ_S下信号重构误差有波动，其原因为Φ_S下要求N＝M×L,M＝N×R，若不满足条件，例如N×R不为整数，则Φ_S最后一行元素个数小于L，捕捉原信号的信息量减少，从而影响了重构效果。为了避免此问题，构造Φ_S时需根据R选择合适的N、M和L。针对相同的R，比较Φ_S和Φ_R下SegSNR可知：即使Φ_S下SegSNR有波动，但是均为重构质量好(20～50dB)，且大于Φ_R下的SegSNR。当R≤1/4时，Φ_R下重构质量从中等(10～20dB)迅速下降为差(0～10dB)，这反映了Φ_R下重构质量对采样个数M较敏感，而经过M和L的配合，Φ_S可通过局部采样捕捉信号的全局信息，使得重构质量受M的影响小，方便应用于实际CS采样。图5(c)(d)给出了R＝1/4时，上述信号在Φ_S和Φ_R下针对不同信号长度N(即不同帧长)的SegSNR。根据R值，N均取4的倍数，图5(d)显示Φ_S下SegSNR趋于稳定，避免了图3中重构误差波动的问题。图4和图5的实验结果表明基于Φ_S的CS采样是切实可行的，能够保证较好的重构精度和性能，优于基于Φ_R的CS采样。

本发明针对压缩感知中随机采样硬件实现、存储和计算困难的问题，依据有限脉冲响应滤波器原理，构造了一种新的基于滤波的压缩感知信号采集方案，并构造了相应的托普利兹测量矩阵，新方案降低了CS采样的硬件实现难度和计算开销。本发明理论证明了托普利兹测量矩阵以高概率满足RIP，易于通过l₁最优化问题求得最优稀疏解；同时，实验验证了基于托普利兹测量矩阵的CS采样在采样计算开销和重构性能方面都优于随机矩阵，具有有效性、可行性和应用性。

参考文献

[1]DONOHO D L.Compressed sensing[J].IEEE Trans.on Information Theory,2006,52(4):1289-1306.

[2]E.Compressive sampling[C].International Congress of Mathematicians.European Mathematical Society PublishingHouse,Madrid,Spain,2006:1433-1452.

[3]余恺,李元实.基于压缩感知的声音信号采集方法.仪器仪表学报,2012,33(1):105-112.

YU K,LI Y SH.New method for acoustic signal collection based on compressed sampling[J].Chinese Journal of ScientificInstrument,2012,33(1):105-112.

[4]CHEN W,WASSELL I J.Energy-efficient signal acquisition in wireless sensor networks:a compressive sensing framework[J].IET Wireless Sensor Systems,2012,2(1):1-8.

[5]FANG LY,LI S T,NIE Q.Sparsity based denoising of spectral domain optical coherence tomography images[J].BiomedicalOptics Express,2012,3(5):927-942.

[6]SUN X L,YU A X,DONG Z.Three-dimensional SAR focusing via compressive sensing:the case study of angel stadium[J].IEEE Geoscience and Remote Sensing Letters,2012,9(4):759-763.

[7]顾彬,杨震,胡海峰.基于序贯压缩感知的自适应宽带频谱检测.仪器仪表学报,2011,32(6):1272-1277.

GU B,YANG Z,HU H F.Adaptive wide-band spectrum detection based on sequential compressed sensing[J].Chinese Journal ofScientific Instrument,2011,32(6):1272-1277.

[8]CANDES E,ROMBERG J,TAO T.Robust uncertainty principles:exact signal reconstruction from highly in complete frequencyinformation[J].IEEE Trans.on Information Theory,2006,52(4):489-509.

[9]RONALD A,DEVORE.Deterministic constructions of compressed sensing matrices[J].Journal of Complexity,2007,23(4):918-925.

[10]DO T T,TRANY T D,GAN L.Fast compressive sampling with structurally random matrices[J].IEEE International ConferenceAgouties,Speech and Signal Processing,Washington D.C.,USA,2008:3369-3372.

[11]HAUPT J,BAJWA W U,RAZ G.Toeplitz compressed sensing matrices with applications to sparse channel estimation[J].IEEETrans.on Information Theory,2010,56(11):5862-5875.

[12]WANG K,LIU Y L.RIP analysis for quasi-toeplitz CS matrices[C].International Conference on Future Information Technologyand Management Engineering,2010,2:223-226.

[13]MARCO D,MARK D.Single pixel imaging via compressive sampling[J].IEEE Signal Processing Magazine,2008,25(2):53-92.

[14]HEGDE C,WAKIN M,BARANIUK R.Random projections for manifold learning[C].Neural Information Processing Systems(NIPS),Vancouver,Canada,2007:1-8.

[15]MONA S,OLGIEA M.Compressed sensing DNA microarrays,Rice ECE Department Technical Report TREE 0706,2007.

[16]叶蕾,杨震.行阶梯观测矩阵、对偶仿射尺度内点重构算法下的语音压缩感知.电子学报,2012,40(3):429-434.

YE L,YANG Z.Compressed sensing of speech signal based on special row echelon measurement matrix and dual affine scalinginterior point reconstruction method[J].Acta Electronica Sinica,2012,40(3):429-434.

[17]WANG X Z,WEI E I.Random sampling using shannon interpolation and poisson summation formulae.http://arxiv.org/abs/0909.2292.

[18]LASKA J N,KIROLOS S,DUARTE M F.Theory and implementation of an analog-to-information converter using randomdemodulation[C].IEEE International Symposium on Circuits and Systems.New Orleans,USA.2007:1959-1962.

[19]TROPP J A,WAKIN M B,DUARTE M F.Random filters for compressive sampling and reconstruction,IEEE InternationalConference on Acoustics,Speech,and Signal Processing(ICASSP),Toulouse,France,2006:872-875.

[20]刘树棠译.信号与系统.西安交通大学出版社,1998.

LIU S T.Signals$Systems[M].Xianjiaotong University Press,1998.

[21]李新，何传江.矩阵理论及其应用.重庆大学出版社,2005.

LI X,HE C J.Matrix theory and its application[M].Chongqing University Press,2005.

Claims (1)

1.一种基于滤波的压缩感知信号采集方法，步骤包括：首先，用感知设备在独立的采样周期内对目标信号x(t)进行采集，并用A/D方式对信号进行数字量化；然后，对量化后的信号x(i)进行降维；最后，对降维后的信号进行重构；其中t为采样时刻，i为量化后的信号排序；

其特征是对量化后的信号进行降维，具体是对量化后的信号通过有限脉冲响应滤波器的差分方程其中h(0),…,h(L-1)为滤波器系数，设计基于滤波的压缩感知信号采集框架，构造如下托普利兹测量矩阵：

则观测 $y\left(i\right)=\underset{l=1}{\overset{L}{\&Sigma;}}{b}_{l}x(i+l-1),i=1,...,M,$ 其中b₁,…,b_L看作滤波器系数；

子矩阵Φ_FT的奇异值是格拉姆矩阵G(Φ_F,T)＝Φ′_FTΦ_FT特征值的算术根，验证G(Φ_F,T)的所有特征值λ_i∈(1-δ_K,1+δ_K),i＝1,…,T，则Φ_F满足RIP，并通过求解如下l₁最优化问题来重构原信号：

$\begin{array}{ccc}\underset{x}{min}\left|\right|x|{|}_1& s.t.& y=\mathrm{\&Phi;}x,\end{array}$

即通过线性规划方法来重构原信号，亦即BP算法；

针对实际可压缩信号，则修改Φ_F为如下形式：

${\mathrm{\&Phi;}}_S=\left[\begin{array}{cccccccccc}{b}_1& \mathrm{...}& b_L& 0& \mathrm{...}& \mathrm{...}& \mathrm{...}& \mathrm{...}& \mathrm{...}& 0\\ 0& \mathrm{...}& 0& b_1& \mathrm{...}& b_L& 0& \mathrm{...}& \mathrm{...}& 0\\ \mathrm{...}& & \mathrm{...}& & \mathrm{...}& & \mathrm{...}& & \mathrm{...}& \\ 0& \mathrm{...}& 0& \mathrm{...}& 0& \mathrm{...}& 0& b_1& \mathrm{...}& b_L\end{array}\right];$

如果信号在变换基矩阵Ψ上具有稀疏性，则通过求解如下l₁最优化问题，精确重构出原信号：

$\begin{array}{ccc}\underset{\mathrm{\&alpha;}}{\min }\left|\right|\mathrm{\&alpha;}|{|}_1& s.t.& y=\mathrm{\&Phi;}x=\mathrm{\&Phi;}\mathrm{\&Psi;}\mathrm{\&alpha;}=\mathrm{\&Xi;}\mathrm{\&alpha;},\end{array}$ 其中Φ与Ψ不相关，Ξ称为CS矩阵；

所述RIP为有限等距特性(Restricted Isometry Principle，RIP)；

上面式中，Φ为测量矩阵，Ψ为基矩阵，y为观测向量，α为稀疏表示系数，x为量化后的稀疏信号，δ_K为常数且δ_K∈(0,1)。