Sparse Signal Recovery from Fixed Low-rank Subspace via Compressive Measurement

This paper designs and evaluates a variant of CoSaMP algorithm, for recovering the sparse signal s from the compressive measurement () v Uw s  given a fixed low-rank subspace spanned by U. Instead of firstly recovering the full vector then separating the sparse part from the structured dense part, the proposed algorithm directly works on the compressive measurement to do the separation. We investigate the performance of the algorithm on both simulated data and video compressive sensing. The results show that for a fixed low-rank subspace and truly sparse signal the proposed algorithm could successfully recover the signal only from a few compressive sensing (CS) measurements, and it performs better than ordinary CoSaMP when the sparse signal is corrupted by additional Gaussian noise.


Introduction
In the last decade, from the pioneering theoretic foundation work [1][2][3][4], compressive sensing (CS) has been attracting many researchers' interests from various fields, such as signal processing [5], OPEN ACCESS medical imaging [6], sensor network [7,8], machine learning [9,10], etc. Regarding the algorithmic aspect of CS, the main problem is to solve the underdetermined linear equation system vx  efficiently with the constraint that x is sparse.There are several well-known approaches to solve the sparse recovery problem, such as the greedy solver CoSaMP [11], the classic convex optimization approach [12], etc.
In this paper, we consider a special noisy sparse recovery problem: Here  is the Gaussian noise with small variance added into data itself, and is the small Gaussian measurement noise.The sparse signal s is buried in structured dense signal l Uw  [13] generated by a fixed low dimensional subspace spanned by the column space of U. In the CS field, this special problem represents a large body of applications in which the measured object can be regarded as the superposition of a sparse part and a dense but structured part.For example, in video surveillance, each frame contains a slowly-changing background or the background is static, which can be faithfully modeled by a low-dimensional subspace; and the remaining part-foreground, which is always of interest, is sparse.
From the highly compressive measurement can we exactly recover the sparse signal s?In statistics and optimization, this problem can be regarded as a compressive sensing version of the classic least absolute deviations (LAD) problem [14] or 1 regression problem.Instead of firstly recovering the full vector then separating the sparse part from the structured dense part, we propose a variant of CoSaMP approach named CoSaMP_subspace, which directly works on the compressive measurement to separate the sparse signal from low-rank background signal.
This paper is organized as follows.Section 2 describes the model of the special sparse recovery problem, introduces the proposed CoSaMP_subspace algorithm, and discusses its relation to ordinary CS technique.In Section 3, we conduct extensive experiments on both simulated data and the real surveillance video.Section 4 concludes this paper and points out future work.

The Model
We denote the fixed d-dimensional subspace of n as .In application of interest we always suppose dn , say is a fixed low-rank subspace.Let the columns of an nd  matrix U be orthonormal and span .{ | , } nm mn  is a linear compressive measurement operator which should satisfy the restricted isometry property (RIP) [2] We want to recover the K-sparse signal n s  buried in the low-rank background signal l Uw  and small Gaussian noise corruption  from the compressive measurement v with small Gaussian measurement noise as Equation ( 1), here, K-sparse means 0 sK  ‖‖ .

Variant of CoSaMP for Fixed Subspace
To tackle this special sparse recovery problem, we adopt the basic idea of CoSaMP [11].Ordinary CoSaMP could efficiently solve the following underdetermined linear equation system if x is a K-sparse signal and the measurement operator satisfies RIP with constant K  .
() yx  Comparing the problems between Equations ( 1) and ( 3), for a fixed low-rank subspace, if we can remove the dense signal generated by the subspace, then problem Equation (1) will be converted to the classic sparse recovery problem Equation (3).We consider the greedy approach to iteratively remove the effect of the dense signal then estimate the sparse signal in the framework of CoSaMP.The main idea is as follows: we iteratively estimate the weights k w which constructs the low-rank background signal k l Uw  ; then form the signal proxy from the current residue () identify the support of the sparse signal just as ordinary CoSaMP; finally estimate the signal by least squares over the identified support.
We summarize the proposed CoSaMP_subspace as Algorithm 1.Here * denotes the adjoint of the operator .The notation supp(y;K) denotes the largest K-term support of the vector y.Let |T denote the restriction of the operator to the support T. Algorithm 1.Given a fixed subspace spanned by the column space of an nd  orthonormal matrix U, the variant of CoSaMP solver for the sparse recovery problem Equation ( 1).

3:
Estimate weights w: Signal estimation by least squares: It can be easily seen from Algorithm 1 that, compared with CoSaMP [11], the variant tries to estimate the weights of the dense signal in step 3 and remove its effects from the CS measurement in Step 8.The main framework follows CoSaMP.Thus, the proposed algorithm inherits the merits of CoSaMP, for example the total cost per iteration is O(nlogn) [11], and the total storage is, at worst, O(nd) because of holding the fixed low-rank subspace.

Relation to Ordinary CS
For this special form of sparse recovery problem Equation ( 1), it could be also tackled via ordinary CS algorithms.We discuss its relation to our proposed algorithm in this subsection.
We rewrite problem Equation (1) as: Here I is the nn  identity matrix.Let Then problem (1) converts to the classic CS model regardless of the effect of data noise  .
  As is a CS operator satisfying RIP, for example, a random Gaussian matrix, and fixed basis, it is well-known in CS theory that is largely incoherent with B [2].Then, x can be exactly recovered via CS techniques given enough linear measurements, for example CoSaMP, so does s.
However, we want to point out that though ordinary CS could also estimate the sparse signal s if we make the proper linear algebraic transform as above, the proposed algorithm works better when the sparse signal is also corrupted by small noise  and not only the measurement noise .It is because in step 3 of Algorithm 1 we explicitly make the least square estimation to cancel the noise effect.We put the comparisons in the experiments section.

Experiments Evaluation
In this section, we evaluate the performance of Algorithm 1 extensively on both simulated data and video CS application.For simulated data, we use Equation (1) to generate a series of CS vectors v. Like most CS literature, we use the noiselet operator [15,16] as the measurement operator , which compresses an 1 n  vector to a 1 p  vector, and we denote cs p n   as the CS measurement ratio.
U is an nd  matrix whose d columns are realizations of i.i.d.(0, ) n I random variables that are then orthornomalized.The weight vector w is a 1 d  vector whose entries are realizations of i.i.d.(0, ) n I random variables, which are Gaussian distributed with mean zero and variance 1.The K-sparse signal s is an 1 n  vector whose supports are chosen uniformly at random without replacement and we denote s K n   as the sparsity of the signal.We use relative error to quantify the sparse recovery performance as follows: Note that we manually add -1‖ to the denominator of RelErr because s may be 0-sparse in our following experiments.In all the following experiments, we use Matlab R2010b on a Macbook Pro laptop with a 2.3 GHz Intel Core i5 CPU and 8 GB RAM.We always set the rank of the fixed subspace as d = 5, and the ambient dimension of signal is n = 512, unless otherwise noted.

Algorithm Behavior on Simulated Data
As a special CS algorithm, given a fixed subspace the performance of Algorithm 1 mainly relates to the sparsity of the signal s  and the CS measurement ratio cs  .We will also investigate how the rank of fixed subspace affects the performance of Algorithm 1.

Recovery on Signal Sparsity
In the first experiment, we investigate how signal sparsity affects the success of recovery.Here we fix the CS measurement ratio cs  as 15%, 30%, and 50%, respectively, and manually tune the signal sparsity from 0% to 30%. Figure 1 shows the fact that for truly sparse signal s, even a very small ratio of CS measurements could recover s from the low-rank background signal.For example, a 25-sparse signal can be successfully recovered from only 150 CS measurements.

Recovery on CS Measurements
In the second experiment, we investigate how CS measurement ratio affects the success of recovery.
Here we fix the signal sparsity s  as 5%, 10%, and 20%, respectively, and manually tune the CS measurement ratio from 10% to 80%. Figure 2 shows the fact that, given enough CS measurements, even a near sparse signal can be precisely estimated, for example a 100-sparse signal can be recovered with high quality from 350 CS measurements

Recovery on the Rank of Fixed Subspace
As the special sparse recovery problem depends on the fixed subspace, in the third simulated experiment we demonstrate how the rank of the fixed subspace affects recovery.To observe the algorithm behavior, we fix the signal sparsity s  as 5% and generate a series of simulated signals with different rank setting from rank = 1 to rank = 100.Figure 3 shows the recovery results with CS measurement ratio cs  as 30%, 50%, and 70% respectively.It can be clearly seen that given moderate enough CS measurements the proposed algorithm can guarantee the recovery with up to a relative high rank of the fixed subspace.For example, with the CS measurement 50% cs   the sparse signal can be stably separated from the fixed subspace with up to rank = 40.

Comparisons with Ordinary CS
From the previous discussion we know that, if we transform the original problem Equation (1) to Equation ( 5), the sparse signal could also be recovered by ordinary CS.However, as we just pointed out above, if the noise not only occurs in measurement phase but also lies in data itself, Algorithm 1 can give a better estimation.Here, we address this comparison by adding relative strong noise corruption to the data in which the variance of  in Equation ( 1) is 3 10  .In addition to this noise setting, we fix the sparsity 5% s   and vary the CS measurement ratio cs  from 10% to 80% to observe the performance comparison between Algorithm 1 and ordinary CS.From Figure 4 it can be easily seen that Algorithm 1 outperforms ordinary CS (CoSaMP) when data itself is also corrupted by relative strong Gaussian noise, especially when the CS measurement ratio is small.We make a more clear comparison in Figure 5. Here, we simulate a sparse signal s with sparsity 5%

Video Compressive Sensing
In video surveillance, background images are always modeled as a low-rank subspace and foreground moving objects are regarded as sparse signals.Then, the following -Low-rank + Sparse‖ model [17,18] can well represent the video frame captured by stationary camera, which has been studied extensively in recent literature.
Here, we consider applying Algorithm 1 to video compressive sensing application [19], for example, the video frames are captured by single-pixel camera [20] in CS literature.We simulate video CS by performing noiselet operator [15] on each normal surveillance video frame.Then recovering the foreground moving objects is just the problem Equation (1) in this paper.We test our simulation on the two well-known datasets (Dataset can be downloaded from http://perception.i2r.a-star.edu.sg/bk_model/bk_index.html.)-airport‖ and -lobby‖ in which the image dimension of the two datasets is 144 × 176.In order to perform CS measurement on each video frame, conforming to the noiselet convention, we first manually resize each frame as 128 × 128, which is the power of 2, thus, the ambient dimension of each video frame t I in our experiment is 128 × 128.Then, for each video frame t I , we take the CS measurement () tt vI  with CS measurement ratio cs  .As the two surveillance videos are static, our fixed low-rank subspace assumption holds.We obtain the fixed low-rank orthonormal matrix U as follows.For each dataset, we choose 200 video frames and perform the Robust PCA algorithm Inexact ALM (IALM) [21] on those frames to get the clean background images, which can be regarded as lying in a low-rank subspace.

, and 1%
cs   , respectively.For -airport‖, Figure 6 shows that the moving foreground objects can be well estimated though there are some ghost-effects in the recovery images, for example the first column.This effect is because IALM can not totally remove foreground from background due to the objects moving slowly in the scene which does not strictly follow the PCP recovery theorem [17].The second row of Figure 6 shows the background estimated by IALM in which there does exist some slight shadows.

Conclusions and Future Works
This paper presents a variant of CoSaMP algorithm, which can tackle a special sparse signal recovery problem () v Uw s  .The algorithm could separate the sparse signal and low-rank background signal from their sum given very few CS measurements.From experiments we show that the recovery performance is similar to ordinary CoSaMP and performs better if the sparse signal also corrupted by Gaussian noise not only the measurement noise.This proposed algorithm assumes that the low-rank subspace is known as a prior.However, what if the subspace is not known or even the subspace is time-varying?It is a challenging problem, to simultaneously estimate the subspace and recover the sparse signal only from a few CS measurements, which is attracting the attention of CS community [22][23][24].We are very interested in incorporating the recent online subspace learning technique [18,25] into this CS framework.We put this endeavor for future work.

Figure 1 .
Figure 1.Sparse recovery performance with respect to signal sparsity.

Figure 2 .
Figure 2. Sparse recovery performance with respect to CS measurement ratio.

Figure 3 .
Figure 3. Sparse recovery performance with respect to the rank of fixed subspace.

Figure 4 .
Figure 4. Performance comparison between CoSaMP_subspace and ordinary CS (CoSaMP) for noisy sparse recovery problem.

. 10 
The variance of additional Gaussian noise  and measurement noise are all set as 3 .From this noisy CS measurement setting, Figure5shows that the estimated sparse signal by our algorithm can exactly match the original signal strong data noise  would introduce more uncertainty when the algorithm tries to identify the sparse support  .

Figure 5 .
Figure 5. Noisy sparse signal recovery comparison between CoSaMP_subspace and ordinary CS (CoSaMP).Top-left is the compressive measured signal v; top-right is the uncompressed signal--dense structured signal + sparse signal‖; bottom-left is the recovered sparse signal by our algorithm; bottom-right is the recovered sparse signal by ordinary CoSaMP.

Figure 6 .
Figure 6.Foreground recovery from video compressive sensing on -airport‖ dataset.The 1 st row is the original video frames; the 2 nd row shows the estimated background images by IALM[21] that are used to train our fixed low-rank subspace; the 3 rd row shows the recovered foreground with CS measurement ratio 30% cs  

Figure 7 .
Figure 7. Foreground recovery from video compressive sensing on -lobby‖ dataset.The 1 st row is original video frames; the 2 nd row shows the estimated background images by IALM [21] that are used to train our fixed low-rank subspace; the 3 rd row shows the recovered foreground with CS measurement ratio 20% cs   Then, the rank d orthonormal matrix U is obtained by performing SVD on those background images and keeping the columns corresponding to the largest d singular values.In our experiments we set d = 5.In order to show how CS measurement ratio affects the recovery results, for the moderate sparse dataset -airport‖ we set the sparsity 20%