The Effect of Spatial Coupling on Compressive Sensing

Shrinivas Kudekar, Henry D. Pfister

I Introduction

This work investigates the effect of spatial coupling in compressed sensing. Spatially-coupled codes are a class of protograph-based low-density parity-check (LDPC) codes capable of achieving near capacity performance, under low-complexity belief propagation (BP) decoding, when transmitting over binary-input memoryless symmetric (BMS) channels. The history of these codes can be traced back to the work of Felström and Zigangirov , where they were introduced as convolutional LDPC code ensembles. There is a considerable literature on convolutional-like LDPC ensembles. Variations in their constructions as well as some analysis can be found in .

The fundamental reason underlying the remarkable performance was recently discussed in detail in for the case where the transmission takes place over the binary erasure channel. Before we go further, we briefly explain the construction of spatially-coupled codes; for more details, see .

Recall that a regular (dl,dr)(d_{l},d_{r}) LDPC code ensemble can be represented by the protograph (or base graph) as shown in Fig. 1. Spatially coupled code ensemble, denoted by (dl,dr,L)(d_{l},d_{r},L), is constructed by considering a protograph created by taking multiple copies of the (dl,dr)(d_{l},d_{r}) protograph (see the figure on the right-hand-side in Fig. 1) and connecting them as shown in the Fig. 2. We stress here that this is only the protograph. The code is constructed by taking multiple copies of this base graph and interconnecting them using a random permutation.

There are two main effects resulting from this coupling (see for details): (i) Rate Reduction: Recall that the design rate of the underlying standard (dl,dr=kdl)(d_{l},d_{r}=kd_{l})-regular ensemble is 1−dldr=k−1k1-\frac{d_{l}}{d_{r}}=\frac{k-1}{k}. The design rate of the corresponding (dl,dr=kdl,L)(d_{l},d_{r}=kd_{l},L) ensemble, due to boundary effects, is reduced to

(ii) Threshold Increase: Let ϵBP(dl,dr,L)\epsilon^{\text{\tiny BP}}(d_{l},d_{r},L), ϵMAP(dl,dr,L)\epsilon^{\text{\tiny MAP}}(d_{l},d_{r},L) and ϵMAP(dl,dr)\epsilon^{\text{\tiny MAP}}(d_{l},d_{r}) denote the thresholdThe BP(MAP) threshold, of a fixed ensemble of codes, denotes the channel parameter value below which the BP(MAP) decoder succeeds and fails above it. of the BP decoder for the (dl,dr,L)(d_{l},d_{r},L) ensemble, MAP threshold of the (dl,dr,L)(d_{l},d_{r},L) ensemble and the MAP threshold of the underlying (dl,dr)(d_{l},d_{r}) ensemble, respectively. Then the main result of is that, when we transmit spatially coupled codes over the BEC we have,

See for a precise statement of the main theorem. The effect of coupling can been nicely seen by plotting the EXIT curves for the uncoupled and coupled codes. This is shown in Fig. 3. Similar phenomena can also be empirically observed when transmitting over more general BMS channels .

I-B Outline

In this work, we study the effect of spatial coupling in the problem of compressed sensing. We begin with our compressive sensing setup and explain our decoders in the next section. In the same section, we introduce spatially-coupled measurement matrices. We then develop the density evolution (DE) equations for the class of decoders which we consider. In Section III we perform experiments depicting the effect of spatial coupling. We conclude with a short discussion of interesting open questions.

II Compressed Sensing

Compressed sensing (CS) is now one of the most exciting new areas in signal processing. It is based on the idea that many real-world signals (e.g., those sparse in some transform domain) can be reconstructed from a small number of linear measurements. This idea originated in the areas of statistics and signal processing , but is also quite related to previous work in computer science and applied mathematics . CS is also very closely related to error-correcting codes, and can be seen as source coding using linear codes over real numbers .

The main goal of this paper is to study the effect of spatial coupling on the performance of CS systems. Spatial coupling has been shown to drastically increase the threshold LDPC codes with an asymptotically negligible decrease in the code rate. The tight connection between low-density parity-check (LDPC) codes and CS means that one would expect spatial coupling to improve message-passing decoders for CS systems. A more subtle question is whether or not traditional CS decoders based on convex relaxations (e.g., basis pursuit and LASSO) will also benefit from spatial coupling. In this paper, we compare various standard constructions of measurement matrices with constructions based on spatial coupling. These tests are performed using a variety of CS decoders and a few signal models.

All our CS results rely on sparse measurement matrices. A variety of such matrices have been considered before in . Their main advantage is that they enable a variety of low-complexity reconstruction techniques. In general, the entries of the measurement matrix are either chosen to be plus/minus one with equal probability, or drawn from a continuous distribution. The latter provides some benefit when the signal contains only a small set of non-zero values.

Throughout the paper we will consider two kinds of measurement matrices. The first type of measurement matrix is generated from the parity-check matrix sampled uniformly at random from the ensemble of a regular (dl,dr=kdl)(d_{l},d_{r}=kd_{l}) LDPC codes. The sampling ratio, δ\delta, is given by 1/k1/k. The second type of measurement matrix, which we call a spatially coupled measurement matrix, will be generated from the parity-check matrix sampled from the (dl,dr,L)(d_{l},d_{r},L) ensemble. In this case, the sampling ratio can be obtained from (1) as δ=1k+dl−1kL\delta=\frac{1}{k}+\frac{d_{l}-1}{kL}. For finite LL, the number of measurements of the spatially coupled measurement matrix is larger than the corresponding uncoupled matrix. But as L→∞L\to\infty, the two measurement matrices have an equal number of measurements.

II-B Message Passing Reconstruction

Message-passing (MP) reconstruction for compressed-sensing systems based on sparse measurement matrices was introduced by Sarvotham, Baron, and Baraniuk in . The tight connection between CS and error-correcting codes enabled researchers to quickly analyze other MP reconstruction schemes based on tools from modern coding theory . Recently, there has also been some progress in analyzing the performance of these schemes for approximately sparse signals and noisy observations .

For sparse measurement matrices, the asymptotic performance of MP reconstruction can be analyzed (in theory) using density evolution (DE) . Indeed, this works well for simplified suboptimal reconstruction algorithms like the ”sudocodes” reconstruction . For true belief-propagation (BP) reconstruction, however, numerical evaluation of the DE equations is intractable. This means that it is difficult to determine the asymptotic behavior of a particular sparse ensemble with BP reconstruction. Recently, Donoho, Maleki and Montanari , proposed an approximate message passing (AMP) algorithm for compressed sensing with dense Gaussian measurement matrices. For this algorithm, they introduced a variant of DE (known as state evolution) that provides a precise characterization of its performance.

II-C Analytical Setup

Our analytical results rely on the suboptimal reconstruction technique introduced for ”sudocodes” by Sarvotham, Baron, and Baraniuk . The ”sudocodes” reconstruction technique falls into the class of verification decoders that was introduced and analyzed by Luby and Mitzenmacher for LDPC codes over large alphabets . In this paper, we use the message-passing based implementation of the second (more powerful) algorithm from their paper and refer to it as LM2 . The main drawback of this choice is that the analysis only works for strictly sparse vectors where the measurements are observed without noise. The main benefit is that one can analyze its performance precisely using density evolution (DE) and construct EXIT-like curves to illustrate the benefits of spatial coupling. For example, in the large system limit, this leads to provable sparsity thresholds where reconstruction succeeds with probability one .

Although the LM2 decoder is unstable in the presence of noise, this does not mean that its threshold is meaningless in practice. The LM2 decoder can be seen as a suboptimal version of list-message passing (LMP) which itself can be seen as a limiting case of the full belief-propagation (BP) decoder for CS . Ideally, one would analyze the BP decoder directly, but performing a DE analysis for decoders which pass real functions as messages remains intractable. Still, we expect that a complete analysis of the BP decoder would show that its expected performance is always better than the LM2 decoder and that the BP decoder allows stable reconstruction below its sparsity threshold.

Verification decoding rules for message-based LM2 decoding in a CS system: • (Check Node) The output message on an edge: – equals the unique value which satisfies the observation constraint based on all other input edges; – this output message is verified if and only if all other input messages are verified. • (Symbol Node) The output message is: – verified and equal to the input message if any input message is verified, – verified and equal to 0 if any input message on another edge is equal to 0, – verified and equal to the matching value if any two other input messages match, – and unverified with a value of 0 otherwise.

Initially, the check nodes with measurement equal to zero, transmit zero on all their outgoing edges. This arises from the basic property of verification decoders when we consider the non-zero values to come from a continuous distribution. The scheme described above does not guarantee that all verified symbols are actually correct. The event that a symbol is verified but incorrect is called false verification (FV). If either the non-zero matrix entries, or the signal values, are drawn from a continuous distribution, then the FV event has probability zero.

For generality, we consider a verification-based decoding rules for channels with erasures and errors. The variables w,x,y,zw,x,y,z will be used to denote, respectively, the probability that the message type is erasure (E), incorrect (I), correct and unverified (C), and verified (V). The DE equations for LM2 decoding of the standard irregular LDPC code ensemble are given in [43, p. 125]. For the protograph setup (including erasures), we derive the DE equations below.

For a check node of degree dd, let wi,xi,yi,ziw_{i},x_{i},y_{i},z_{i} the message-type probability for the iith input edge and wi′,xi′,yi′,zi′w_{i}^{\prime},x_{i}^{\prime},y_{i}^{\prime},z_{i}^{\prime} be the message-type probability for the iith output edge. Then, we have the update rules

In words, these four disjoint probabilities are for the events: “at least one E”, “at least one E or I minus at least E”, “all C or V minus all V”, and “all V”.

For a bit node of degree dd, let wi,xi,yi,ziw_{i},x_{i},y_{i},z_{i} the message-type probability for the iith input edge and wi′,xi′,yi′,zi′w_{i}^{\prime},x_{i}^{\prime},y_{i}^{\prime},z_{i}^{\prime} be the message-type probability for the iith output edge. If ϵ\epsilon is the probability of channel error and pp is the probability of channel erasure, then we have the update rules

Let AA be the event that “all input edges are E or I except for at most one C”. In words, these four disjoint probabilities are for the events: “channel erased and AA”, “channel error and AA”, “channel correct and all input edges are E or I”, and “not AA or channel correct and AA”.

We can specialize the above DE equations to the case when we do not have spatial coupling. More precisely, consider the measurement matrix generated by the regular (dl,dr)(d_{l},d_{r}) LDPC ensemble. Then for the check node we have,

These equations can be compared to [43, p. 125] by setting p=0p=0, w=0w=0, mapping x→bx\rightarrow b , and mapping y→ay\rightarrow a. We remark here that the channel error probability, ϵ\epsilon, is equal to the fraction of non-zero symbols in the signal.

II-C2 EXIT-like curves

In theory of iterative codes, the EXIT function is defined by H(Xi∣Y∼i)H(X_{i}|Y_{\sim i}) , where Y∼iY_{\sim i} denotes the vector of all observations except YiY_{i}. In words, it is the uncertainty in decoding a bit, when its channel observation is discarded. In general one can replace Y∼iY_{\sim i} by ϕi,dec(Y∼i)\phi_{i,\text{dec}}(Y_{\sim i}), where ϕi,dec(Y∼i)\phi_{i,\text{dec}}(Y_{\sim i}) is the extrinsic estimate of any other decoder (e.g., BP decoder). Such EXIT-like curves are useful visualization tool in iterative coding theory. They have also been used in iterative coding theory to provide deep results relating the BP decoder and the optimal MAP decoder . It is hence of great interest to visualize our results by plotting EXIT-like curves even in the case of compressed sensing.

According to the rules of LM2 decoder, a variable node is verified, or its value is perfectly known, when the variable node receives either two or more CC (correct but unverified) messages or one or more VV (verified) message. In this case the EXIT value of the variable node is zero. Thus the EXIT value of a variable node is proportional to the probability that a variable node is unverified. In other words, it receives “all I (incorrect) except at most one C (correct but unverified)”. This is exactly the probability of the event A (cf. Section II-C). To summarize: in our experiments we plot an EXIT-like curve which is the probability of a variable node being unverified when we change the sparsity ratio continuously.

III Experiments

We perform DE analysis of the verification decoder. More precisely, we consider the DE equations (II-C1) to (II-C1). In particular we consider the case when p=0p=0, i.e., we have corruption of symbols only via errors in the channel.

We first consider measurement matrices generated by the regular (4,8)(4,8) LDPC code. This fixes the sampling ratio to 1/21/2. We then fix the sparsity ratio 0<ϵ<10<\epsilon<1 and run DE equations (II-C1) and (II-C1), till a fixed point is reached. We then use the fixed point densities to evaluate the probability that a particular node is unverified. We plot this value on the vertical axis in Fig. 4, for different values of ϵ\epsilon. In this case, the EXIT-like curve is illustrated by the light gray curve (leftmost) in Fig. 4. We observe that below ϵ≈0.208\epsilon\approx 0.208, the fixed point is trivial. More precisely, for sparsity ratio ϵ≤0.208\epsilon\leq 0.208, the probability that a node is unverified goes to zero. This means that, with high probability, the LM2 decoder is able to reconstruct the signal exactly. For ϵ>0.208\epsilon>0.208, we see that the LM2 decoder fails to reconstruct the signal and the probability of a variable node being unverified is non-zero.

Similar experiment is now done with spatially coupled (4,8,L)(4,8,L) measurement matrices. We run the DE equations given by (II-C1) and (II-C1) for different lengths, L=2,4,8,16,32,64,128,256L=2,4,8,16,32,64,128,256. For L=2L=2, the threshold is ≈0.837\approx 0.837. The reason for such a large value of the threshold is the because the sampling ratio is much larger than 0.50.5 (cf. equation (1)). As LL increases the curves move to the left, which is similar to the effect observed in channel coding (cf. Fig. 3). As LL increases, the resultant measurement matrix resembles more and more like the uncoupled (4,8)(4,8) measurement matrix (locally) and the sampling ratio also approaches 1/21/2. However, the curves seem to get stuck at ϵ≈0.287\epsilon\approx 0.287. To summarize: for large LL we have the sampling ratio δ\delta very close to 0.50.5 and we observe that the sparsity threshold is much larger, ϵ≈0.287\epsilon\approx 0.287.

III-A2 Phase transition

We run the DE equations for different values of the sampling ratio. Let us explain this more precisely. As usual, we consider an uncoupled (dl,drd_{l},d_{r}) measurement matrix generated and its coupled version, (dl,dr,L)(d_{l},d_{r},L). We fix dl=4d_{l}=4 and consider dr=kdld_{r}=kd_{l} for different values of kk. As a consequence, we obtain different values of δ=1/k\delta=1/k. We run DE equations (II-C1) and (II-C1) for k=2,3,…,10k=2,3,\ldots,10 and for L=1000L=1000 fixed for all the experiments.

Figure 5 shows the results. The light circles correspond to the regular case and the dark circles correspond to the coupled case. We see that for the sampling ratios considered, the coupled measurement matrices have an improved sparsity-sampling trade-off. At the first, the trend of the circles may seem a bit strange. But, the normalized coordinates of the plot imply that a family of systems achieving a fixed oversampling ratio family would give a horizontal line.

The continuous curve shown is Donoho’s phase transition for weak CS reconstruction using LP . This was later identified in as the LP threshold for AMP reconstruction of a sparse signal with a worst-case amplitude distribution. For this reason, the comparison is not really fair. The LM2 thresholds are for the noiseless measurement of a special class of signals, while LP reconstruction thresholds are for the stable reconstruction of a random signal drawn from any i.i.d. distribution. Still, this result does highlight the fact that better performance is possible (and might even be achieved by LP reconstruction) in some special cases.

III-B BP and LP Reconstruction

Simulations were also performed for belief-propagation (BP) and linear-programming (LP) reconstruction of coupled and uncoupled measurement matrices. The initial goal was to choose system parameters that would allow direct comparison with [39, Fig. 3]. Unfortunately, the rate loss associated with coupled codes made it difficult to achieve fair comparisons in this case. Therefore, the tested system uses noiseless measurements to reconstruct signals from the two-Gaussian model with N=4032N=4032, ϵ=0.1\epsilon=0.1, σ0=0.5\sigma_{0}=0.5, and σ1=10\sigma_{1}=10. For the coupled and uncoupled systems, measurement matrices are based on (dl,dr)(d_{l},d_{r}) regular codes with dl=5,7,9d_{l}=5,7,9 and dr=2dl,3dl,4dld_{r}=2d_{l},3d_{l},4d_{l}. Non-zero entries in these matrices are randomly chosen to be either +1+1 or −1-1. These parameters give asymptotic sampling ratios of δ=0.25,0.5,0.75\delta=0.25,0.5,0.75, but the coupled systems have slightly more measurements due to finite-LL rate loss.

For each setup, the experiment tested a single (randomly constructed) measurement matrix on the same set of 10 randomly generated signal vectors with N=4032N=4032 and K=403K=403 (the number of non-zero symbols in the signal). During each trial, the BP decoder is run for 30 iterations and the root-mean-square error (RMSE) is calculated. The curves in Figure 6 show the median RMSE (over 10 trials) for the particular parameters.

In Fig. 6 we compare the coupled and uncoupled measurement matrix for increasing average degrees. In Fig. 7 we compare the performance of CSBP and LP decoders when we consider coupled and uncoupled measurement matrix with variable node degree fixed to 7. The figure also shows the performance curve when we use the Sarvotham, Barron, Baraniuk (SBB) (see ) ensemble for measuring.

We remark here that the simulation setup for this section uses LL between 24 and 48. The number of coupling stages must be kept somewhat low for two reasons: (i) to prevent short cycles in the decoding graph and (ii) to reduce the decoding time as the number of required decoding iterations increases with LL. These small values of LL result in a slightly larger sampling ratio (i.e., more measurements) for the coupled measurement matrices as opposed to uncoupled measurement matrices. This effect is handled correctly on the simulation curves, but makes direct comparison somewhat difficult. As a consequence, we observe that there is no appreciable gain by using couplng in the case when we are using the CSBP decoder. The performance changes are very small and do not support a conclusion that spatial-coupling provides large gains for CS with BP reconstruction.

The basis-pursuit LP reconstruction technique was also tested with each measurement matrix and signal vector. The results for symbol degree dl=7d_{l}=7 are shown in Figure 7. From this, one observes that LP reconstruction benefits even less from spatial coupling.

However, results from the channel coding setup imply that the performance for moderate LL is very close to the performance when LL is very large. For example, Figures 3 and 4 show that the threshold for L=16L=16 has already saturated to the L→∞L\rightarrow\infty threshold. Nevertheless, experiments with larger values of LL merit further investigation.

IV Discussion

Recently, it was shown in that the effect of coupling is also observed in many other problems, like the kk-satisfiability, graph coloring and Curie-Weiss model of statistical physics. In this paper, we found that spatially-coupled measurement matrices have an improved sparsity-sampling ratio for model used by verification decoding. This provides evidence that the phenomena of spatial coupling is quite general. On the other hand, we also observed that spatially-coupled measurement matrices provide little (if any) gains for the compressed-sensing problem with moderate blocklengths and belief-propagation reconstruction.

We conclude with possible future research directions.

There are also interesting open questions regarding EXIT-like curves for compressed-sensing reconstruction. What is the meaning, if there is any, of the limiting curve in Figure 4? Can we define an EXIT-like curve for CS reconstruction that obeys an Area Theorem like in the case of channel coding .

V Acknowledgments

SK acknowledges support of NMC via the NSF collaborative grant CCF-0829945 on “Harnessing Statistical Physics for Computing and Communications”. The work of HP was supported by the National Science Foundation under Grant No. CCF-0747470. The authors would also like to thank Misha Chertkov, Andrea Montanari, Rüdiger Urbanke for many useful discussions and encouragement.

References