Generalised Approximate Message Passing for Non-I.I.D. Sparse Signals

TL;DR

This paper extends the GAMP algorithm to handle non-i.i.d. sparse signals by introducing non-uniform priors, enabling more flexible and accurate compressed sensing of structured signals.

Contribution

It introduces new sparse signal priors for GAMP that support non-i.i.d. signals, advancing model-based compressed sensing capabilities.

Findings

Supports non-i.i.d. sparse signal reconstruction

Enables model-based compressed sensing with GAMP

Improves flexibility over traditional i.i.d. assumptions

Abstract

Generalised approximate message passing (GAMP) is an approximate Bayesian estimation algorithm for signals observed through a linear transform with a possibly non-linear subsequent measurement model. By leveraging prior information about the observed signal, such as sparsity in a known dictionary, GAMP can for example reconstruct signals from under-determined measurements - known as compressed sensing. In the sparse signal setting, most existing signal priors for GAMP assume the input signal to have i.i.d. entries. Here we present sparse signal priors for GAMP to estimate non-i.d.d. signals through a non-uniform weighting of the input prior, for example allowing GAMP to support model-based compressed sensing.

Equations

$$\mathbf{y}=\mathbf{Ax}+\mathbf{e}$$

$$\mathbf{z}=\mathbf{Ax},\mathbf{y}=\mathbf{z}+\mathbf{e}$$

$$\Vert \mathbf{x}{\Vert }_0=k\ll n,k<m,$$

$$p(x_j;\tau ,{\mathbf{\theta }}_I)=(1-\tau )\delta (x_j)+\tau \varphi (x_j;{\mathbf{\theta }}_I)$$

$$p(x_j;\tau ,\mu ,\sigma )=(1-\tau )\delta (x_j)+\tau \mathcal{N}(x_j;\mu ,\sigma ),$$

$$p(x_j;\tau ,{\mathbf{\theta }}_I)=(1-w_j\tau )\delta (x_j)+w_j\tau \varphi (x_j;{\mathbf{\theta }}_I)$$

$$f_{{\bar{x}}_j}(s_j,r_j;{\mathbf{\theta }}_I)={\pi }_j^{\text{w}}(r_j,s_j,{\mathbf{\theta }}_I)\left(\frac{\frac{\bar{\theta }}{\tilde{\theta }}+\frac{r_j}{s_j}}{\frac{1}{\tilde{\theta }}+\frac{1}{s_j}}\right)$$

$$f_{{\tilde{x}}_j}(s_j,r_j;{\mathbf{\theta }}_I)={\pi }_j^{\text{w}}(r_j,s_j,{\mathbf{\theta }}_I)$$

$$\cdot \left(\frac{1}{\frac{1}{\tilde{\theta }}+\frac{1}{s_j}}+{\left(\frac{\frac{\bar{\theta }}{\tilde{\theta }}+\frac{r_j}{s_j}}{\frac{1}{\tilde{\theta }}+\frac{1}{s_j}}\right)}^2\right)-f_{{\bar{x}}_j}(s_j,r_j;{\mathbf{\theta }}_I{)}^2$$

$$w_j=e^{-25{\left(\frac{j}{n}-\frac{1}{2}\right)}^2},j=0,1,\dots ,n$$

Full text

1 Introduction

Generalised approximate message passing (GAMP) was introduced by [1] in [1, 2]. GAMP addresses the estimation of signals $\mathbf{x}$ observed through a linear transform as follows

[TABLE]

where $\mathbf{y}\in\mathbb{C}^{m},\mathbf{A}\in\mathbb{C}^{m\times n},\mathbf{x}\in\mathbb{C}^{n},\mathbf{e}\in\mathbb{C}^{m}$. Here we express $\mathbf{e}$ as an additive noise which is classical in linear measurement models. We define the intermediate measurements without noise:

[TABLE]

The noise term $\mathbf{e}$ does not have to be a strictly additive term independent of $\mathbf{x}$. It can more generally be a (possibly non-linear) separable measurement channel expressed through a p.d.f. of $\mathbf{y}$; $p(\mathbf{y}|\mathbf{z};\bm{\theta}_{O})=\prod_{l=0}^{m-1}p(y_{l}|[\mathbf{z}]_{l};[\bm{\theta}_{O}]_{l})$ where $[\bm{\theta}_{O}]_{l}$ represents the parameters of the channel.

GAMP can be used quite generally for estimation of signals from many different distributions of $\mathbf{x}$. Here we consider the compressed sensing setting (see [3]) where $m<n$, i.e. an under-determined system which may be solved if

[TABLE]

where the $\|\mathbf{x}\|_{0}$ operator counts the number of non-zero entries in $\mathbf{x}$. GAMP can solve (1) in the compressed sensing setting when a sparse prior can be imposed on $\mathbf{x}$ to model (3).

For the compressed sensing setting, the algorithm approximate message passing (AMP) was proposed to estimate $\mathbf{x}$ with an i.i.d. Laplacian prior and i.i.d additive white Gaussian noise $\mathbf{e}$ [4]. GAMP can be seen as a generalisation of AMP that allows for a wider range of probability distributions on the signals $\mathbf{x}$ and on the measurements $\mathbf{y}$ given $\mathbf{x}$.

Another prior that can be used to model sparse signals in Bayesian estimators such as GAMP is the so-called spike-and-slab model [5]. According to this model, each entry $x_{j}$ in $\mathbf{x},\ j=0,1,\ldots n$ is distributed according to a linear combination of a Dirac delta p.d.f $\delta(x_{j})$ and another p.d.f. $\phi(x_{j};\bm{\theta}_{I})$:

[TABLE]

The $\delta(x)$ models the fact that many of the entries $x_{i}$ are zero (by the sparsity of $\mathbf{x}$). The variable $0\leq\tau\leq 1$ controls the sparsity of $\mathbf{x}$, i.e. how likely the entries are to be zero. The function $\phi(x_{j};\bm{\theta}_{I})$ can be chosen to represent the p.d.f. of entries of $\mathbf{x}$ that are not zero; $\bm{\theta}_{I}$ represents the parameters of $\phi$.

One example of such a spike-and-slab prior is the Bernoulli-Gauss distribution [6]:

[TABLE]

where the p.d.f. representing the non-zero entries is the Gaussian distribution with mean $\mu$ and std. deviation $\sigma$ (corresponding to the parameters $\bm{\theta}_{I}$). The term Bernoulli-Gauss (BG) reflects the fact that each entry $x_{i}$ can be seen as the product of a Bernoulli random variable (values 0 or 1) and a Gaussian random variable.

Applying GAMP with the signal prior (5) assumes that the entries of $\mathbf{x}$ are i.i.d. In many cases, signals of interest exhibit additional structure that can be exploited to estimate them more accurately [7].

Many different algorithmic approaches to modelling and leveraging such signal structure can be taken. Schniter et al. have for example produced substantial results on incorporating the ability to learn the structure of non-i.i.d. priors into the GAMP framework, see e.g. [8]. We take a different approach here and propose a weighted spike-and-slab model that can be used to model non-i.i.d. signals. We present the corresponding GAMP equations derived for this prior model on $\mathbf{x}$ as well as results from numerical simulations that show the improvements in reconstruction capabilities that are achievable using such a structured prior.

2 Weighted-Prior GAMP

We re-state the uniform variance MMSE GAMP algorithm [2] in Algorithm 1 as given in [9] with some variable changes to match (1).

We stress that our proposed prior may as well be used with the non-uniform variants of GAMP that track individual variances. Due to limited space, a detailed explanation of the algorithms can be found in [9] and references therein.

In order to apply the GAMP algorithm with a specific entry-wise input prior $p(x_{j};\tau,\bm{\theta}_{I})$ and measurement channel function $p(\mathbf{y}|\mathbf{z};\bm{\theta}_{O})$, it is necessary to derive the posterior p.d.f. $p(x_{j}|\mathbf{y};s_{j},r_{j},[\bm{\theta}_{I}]_{j})$. The evaluation of input posterior and measurement channel functions is incorporated in the GAMP iterations in the form of $\bar{\mathbf{x}}_{t}$ and $\breve{x}_{t}$ in Algorithm 1, lines 11-12, respectively $\bar{\mathbf{z}}_{t}$ and $\tilde{\mathbf{z}}_{t}$ in lines 5-6. We refer to these functions, $f_{\bar{x}},f_{\tilde{x}}$ and $f_{\bar{z}},f_{\tilde{z}}$ as the GAMP input- and output channnel functions, respectively.

Channel functions for the BG input channel (5) can be found in [10, 11, 6] and for the additive white Gaussian noise output channel in [1, 2].

Here we propose a modified Bernoulli-Gauss input channel that supports non-uniform sparsity over the signal $\mathbf{x}$, i.e. different entries $x_{j}$ can have different probabilities of being zero. We extend the spike-and slab model (4) as follows [9, p. 15]:

[TABLE]

General input channel functions have been derived for this model [9, p. 16]. As one example of such a weighted spike-and-slab input channel, we have derived the following closed-form expressions for the weighted BG input channel functions [9, p. 26], i.e. where $\phi(x_{j};\bm{\theta}_{I})=\mathcal{N}(x_{j};\bar{\theta},\tilde{\theta})$:

[TABLE]

where the function $\pi_{j}^{\text{w}}$ is given in eqs. (3.46)-(3.50) in [9, p. 15] and $\bar{\theta}$ and $\tilde{\theta}$ are the mean and the variance, respectively, of the Gaussian term in each entry of the GAMP input channel.

With the proposed weighted input prior (6) it is possible to model signals with varying probability of zero entries across the signal $\mathbf{x}$. Additionally, similar to the mechanism for the i.i.d. BG input channel in [6], we have derived formulas for updating the parameters $\tau,\bm{\theta}_{I},$ and $\bm{\theta}_{O}$ using expectation-maximization (EM) as part of the GAMP algorithm [9, sec. 6]. We retain separate parameters $\tau$ and $w_{j}$ to enable estimating the overall sparsity via EM without too many free parameters.

Next, we demonstrate by numerical experiments how the proposed model can improve estimation of compressively sensed signals with a known non-uniform sparsity structure.

3 Numerical Example

We perform a set of numerical experiments to demonstrate the benefits of the proposed weighted sparse prior for GAMP when such a model matches the signal of interest. We simulate the algorithm’s reconstruction capabilities in the form of a phase transition diagram over the full sparsity / under-sampling parameter space. See [12] for an introduction to the phase transition in compressed sensing problems.

We simulate sparse BG signals according to the proposed weighted model (6) where the weights $w_{j}$ are selected to have a Gaussian shape over the support of the signal vector $\mathbf{x}$:

[TABLE]

We partition the parameter space of $\delta=m/n$ and $\rho=k/m$ into a grid of points in each of which we reconstruct 10 random signals $\mathbf{x}$ with $k$ non-zero entries on average, measured with a randomly generated matrix $\mathbf{A}$ with i.i.d. random Gaussian entries. The reconstruction success (exact to within a small tolerance) of each signal is evaluated and the success rates in the $(\delta,\rho)$ points are used to estimate the phase transitions shown in Figure 1. The results shown here are a small subset of the simulation results available along with source code in [13].

For reference, the figure includes the theoretical $\ell_{1}$-optimisation phase transition curve [12], simulation results for a non-weighted BG prior (5) as well as simulation results for the implied Laplacian model used by AMP, “DMM AMP” [4]. We present results for a weighted prior that matches our Gaussian weights in (9) for: 1. The algorithm knows the true model parameters (genie) and 2. EM is used with a re-weighting scheme to estimate the parameters [9, sec. 6]. Comparing the two, we can see that knowing the underlying sparsity structure of the signal of interest improves the reconstruction capabilities substantially compared to simply assuming no weighting. Note, however, that it is important to allow the GAMP algorithm some “slack” in the form of estimating the parameters $\tau,\bar{\theta}_{j}$, and $\tilde{\theta}_{j}$ using EM to get the best performance. The error bars shown on the result curves correspond to the 10% to 90% percentile range of the logistic sigmoid functions fitted to the reconstruction outcomes to estimate the phase transition.

4 Conclusion

We have proposed a model for a class of non-uniformly structured sparse signals for use in the generalised approximate message passing (GAMP) algorithm. The proposed approach models sparse signals where the probability of zero vs. non-zero entries can vary across the support of the signal of interest. We have demonstrated through numerical examples how exploiting such structure when present in signals can substantially improve reconstruction of compressively sensed signals.