Compressive Sensing with Optimal Sparsifying ... - Harvard University

Report 4 Downloads 119 Views
Compressive Sensing with Optimal Sparsifying Basis and Applications in Spectrum Sensing Youngjune Gwon, H. T. Kung, and Dario Vlah Harvard University Abstract—We describe a method of integrating KarhunenLo`eve Transform (KLT) into compressive sensing, which can as a result improve the compression ratio without affecting the accuracy of decoding. We present two complementary results: 1) by using KLT to find an optimal basis for decoding we can drastically reduce the number of measurements for compressive sensing used in applications such as radio spectrum analysis; 2) by using compressive sensing we can estimate and recover the KLT basis from compressive measurements of an input signal. In particular, we propose CS-KLT, an online estimation algorithm to cope with nonstationarity of wireless channels in reality. We validate our results with empirical data collected from a wideband UHF spectrum and field experiments to detect multiple radio transmitters, using software-defined radios.

I. I NTRODUCTION Sparse recovery techniques provide a new opportunity to build communication and data-driven systems. We consider compressive sensing (see, e.g., Cand`es and Tao [1]), a recent theory to extract critical information of data and represent it with substantially reduced measurements of the original data. The simple but effective framework of compressive sensing has enabled a broad range of applications in imaging, data processing in wireless sensor networks [2], cloud computing [3], cognitive radios [4], spectrum sensing [5], and IC Trojan detection [6]. Compressive sensing has three important properties. First, the encoding is blind to the content of a signal (or data) and has low computational complexity suitable for fast, real-time usage. Secondly, the number of measurements required for exact recovery is approximately proportional to sparsity of the signal, not its size. Lastly, the decoding is adaptive in the sense that the quality of recovered data can improve under a fixed number measurements—or equivalently, the required number of measurements that achieves the same quality can decrease— when a more effective sparsifying basis becomes available. Therefore, an algorithm to find an optimal basis would be critical for performance improvement. In this paper we propose a method to integrate the finding and the use of optimal sparsifying bases into compressive sensing to lower the number of measurements without affecting the accuracy of decoding. We demonstrate the benefit of our approach with an application system for spectrum sensing. The principle behind our method is Karhunen-Lo`eve Transform (KLT) [7], a classical procedure that reveals the correlation structure of a signal. KLT is optimal such that it decorrelates the signal into a representation comprising only statistically non-redundant coefficients. Thus we can use KLT to derive

sparsifying basis that can improve the sparse recovery process of compressive sensing. Despite the optimal sparsification, KLT has a well-known drawback for being data-dependent. That is, when the correlation structure of a signal changes, the KLT basis also has to be updated. This has been a prohibitive limitation for KLT in realistic settings, where the KLT basis is computed numerically from a sample covariance matrix of training data sets. It is contrasted to an approach based on a fixed support set such as Discrete Fourier Transform (DFT). The recovery with DFT as the sparsifying basis is suboptimal for a given signal, but it is popularly used in compressive sensing because it does not require any data-dependent adaptation. A novelty in our approach is to recover KLT bases from random projections of a signal, as we normally do to recover the signal with compressive sensing. This allows us to update the KLT basis incrementally from the measurements on the original signal with relatively few measurements. We summarize the main contributions of this paper concerning the use of optimal bases in compressive sensing: 1) we show a drastic reduction in required measurements when an optimal KLT basis, as opposed to DFT or DCT, is used for decoding; 2) we formulate an online algorithm to update the KLT basis directly from the measured input signal via compressive sensing, which only requires a nonadaptive, fixed support set for recovery (e.g., DFT basis); 3) we show how the previous two work together, using the KLT basis to recover the signal and the DFT to update the KLT basis. Although the latter requires more measurements, it needs not take place frequently; and 4) we demonstrate our approach lead to an overall gain in performance of spectrum analysis using multiple software-defined radios (USRPs) [8] and empirical data measured in the field. The rest of this paper is organized as follows. In Section II, we review compressive sensing and KLT. Section III introduces a compressive spectrum analysis system. In Section IV, we explain how to compute KLT basis with compressive sensing. Section V presents the CS-KLT online algorithm and discusses our implementation with USRPs. We evaluate the performance of our approach in Section VI and discuss related work in Section VII. Section VIII concludes the paper. II. BACKGROUND A. Compressive Sensing Compressive sensing exploits a sparse structure (either exposed or hidden) in data. Some data naturally exist as a

Secondly, the in-network measurement collection we consider depends heavily on computational costs. Conventional compression schemes likely require complex comp on platform. We aim simple compressive encoding schemes that can be performed wi nodes. 2

3.2

sparse form, and most of the others can be made sparse so that only several largest coefficients in some basis—which we call sparsifying basis—suffice to represent the data (e.g., Fourier transform of sinusoids). The fundamental premise of compressive sensing is that we can directly encode all significant coefficients without completely analyzing the data. More importantly, the encoding does not require any prior knowledge on sparsification of the data. Traditional compression schemes such as entropy coding, however, require statistical analysis of the data. Compressive sensing combines the measurement and compression of the data into one—non-analytical, lowcomplexity encoding based on matrix-vector multiplications using a randomly generated measurement matrix. Consider a vector x ∈ CN and an M × N measurement matrix Φ for M  N . Compressive sensing encodes x into compressive measurements y = Φx, an N to M reduction in size. The encoding takes place without the knowledge of a transform Ψ (an N × N matrix), which sparsifies x such that in s = Ψx there are only K  N nonzero elements (i.e., x is K-sparse). Compressive sensing forms an underdetermined system of equations susceptible to many solutions for x: min ksk`1

s∈RN

subject to y = Φx = (ΦΨ−1 ) s

(1)

The Robust Uncertainty Principle [9] states that M needs to N for some small constant c > 0 for exact be at least cK log K recovery. Also, the Restricted Isometry Property (RIP) [10] of Φ can guarantee a unique solution with high probability by the `1 -norm minimization decoding of s. Linear programming can be used to solve Eq. (1). Once s is recovered, x can be determined from x = Ψ−1 s. The quality of decoding depends on M . More interestingly, the data recovery is incremental—using a smaller M (than required for exact solution) does not disqualify the decoded result entirely. The corresponding largest components recovered (specific to M ) would still be approximately accurate. In other words, the majority of decoding error is contributed by the unrecovered components, and the sum of their magnitudes quantifies the approximate error. If desired accuracy is not met, one can increase M accordingly. Fixing M , however, by no means fixes the decoding accuracy because the use of a better sparsifying basis for the data can still improve the performance. Estimating an optimal basis and using it for decoding is the approach we take in this paper. B. Karhunen-Lo`eve Transform (KLT) But how can we systematically find an optimal basis? Consider x ∈ CN , a complex-valued, wide-sense stationary signal with mean zero (for simplicity). The covariance matrix of x can be computed numerically: Rx = E[xxH ], where the superscripted H denotes Hermitian transpose (i.e., xH = x∗T ). Rx is real and symmetric, and the eigen-decomposition, Rx = Q Λ QH , gives columns of Q the eigenvectors of Rx and Λ a diagonal matrix of the eigenvalues. Q is an orthogonal matrix, thus Q−1 = QH . The representation s = QH x is known as KarhunenLo`eve Transform (KLT) [7] of x, and we call Q the KLT

Mixing Multiple Compressive Measurements from Heterogeneous Nod =;51,/05!1-46+?!>/6@@3015!;0!%! :;4.!1/-.!1&-"(4) 0.2index Timesample (sample index ##) KLT !"#$%#&'()*!!+),-./0#)1&2#34) Time (sample index #) !"#$%#&'()*!!+),-./0#)1&2#34)

Fig. 2. Time domain representation of an actual wireless channel (source: 512 5/6.-0)78+)9-,1,)"#/"#,#&:-6;&)*=.->1&-"(4) 5/6.-0)78+)9-,1,)"#/"#,#&:-6;&)*=.->1&-"(4) Fig. 2. domain Timeof domain representation of an actual wireless (source: 512 domain representation an actual wireless channel (source: samples measured using USRP2 from UHF Ch. 21 channel (fc512 =channel 515 MHz, 6-MHz Fig. 2. Time representation of an actual wireless (source: 512 Fig. 4. fromThe samples from Fig. 2 represented using optim samples measured using USRP2 from UHF Ch. 21 (fc = 515 MHz, 6-MHz Fig. 4. The samples Fig. 2 represented using optimal KLT basis uredsamples using USRP2 from UHF Ch. 21 (f = 515 MHz, 6-MHz bandwidth) in Cambridge, MA on Saturday May 21, 2011 23:03 EST) c measured using USRP2 from UHF Ch. 21 (f = 515 MHz, 6-MHz Fig. samples 2 represented KLT using basis optim bandwidth) in Cambridge, MA onfor Saturday 21, 2011 Fig.from 4. The from using Fig. 2optimal represented Why Search BetterMay Basis? (2)c23:03 EST) Fig. 4. The samples Cambridge, MA Saturday May 21,Saturday 2011 23:03 bandwidth) in on Cambridge, MA on MayEST) 21, 2011 23:03 EST) Why Search for Better Basis? (2) 0

−0.05

−0.1

50

100

150

0 0.05

−0.05 0

0.2

−0.1 −0.05

0.1

50

200−0.1 250 50 300 Time (sample index #)

100

150

100350

150400

200 250 300 Time (sample index #)

350

200450 250500 300 Time (sample index #)

350

400

400

450

450

Why Search for Better Basis? (2) Frequency response (real)

500

500

Received signal represented in Fourier−basis (real)

0

−0.1

50

100

150

Optimal basis representation via KLT using averaging method (imaginary) Optimal basis representation via KLT using finite averaging method (imaginary) 0.1 finite 0.2

Magnitude Magnitude

0.05

Magnitude

Magnitude Magnitude

Magnitude

0.1

0 0.1

−0.1 0 −0.2 −0.1

50 100 150 200 the 250 smallest 300 350 400 450 500 We note that among all Ψ, Q−0.2makes s with Time (sample index #) !"#$%#&'()*!!+),-./0#)1&2#34) −0.3 a sample covariance matrix with 10 previous m −0.3 300 50 100 150 200 250 350 100 400 150 450 200 500 250 50 300 350 400 450 500 number of nonzeros. Fig.!"#$%#&'()*!!+),-./0#)1&2#34) 4Timedepicts the KLT Time representation (sample index #) (sample measurement index #) !"#$%#&'()*!!+),-./0#)1&2#34) a sample covariance matrix with 10 previous a sample covariance matrix with 10 F, previous sets. Compared to 2.the fixed support we obs of the same samples from Fig. We estimated Q from sets. Compared to the fixed support F, we observed a 14sets. Compared to the fixed support F, we ob a sample covariance with 10for previous measurementSuch fold gainmatrix in sparsity this example. spar !"#$%#&'()*!!+),-./0#)1&2#34) Frequency (FFT sample index #) sets.inCompared the this fixed support for F, Such we observed a 14fold gain sparsity for example. sparsity gain canspar fold to gain in sparsity this example. Such !"#$%#&'()*!!+),-./0#)1&2#34) significantly improve the accuracy ofcan recovery w !"#$%#&'()*!!+),-./0#)1&2#34) fold gain in sparsity for this example. Such sparsity gain significantly improve the accuracy of recovery with the same Fig. 3. The samples from Fig. 2 represented using fixed support F significantly improve the accuracy of recovery w Fig. 3. The samples from Fig. 2 represented using fixed support significantly F (orimprove even the fewer) number of compressive accuracy of recovery with the same measurem The samples from Fig. 2 represented using fixed support F (or even fewer) number of compressive measurements. Fig. 3. The samples from Fig. 2 represented using fixed support (or F even fewer) (or even numbermeasurements. of compressive measure numberfewer) of compressive IV. RKLT ECOVERING KLT BASIS WITH C OMPR fnk = e−j2πnk/N for n, k!"#$%#&'()*!!+),-./0#)1&2#34) = 0, 1, . . . , N − 1. F is invertible, R ECOVERINGKLT BASIS C OMPRESSIVE !"#$%#&'()*!!+),-./0#)1&2#34) IV. RIV. ECOVERING B ASISWITH WITH C OMPRESSIVE −1 −j2πnk/N !"#$%#&'()*!!+),-./0#)1&2#34) IV. R ECOVERING KLT BASIS WITH C OMP so we have x = F X, the Inverse DFT (IDFT). M EASUREMENTS f = e for n, k = 0, 1, . . . , N − 1. F is invertible, M EASUREMENTS nk πnk/N −j2πnk/N M EASUREMENTS for n, k = 0, 1, . . . , N − 1. F is invertible, −1 Usingx the DFT F, we can. .rewrite the M EASUREMENTS fso = ehave for basis n, kthe = Inverse 0, 1, .DFT ,N − 1.compressive F is invertible, nkwe =F X, (IDFT). We have mentioned earlier that compressive We have mentioned earliermeasurements that compressive m measurements in Eq.DFT (1): (IDFT). −1 x so = we F −1have X, the Inverse x = F X, the Inverse DFT (IDFT). can be used to recover the KLT matrix. To explain this,measurements we start We have mentioned earlier that compressive Using the DFT basis F, we can rewrite the compressive canWe have mentioned earlier compressive be usedthe to recover theencoding KLT that matrix. To explainm −1the compressive −1 e DFT basis F, we can rewrite with the tie-in between compressive at sensor =basis Φx = F, Φ(FweX)can = (ΦF )X the compressive (2)be used to recover the KLT matrix. To explain this, Using the DFT rewrite can we start measurements inyEq. (1): can used to recover KLT matrix. the tie-in between the encodi nodes and with KLT be basis estimation. Recallthe that thecompressive KLT matrixTo explain ntsmeasurements in Eq. (1): in Eq. (1): with the tie-in with between thecovariance compressive encoding at sensor The key intuition here is that from the matrix of input signal the tie-in between the compressive encod −1 a fixed support −1 F is used Q is computed nodesH and KLT basis estimation. Recall that the y −1 = Φxbasis = Φ(F X) = (ΦF )X (2) −1 nodes and KLT basis estimation. Recall that the KLT matrix as a sparsifying Ψ for compressive sensing decoding. x, R = E[xx ]. Similarly, the covariance matrix of the x −1)X −1 y = Φx = Φ(F X) = (ΦF (2) nodes and KLT basis estimation. Recall that the is computed the Hcovariance matrix of = Φx = Φ(F (ΦF )X can compressive (2) Qmeasurements With theymatrix product ΦF −1 , X) the `= y is from Ry = E[yy ]. Byofcompres1 -minimum decoder Q is computed from the covariance matrix input signal H QHR is from the covariance H Tcovariance The key intuition that a fixed support F is sive used x, = E[xx ]. R Similarly, m recover thenahere xfixed fromis X by IDFT. compressive encoding yx =computed Φx, we know: Φ ] =matrix of y = E[Φxxthe ntuition here is X, that support FWeis use used x, R = E[xx ]. Similarly, the H T covariance matrix of the H T x The key intuition here is that a fixed support F is used x,. R E[xx covariance sensing with basis F as the approach against which our ΦE[xx ]Φ So,x R= ΦRx Φ]. . Similarly, Note that Φ the is ynot aE[yyH ]. m as a sparsifying Ψbaseline for compressive sensing decoding. y = compressive measurements y is R = fying for compressive sensing decoding. T † H ]. By comprescompressive measurements ymeasurements is Ry =(Φ E[yy H newΨ approach will beΨ compared. square Using the pseudo-inverse )y , we as abasis sparsifying basis for−1 compressive sensing decoding. compressive is can Ryhave =yE[yy ]. With the matrix , the �1 -minimum decoder can matrix. sive encoding y = Φx, we know: R = E[Φ −1product ΦF H T atrix product ΦF product , the �1 -minimum decoder can decoder thecan following expression: −1 sive encoding y = Φx, we know: R = E[Φxx Φ ] = y With the matrix ΦF , the � -minimum H T T 1 sive encoding y =R Φx, we know: Ry = that E[Φ recoverC. X, then x Fourier from X by IDFT. We use compressive ΦE[xx ]Φ T. So, HowX Good Basis? Ty = ΦRx Φ . Note † then x from byxIsIDFT. We useIDFT. compressive ΦE[xxH ]ΦT . ΦE[xx So, RyHR]Φ = . Note that ΦΦT(4) isNote notTathat TΦR (Φ ) = ΦR xΦ y x recover X, then from X by We use compressive . So, R = ΦR . sensing with as the the baseline approach againstof which our square matrix. Using the y pseudo-inverse (Φ )† , Fig. 2 Fillustrates time-domain representation the T †x h sensing F as thewith baseline approach against which our square matrix. Using the pseudo-inverse (Φ ) , we can have T † Fwill as be the baseline approach which Here, our wethe findfollowing thatmatrix. we have been compressively measuring (Φ square Using the pseudo-inverse ) , new approach compared. measured from an actual UHF channelagainst at its Nyquist expression: T † chnew will approach besamples compared. the following expression: R in R (Φ ) , which can be approximated from y = Φx x y rate with will N =be512, using GNU Radio-USRP2 [8], [11]. compared. the following expression: x. Thus compressive meaT † Evidently, we Fourier cannot conclude any sparse structure in the that we used to encode ourT data C. How Good Is Basis? Ryto(Φ ) = ΦR † ood Is Fourier Basis?We take FFT and depict the frequency domain surements y have sufficient information recover R fromx R (Φ ) = ΦR (4) T † x y x time domain. C. How Good Is Fourier Basis? R (Φ ) = ΦR y (4).the Below is a procedure to estimate KLT basis Q withx Fig. 2 illustrates the time-domain representation of representation of the same samples in Fig. 3. About 23% of llustrates the time-domain representation of the Here,wewehave find that we have been measuring compressivel in four steps: Fig. the 2measured illustrates the time-domain ofcompressive the find measurements we that been compressively samples (≈ 117 outan of 512 samples) arerepresentation observed at to have samples from actual UHF channel itsHere, Nyquist T †−1 Here, we find that we have been compressive easured from an actual UHF channel at its Nyquist Ryy(Φ ) , which approximated TR†xXin Decode from = (ΦF )X usingcan fixedbe support F; relatively large from magnitudes that canUHF be saidchannel nonzero. atR samples measured an actual itsx[8], Nyquist in [11]. R1)y (Φ ) , which can from y = Φx fr rate with N = 512, using GNU Radio-USRP2 Tbe† approximated −1 R in R (Φ ) , which can be approximated N = 512, using GNU Radio-USRP2 [8], [11]. 2) Recover x by computing x = F X; x we used y that to encode our compressive data x. Thusmeacompf rate with Nwe=cannot 512, using GNU Radio-USRP2 [8],we [11]. that used to encode oursteps data x. Thus Evidently, conclude any sparse structure in the D. Compressive Sensing with Optimal Sparsifying Basis 3) Repeat the previous l times to numerically comwe cannot conclude any sparse structure in the that we used to our information data x. Thus Pencode surements y have to comp recov l sufficient 1 Evidently, wesupport cannot conclude any sparse structure in theypute: surements have xi xi H ;to recover Rx from Rx sufficient = E[xxH ] =information time domain. We take FFT and depict the frequency domain i=1 Fixed sets provide signal-independent sparsifying l n. We take FFT and depict the frequency domain surements y have sufficient information to reco H Below istoa estimate procedure to estimate 4)ofObtain by eigen-decomposition Rx KLT = QΛQ . time domain. takesame FFT and isdepict the Sparsity frequency domain (4). Below is (4). a Qprocedure basis Q KLT with b matrix Ψ We for but there a drawback. under23% representation ofdecoding, the samples in Fig. 3. About 6

Magnitude Magnitude Magnitude

Magnitude

4 2 0

−2 −4

50

100

2−2 0−4

−6 −2

2 0

−2

100

100

150

200 250 300 Time (sample index #)

350

400

150

450

500

450

500

4 0

2 −2 0 −4

50 100 150 200 250 300 Time (sample index #) 200 −4 250 300 350 400 450 500 50 100 150 200 250 300 Time (sample index #) Time (sample index #) −2

50

50

Received400 signal represented in Fourier−basis (imaginary) 200 −6 250 300 350 450 500 100 150 200 250 300 350 400 6 (sample 50 Time index #) Time (sample index #) Received signal represented in Fourier−basis (imaginary) 4 Received signal represented in Fourier−basis (imaginary) 62 150

Magnitude Magnitude

4

Magnitude

40

−4

6

−4

−0.3

4 Received signal represented in Fourier−basis (real) Received signal represented in Fourier−basis (real) 62

6

−6

−0.2

350

400

450

500

350

400

450

500

on of the asame samples inuniform Fig. 3.across About 23%subbands, of (4). Below is a procedure to estimate KLT b compressive measurements fixed (≈ support isout not different and23% compressive in four steps: in four steps: representation of117 the same samples in Fig. 3. About of measurements V. I MPLEMENTATION the samples of 512 samples) are observed to have s (≈ 117 out 512 samples) are take observed have (say, Mi compressive measurements in four −1 steps: as aof consequence, we should differenttonumber −1 1) Subsystems Decode X from = (ΦF )X usingF; fixed the samples (≈ 117 out of 512 samples) are nonzero. observed A.have Workflow 1) toDecode X and from y = (ΦF )Xyusing fixed large magnitudes that can besubband. said for subband i) can of measurements for each This partly −1 support argerelatively magnitudes that be said nonzero. −1 1) Decode X from y = (ΦF )X using 2) computing Recover x xby =F Fig. 5 outlines the workflow between a X; sensor x node and X; fixe relatively largethe magnitudes that KLT can basis be said nonzero. 2) Recover x by =computing F −1 motivates use of an optimal for each subband, 2) faithfully Recover x by computing xl sensing = F −1toX;nume the system. We implement a compressive D. Sensing Compressive Sensing withdifferent Optimal Sparsifying Basis 3) Repeat the previous steps times whichwith makesOptimal sparsity across subbands more uniform ssive Sparsifying Basis 3) Repeat the previous steps l times to numerically com� l times �E[xx any modification. are to 1 nodes H num D. Compressive Sensing for with Optimal Sparsifying and, more importantly decoding accuracy, also smaller.Basisencoder without 3) Repeat the steps l l OurH ]sensor 1= previous H Rx H pute: = x x ; i i support II.B, setsweprovide signal-independent sparsifying pute: Rx = inE[xx ] = Radio xH ; l1 �i=1 i xi framework the GNU software l for pport Fixed sets provide i=1 l H In Section signal-independent derived that the sparsifying KLT basis Q can give implemented pute: R = E[xx ] = x x H xdirect i x ;= Q 4) Obtain by eigen-decomposition Fixedthe sets provide signal-independent sparsifying i=1 l QΛQ USRP2. Q Theby nodes performQ bandpass [12] i .R 4) Obtain eigen-decomposition Rxsampling = matrix Ψsupport for there decoding, but there is drawback. under most compact for xa with s under = QH x.Sparsity In fact, or decoding, but is arepresentation drawback. Sparsity Obtain and convert the4)analog RF/IFQ intoby theeigen-decomposition complex digital in-phase Rx = Q matrix Ψsupport but there is a drawback. Sparsity isfor an decoding, ideal candidate for Ψ.across Similar to Eq. and (2), we write under a is fixed isacross not uniform different subbands, and port notQ uniform different subbands, and quadrature (I-Q) samples. Each sampling instance buffers V. I MPLEMENTATION compressive sensing with optimal KLT basis: V. I MPLEMENTATION aasfixed support not uniform across different subbands, and a we consequence, wedifferent should take different number (say, M i quence, should is take number (say, M uncompressed x containing N I-Q samples, and the encoder i V. I MPLEMENTATION as a consequence, we should take different number (say, M A. Workflow and Subsystems i uses a pre-generated, M × N Bernoulli random matrix Φ to y = Φx = Φ(Qs) = (ΦQ)s (3) for subband i) of measurements for each subband. This partly d i) of measurements for each subband. This partly A. Workflow and Subsystems A. Workflow and Subsystems subband of measurements each subband. partly motivates thei) use of an basis optimal KLT basis for eachThis subband, he for use of an optimal KLT for for each subband, outlines between the workflow between sens Fig. 5 outlinesFig. the 5workflow a sensor node aand motivates the use of an optimal KLT basis for each subband, which makes sparsity across different subbands more uniform Fig. 5 outlines the workflow between a sens es sparsity across different subbands more uniform the system. Wethefaithfully system. implement We faithfully implement asensing compres a compressive which makes sparsity across different subbands more uniform and, moreforimportantly for decoding smaller.without the any system. We faithfully implement a compre mportantly decoding accuracy, also accuracy, smaller. alsoencoder encoder without any Our modification. Our senso modification. sensor nodes are and, more importantly for decoding accuracy, also smaller. In we Section II.B, that the KLTgive basis implemented Q can give in encoder without any modification. Our senso n II.B, derived thatwethederived KLT basis Q can implemented in the software GNU Radio software fra the GNU Radio framework for

4 Compressive Sensing

(Encode) Buffer N x(t) (fc, B) sampled

xNx1

!MxN

Algorithm 1 The CS-KLT Online Optimal Basis Estimation

Transmit yMx1

Compressive Sensing

(Decode)

"=Q Compressive Sensing

(CS-KLT)

Fig. 5. We implement unmodified compressive sensing encoder but place an additional decoder block, namely CS-KLT, that estimates a KLT basis from the procedure described in Section IV.

produce M compressive measurements before transmitted innetwork. We simulated wireless uplinks to the base station for in-network transmissions. A novelty in our implementation is in the decoder. We added a new block for KLT basis estimation. It is labeled CSKLT in Fig. 5. CS-KLT is essentially a compressive sensing decoder that recovers an estimate of the covariance matrix from compressive measurements of the input signal. Based on the procedure described in Section IV, we formulated an online algorithm for CS-KLT, explained next. B. Online Algorithm for CS-KLT Algorithm 1 presents our pseudo-code implementation for the CS-KLT online optimal basis estimation. The algorithm aligns to incremental update or complete update duty cycles (tested by a boolean completeUpdate) to determine whether the sample covariance matrix Rx should be recalculated completely. This reestimate is necessary because in reality wireless channels are time-varying, and stationarity does not persist. Hence, we need to schedule a complete update once every certain number of measurement intervals. In the algorithm, procedure CS KLT invokes CS D ECODE R X, triggering the system to flush the current Rx and reacquire an estimate. On a normal measurement interval, procedure CS KLT invokes DO I NC U PDATE R X . The incremental update is based on weighted averaging (with weight α for 0 < α < 1) between the two successive intervals (see lines 20–22). C. Optimal Interval to Update KLT Basis The complete update of Rx may take up multiple measurement intervals and operate at low compression, costing much more measurements than recovering the input signal on a normal interval. So we should not do complete updates too frequently. An optimal interval for complete updates can be determined from the performance requirement (estimation/decoding errors), stationarity of the channel, and the relative cost of a complete update to an incremental update. VI. E VALUATION We evaluated the performance of our CS-KLT against Fourier basis in the following scenarios: (1) spectrum analysis of 200-MHz UHF whitespace; (2) sensing short radio transmissions. We collected all of the data presented in this section empirically from our field and laboratory experiments, using commercially available software radios.

1: procedure CS KLT(curRx, y, α, ΦCU , Φnormal , Q) 2: if completeU pdate == True then 3: newRx = CS D ECODE R X(y, ΦCU , F −1 ) 4: else 5: newRx = DO I NC U PDATE R X (curRx, y, α, Φnormal , Q) 6: end if 7: newQ = EIG(newRx) 8: return newQ 9: end procedure 10: procedure CS D ECODE R X(y, Φ, Ψ−1 ) 11: a ← MTX M ULTIPLY(Φ, Ψ−1 ) 12: b ← `1 - MIN D ECODE(y, a) 13: return b 14: end procedure 15: procedure DO I NC U PDATE R X(Rx , y, α, Φ, Ψ) 16: x ← `1 - MIN D ECODE(y, Φ, Ψ) 17: a ← CMPLX C ONJ(x) 18: b ← VTRT RANSPOSE(a) 19: c ← VTR M ULTIPLY(x, b) 20: d←α× c 21: e ← (1 − α) × Rx 22: f ←d+e 23: return f 24: end procedure 25: procedure CS E NCODE C H(f , bw, Φ) 26: (M, N ) ← SIZE O F(Φ) 27: TUNE RF(f ) 1 28: ts ← 2·bw 29: a ← BANDPASS S AMPLE(ts ) 30: b ← MTX M ULTIPLY(Φ, a) 31: return b 32: end procedure

A. Description of Experiments 1) UHF spectrum analysis: We prepared 4 USRP2 radios in an indoor lab. The spectrum is partitioned into J = 8 subbands (i.e., each with 25-MHz bandwidth) with center frequencies fc ∈ {512.5, 537.5, 562.5, 587.5, 612.5, 637.5, 662.5, 687.5} MHz. This spans UHF channels 19 to 51. We configured each USRP2 node to alternate between two subbands and measure. Duty cycles were counted in 1-msec unit measurement intervals. The Nyquist rate of each subband is 25 MHz × 2 = 50 million samples/sec (or 50,000 samples per unit measurement interval). We ran multiple sessions over different days at various hours, with each session lasting up to a minute. 2) Detecting Harris radios: This was conducted in an outdoor wireless test field. We used 4 USRP2 radios for compressive sensing of 100-MHz target spectrum partitioned in J = 4 subbands. The center frequencies for this experiment were {363, 386, 409, 432} MHz, meaning there is a 2-MHz overlap between two adjacent sensors. We placed 4 Harris military radios [13], two of which were engaged in live voice communications in the first subband. The other two were turned on and off randomly in the second and the fourth subbands. No radios transmitted in the third subband. The sessions lasted 10–30 seconds, and we kept other measurement parameters the same as before.

5 1

0.2 0.18

CS with Fourier basis CS−KLT

0.8

Detection error rate

Decoding error (L2 norm)

0.9

0.7 0.6 0.5 0.4 0.3

0.14 0.12 0.1 0.08 0.06 0.04

0.2

0.02

0.1 0

0.16

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

0 10

1

Fig. 6.

Error performance on recovering 200-MHz UHF spectrum

20

30

40

50

60

70

80

90

100

110

120

Number of compressive measurements (M)

Compression ratio (M/N)

Fig. 8.

Error performance of CS-KLT in a real-world experiment

Decoding errors (L2 norm)

0.8 0.7

CS Encode−Fourier Decode CS Encode−KLT Decode

0.6 0.5 0.4 0.3 0.2 0.1 0 0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

Compression ratio (M/N)

Fig. 7.

Error performance on sensing Harris radio transmissions

B. Spectrum Recovery Performance 1) UHF spectrum: Fig. 6 compares decoding errors of our CS-KLT against Fourier basis under various compression ratios (M/N ). For this analysis, we used N = 256 per subband (resulting 2048-point frequency response for the 200MHz spectrum). CS-KLT had 100-msec incremental and 1sec complete update duty cycles. We adopted a normalized `2 -norm error metric computed in the frequency domain, ˆ ` /||X||` (X ˆ indicates recovered frequency e`2 = ||X − X|| 2 2 response). Given a fixed error budget, our method resulted substantial, 2x to 4.3x gains in compression. For illustrative purposes, we plot the original frequency response and the recovered frequency responses from compressive measurements in Fig. 9 (Fourier) and in Fig. 10 (KLT). 2) Harris radio sensing: We also observe a clear advantage of CS-KLT in sensing Harris radio transmissions as indicated by the error performances in Fig. 7. The use of KLT basis for this case resulted about 2x to 3x gains in compression compared to using Fourier basis for decoding. C. Dynamic Identification of Trained Patterns via KLT We examined a preliminary application of detecting the ON/OFF state combinations of multiple transmitters using previously trained KLT matrices. For this evaluation, we reused the data collected from our Harris radio experiments. Suppose we wish to detect which of the target transmitters is in the ON state. To do this, we propose the following target detection strategy based on CS-KLT. First, we perform a training phase, where we estimate the KLT bases for individual Harris radios by activating them one at a time. Let us call the 4 KLT matrices obtained in our field experiment Q1 , . . . , Q4 . There could be other ways to obtain the KLT characterizations

of the targets requiring less access—for example, through training on reference signals or from known signal models such as the power-spectral densities of specific modulation schemes. Secondly, we perform a separate CS-KLT decoding operation, using each KLT matrix Qi as the sparse transform. Each decoding based on Qi detects the presence of target i, as desired. Clearly, this decoding process can detect any combination of active transmitters. Furthermore, since the KLT bases of individual transmitters do not change, there is no need for periodic updates. We computed the detection error rate by taking 5-sec long spectrum samples from each of our 4 USRP2s, encompassing 400 blocks of 2048 samples each. We added the samples together into a single stream of measurements and encoded them with compressive sensing. On the decoding, we recovered the combined samples with the KLT matrix Q3 and used a simple threshold to decide whether Harris radio 3 was ON or OFF. We present the detection error rate—the fraction of 400 CS-KLT instances in which the decoding decision was incorrect—in Fig. 8 versus the number of measurements (M ). We observe that even for small M , the error rate is quite low. This is an encouraging outcome, considering that the raw number of samples per measurement interval is close to 104 . VII. R ELATED W ORK Polo et al. [5] presented compressive spectrum sensing with distributed cognitive radios, which assumed a similar system model as ours. Wang et al. [14] proposed a compressive sensing decoder that combined analog-to-information converter (AIC) and the SOMP [15] algorithm. Duarte et al. [15] leveraged model-based joint sparsity of wireless signals to improve compressive sensing decoding. Our work is not just about applying classical transforms in the pursuit of optimal sparsifying bases for compressive spectrum sensing, but adaptation to time-variance and nonstationarity, which pose difficult yet most realistic problems. This differentiates our work from others. VIII. C ONCLUSION In this paper, we described how compressive sensing can benefit from using optimal sparsifying basis. Our approach was based on Karhunen-Lo`eve Transform (KLT) that maximally leverages the correlation structure of a signal. By using an

6

Uncompressed (M/N = 1)

M/N = .1 (Fourier)

1 0.5 0

M/N = .2 (Fourier)

1

1

0.5 550

600

650

700

0

0.5 550

M/N = .3 (Fourier)

600

650

700

0

1

1

1

0.5

0.5

550

600

650

700

0

550

M/N = .6 (Fourier)

600

650

700

0

1

1

0.5

0.5

600

650

Frequency (MHz)

700

0

550

600

650

Frequency (MHz)

700

0

Uncompressed (M/N = 1) 0.5

0.5

0.5

0

0

600

650

700

550

M/N = .3 (KLT)

600

650

700

0

700

600

650

700

1

1

0.5

0.5

650

700

0

550

M/N = .6 (KLT)

600

650

700

0

650

700

1

1

1

0.5

0.5

650

Fig. 10.

700

0

550

600

600

650

700

650

700

M/N = .8 (KLT)

0.5 600

550

M/N = .7 (KLT)

Frequency (MHz)

600

M/N = .5 (KLT)

1

600

550

M/N = .4 (KLT)

0.5

550

650

M/N = .2 (KLT) 1

0

550

M/N = .1 (KLT) 1

550

600

Frequency (MHz)

1

0

700

Recovering 200-MHz UHF spectrum with the Fourier basis under various compression ratios ( M ) N

Fig. 9.

550

650

M/N = .8 (Fourier)

1

550

550

M/N = .7 (Fourier)

0.5 0

600

M/N = .5 (Fourier)

0.5 0

550

M/N = .4 (Fourier)

650

Frequency (MHz)

700

0

550

600

Frequency (MHz)

) Recovering 200-MHz UHF spectrum with the KLT basis under various compression ratios ( M N

optimal KLT basis for decoding, we were able to drastically reduce the number of measurements without affecting the accuracy of data recovery. Our method, highlighted by the CSKLT online algorithm, was designed to work with time-varying wireless channels in reality, updating KLT basis incrementally from compressive measurements. We demonstrated empirical performance gains of our method in spectrum sensing applications. The results of this paper suggested that we can make an effective use of optimal bases in compressive sensing when the data dependency of KLT can be overcome with a proper online updating mechanism. This represents a significant departure from the current practice of compressive sensing applications that almost always uses a fixed support set such as DFT. While encouraged by our results, we realize the need for further work. It would be useful to develop a comprehensive theory on the efficiency limit of optimal basis compressive sensing. For the time being, we recover the covariance matrix with DFT basis, but other estimation method could be used for this purpose. Lastly, we plan to explore applications that leverage the new opportunity enabled by optimal basis compressive sensing. ACKNOWLEDGMENT This material is in part based on research sponsored by Air Force Research Laboratory under agreement number FA8750-10-2-0180. The U.S. Government is authorized to reproduce and distribute reprints for Governmental purposes notwithstanding any copyright notation thereon. The views and conclusions contained herein are those of the authors and should not be interpreted as necessarily representing the official policies or endorsements,

either expressed or implied, of Air Force Research Laboratory or the U.S. Government.

R EFERENCES [1] E. J. Cand`es and T. Tao, “Decoding by Linear Programming,” IEEE Trans. on Information Theory, vol. 51, no. 12, pp. 4203–4215, 2005. [2] C. Luo, F. Wu, J. Sun, and C. W. Chen, “Compressive data gathering for large-scale wireless sensor networks,” in ACM MOBICOM, 2009. [3] H. T. Kung, C.-K. Lin, and D. Vlah, “CloudSense: Continuous Finegrain Cloud Monitoring with Compressive Sensing,” in USENIX HotCloud, 2011. [4] Z. Tian and G. Giannakis, “Compressed Sensing for Wideband Cognitive Radios,” in IEEE ICASSP, 2007. [5] Y. Polo, Y. Wang, A. Pandharipande, and G. Leus, “Compressive Wideband Spectrum Sensing,” in IEEE ICASSP, 2009. [6] Y. Gwon, H. T. Kung, and D. Vlah, “DISTROY: Detecting IC Trojans with Compressive Measurements,” in USENIX HotSec, 2011. ¨ [7] K. Karhunen, “Uber Lineare Methoden in der Wahrscheinlichheitsrechnung,” PhD Dissertation, University of Helsinki, 1947. [8] The Universal Software Radio Peripheral (USRP) Products, Ettus Research LLC. [Online]. Available: http://www.ettus.com/ [9] E. J. Cand`es, J. Romberg, and T. Tao, “Robust Uncertainty Principles,” IEEE Trans. on Information Theory, vol. 52, no. 2, pp. 489–509, 2006. [10] E. J. Cand`es, “Compressive Sampling,” in International Congress of Mathematicians, 2006. [11] GNU Radio. [Online]. Available: http://gnuradio.org [12] W. Waters and B. Jarrett, “Bandpass Signal Sampling and Coherent Detection,” IEEE Trans. on Aerospace and Electronic Systems, vol. 18, no. 6, pp. 731–736, 1982. [13] Harris Radio. [Online]. Available: http://www.harris.com/ [14] Y. Wang, A. Pandharipande, Y. Polo, and G. Leus, “Distributed Compressive Wide-band Spectrum Sensing,” in ITA Workshop, 2009. [15] M. Duarte, S. Sarvotham, D. Baron, M. Wakin, and R. Baraniuk, “Distributed Compressed Sensing of Jointly Sparse Signals,” in Signals, Systems and Computers, 2005.