PaperlayerView exact revision on arXiv

Source-linked AI summary

Statistical physics-based reconstruction in compressed sensing

Florent Krzakala, Marc Mézard, François Sausset, Yifan Sun, Lenka Zdeborová

arXiv:1109.4424v4cond-mat.stat-mechcs.IT

TL;DR

Compressed sensing aims to recover sparse signals from fewer measurements, but conventional practical methods require more than the signal density. The paper combines probabilistic reconstruction, EM-BP message passing, and crystal-nucleation-inspired seeding matrices, and reports reconstruction rates approaching α = ρ0 in large systems. The analysis uses replica and cavity methods, with numerical studies supporting the result.

  • Problem

    Practical compressed-sensing reconstruction methods require measurement rates larger than the signal density, although α = ρ0 is the theoretical limit for exact noiseless recovery.

  • Method

    The s-BP method combines a probabilistic signal model, EM-BP message passing with learned parameters, and structured seeding measurement matrices.

  • Results

    s-BP approaches the optimal measurement rate α = ρ0 in the large-N limit, with statistical-physics analysis and numerical studies supporting the improvement.

  • Takeaways & Limitations

    Seeding enables practical message-passing reconstruction beyond the threshold of unstructured EM-BP and close to the fundamental compression limit.

  • Takeaways & Limitations

    The reported analysis assumes a Gaussian nonzero-value model, and using it for binary signals yields nearly no improvement over ℓ1 reconstruction.

Abstract

from arXiv · show

Compressed sensing is triggering a major evolution in signal acquisition. It consists in sampling a sparse signal at low rate and later using computational power for its exact reconstruction, so that only the necessary information is measured. Currently used reconstruction techniques are, however, limited to acquisition rates larger than the true density of the signal. We design a new procedure which is able to reconstruct exactly the signal with a number of measurements that approaches the theoretical limit in the limit of large systems. It is based on the joint use of three essential ingredients: a probabilistic approach to signal reconstruction, a message-passing algorithm adapted from belief propagation, and a careful design of the measurement matrix inspired from the theory of crystal nucleation. The performance of this new algorithm is analyzed by statistical physics methods. The obtained improvement is confirmed by numerical studies of several cases.

Reconstruction in Compressed Sensing

Compressed sensing seeks exact recovery of sparse signals from fewer measurements, but practical ℓ1 methods require rates above the fundamental density limit. The s-BP approach combines probabilistic reconstruction, belief-propagation message passing, and seeded measurement matrices to approach that limit.

  • Reconstruction in Compressed Sensing: Exact reconstruction is theoretically possible at measurement rate α ≥ ρ0, but exhaustive enumeration is exponentially slow in N.Here ρ0 is the signal density and α = M/N is the measurement rate.
  • Reconstruction in Compressed Sensing: The ℓ1 phase transition satisfies αℓ1(ρ0) > ρ0, so ℓ1 reconstruction requires more measurements than the information-theoretic minimum.For α below αℓ1, exact recovery probability vanishes in the large-N limit.
  • Reconstruction in Compressed Sensing: s-BP combines a probabilistic signal model, belief-propagation message passing, and measurement matrices designed using seeding ideas from crystal nucleation.The paper attributes reaching the α = ρ0 limit to the combined use of all three ingredients.
  • Reconstruction in Compressed Sensing: Phase-diagram experiments show s-BP thresholds close to the optimal line α = ρ0, with error bars shrinking as N increases.The figure compares Gaussian and binary nonzero components and reports thresholds from 50% success, with error-bar endpoints at 10% and 90% success.
  • Reconstruction in Compressed Sensing: The probabilistic reconstruction uses a Gauss-Bernoulli measure restricted to vectors satisfying the measurement constraints, without requiring the true signal distribution.The assumed sparsity density and nonzero-value distribution may differ from the signal’s generating parameters, provided the assumed density is below one.
  • Reconstruction in Compressed Sensing: Statistical-physics analysis describes reconstruction through a free-entropy landscape, where seeding removes a metastable local maximum that traps ordinary EM-BP.For ρ0 = 0.4, increasing the seeding structure approaches α = ρ0 from α1 = 0.7.

Sampling with Expectation Maximization Belief Propagation

EM-BP approximates probabilistic sampling through message passing while learning the signal-model parameters. Its reconstruction threshold depends on the signal distribution, and its performance can substantially differ when the assumed and actual distributions mismatch.

  • Sampling with Expectation Maximization Belief Propagation: EM-BP approximates sampling from the constrained probabilistic measure using belief-propagation messages and iterated parameter updates.The learned parameters are the sparsity density ρ and the Gaussian mean x and variance σ2.
  • Sampling with Expectation Maximization Belief Propagation: The algorithm iterates messages and the three model parameters from random initialization until reaching a fixed point.Perfect reconstruction corresponds to messages satisfying ai→µ = si and vi→µ = 0.
  • Sampling with Expectation Maximization Belief Propagation: EM-BP achieves perfect reconstruction in the large-N limit when α exceeds a transition αEM−BP that depends on the signal distribution.This dependence contrasts with the distribution-independent ℓ1 transition.
  • Sampling with Expectation Maximization Belief Propagation: When the assumed Gaussian nonzero distribution matches the signal, EM-BP improves substantially over ℓ1; for binary signals, the improvement is nearly absent.The analysis considered a Gaussian assumed distribution in both cases.

Designing seeding matrices

Seeding matrices partition variables and measurements into coupled blocks whose measurement rates and variances are designed to initiate reconstruction in one block and propagate it through the rest.

  • Designing seeding matrices: The blockwise reconstruction is designed to eliminate the metastable state that traps EM-BP with unstructured measurement matrices.The mechanism is compared with crystal growth from a nucleated seed.
  • Designing seeding matrices: Seeding matrices divide variables and measurements into L groups and couple each pair of groups through variances Jq,p/N.Standard compressed-sensing matrices correspond to L = 1 and α1 = α.
  • Designing seeding matrices: A high measurement rate α1 in the first block seeds reconstruction, which then propagates as a wave through later blocks with smaller rates.In practice, the later blocks use α2 = ··· = αL = α′, giving α = [α1 +(L−1)α′]/L.

Analysis of the performance of the seeded Belief Propagation procedure

Statistical-physics analysis explains why EM-BP can be trapped for unstructured matrices and how seeded matrices enable exact reconstruction near the optimal measurement rate. Numerical and analytical results show reconstruction nucleating in one block and propagating across the signal.

  • Seeded reconstruction: s-BP combines seeded measurement matrices with EM-BP, whose message dynamics follow gradient ascent toward a maximum of the free entropy.The analytical treatment combines replica and cavity methods, with the cavity method describing the message-passing dynamics.
  • Unstructured matrices: Replica analysis shows that unstructured matrices have a correct free-entropy maximum at D = 0 for α > ρ0, but EM-BP can converge to a secondary maximum at D > 0.The secondary maximum appears for αEM−BP > α > ρ0, preventing exact reconstruction by EM-BP.
  • Unstructured matrices: The EM-BP convergence time diverges as α approaches αEM−BP, and below this threshold it converges to a fixed point with positive mean-squared error.The threshold is analytically identified as the smallest α for which Φ(D) is decreasing.
  • Seeded reconstruction: The seeded matrix uses blockwise Gaussian couplings, with strong backward coupling J1 and weak forward coupling J2 typically giving good performance.The construction divides the signal into L blocks and assigns measurement-dependent variances through Jq,p/N.
  • Seeded reconstruction: For seeded matrices, the cavity-derived dynamics can optimize α1, L, and Jp,q so that exact reconstruction approaches α → ρ0 as L and N/L become large.Finite-size effects slightly degrade this asymptotic threshold saturation.
  • Numerical illustration: In a tested signal with ρ0 = 0.4 and α = 0.5, s-BP nucleates the correct state in the first block and propagates it through the system, while also performing well on real images.The tested regime is described as one where other known algorithms fail, and the non-zero signal components are Gaussian.

Perspectives

The paper outlines extensions aimed at reducing implementation cost, handling noise and prior information, and broadening the measurement model. It also notes independent rigorous confirmation of the asymptotic threshold result.

  • Implementation: Writing EM-BP with N messages instead of M × N parameters connects the approach to AMP and can reduce implementation complexity.This reformulation is especially useful for structured measurement matrices.
  • Implementation: Structured measurement matrices may allow measurements to be obtained in N log N operations rather than M × N operations.The passage presents this as a potential computational advantage when the matrix has suitable structure.
  • Extensions: The approach is reported to be robust to a small amount of measurement noise and extensible to nonlinear measurements and different noise types.These extensions are described as natural directions within the same formalism.
  • Prior information: Using prior information through a better choice of φ can considerably improve performance, especially when non-zero components follow a discrete distribution.The paper identifies continuous non-zero components as the addressed worst case for signals of density ρ0.
  • Related work: A separate work is noted to provide a rigorous proof that α = ρ0 is asymptotically reachable by s-BP when ρ0 = ρ and φ0 = φ.The result is stated for this special matched-prior case.

Appendix A: Proof of the optimality of the probabilistic approach

The appendix proves that, under seeding-matrix conditions and α > ρ0, the probabilistic measure concentrates on the signal in the large-N limit. The proof compares partition-function scaling at zero and positive mean-squared distance.

  • With probability tending to one as N grows, the measure is dominated by the signal when α > ρ0 and α′ > ρ0, assuming finite φ0(0).
  • The proof studies a constrained partition function and its free-entropy density over configurations at fixed distance D from the signal.The construction averages over the measurement matrix and signal, with an indicator enforcing the measurement constraints.
  • For D = 0, the normalized free entropy remains finite as ϵ approaches zero, whereas for D > 0 it vanishes under the stated conditions.This contrast establishes that the measure is dominated by configurations in the neighborhood of the signal.
  • A first-moment bound controls the positive-distance case by comparing the annealed partition function with the quenched quantity.The appendix defines the annealed partition function and uses it to derive the required asymptotic bound.
  • For zero distance, a lower bound is obtained by separating zero and nonzero signal coordinates and substituting the corresponding factors of the prior.The resulting bound contains the explicit term exp[N(ρ0 − α) log ϵ].
  • The positivity of the restricted Wishart-like matrix for α > ρ0 supports the proof that the nonzero signal coordinates yield a strictly positive integral.For seeding matrices, positivity is equivalent to linear independence of the measurement rows restricted to nonzero signal columns.

Appendix B: Derivation of Expectation maximization Belief Propagation

The appendix derives EM-BP by approximating posterior averages with belief propagation and learning the probabilistic prior parameters through Bethe-free-entropy updates. It then reduces the message equations to practical iterative forms and summarizes the implementation.

  • The large-N reduction retains Onsager reaction terms while making messages nearly independent of individual measurements.This yields approximate-message-passing equations asymptotically equivalent to the fully connected belief-propagation formulation.
  • Belief propagation approximates posterior averages that would otherwise require exponential computation, using messages associated with erased measurement constraints.The noiseless model is recovered as the measurement-noise variance tends to zero.
  • For Gaussian priors, the message equations close in terms of means and variances of the component-message distributions.For general φ(xi), these quantities can instead be computed by numerical integration.
  • EM-BP learns the prior density ρ, Gaussian mean x, and variance σ2 while iterating the belief-propagation messages.The parameter updates are based on gradients or stationarity conditions of the Bethe free entropy.
  • The algorithm alternates component and measurement updates, parameter learning, convergence checking, and final signal estimation.The reported implementation initializes messages and parameters, iterates until a criterion or iteration limit, and returns the component estimates.
  • Each iteration costs O(NM) for a general matrix, while recursively computable matrices can reduce the message-passing loop to O(M + N).The observed number of iterations is essentially independent of N, although the constant depends on the parameters and signal.

Appendix D: Replica analysis and density evolution: full measurement matrix

The appendix connects replica calculations with belief-propagation density evolution for a full Gaussian measurement matrix. These equations describe the large-N dynamics and provide the basis for analyzing EM-BP reconstruction.

  • Replica averaging describes both the large-N partition function and the density evolution of belief-propagation messages.
  • The full-matrix free entropy is obtained as a saddle-point function of order parameters, with signal density, signal-component distribution, and noise variance entering the formulation.The appendix states that the noiseless case uses Δ = 0.
  • The order parameters mBP, qBP, and QBP are defined from belief-propagation messages and evolve according to the replica-derived equations.
  • The density-evolution equations for message passing and approximate message passing are identical and close in terms of the mean-squared error E(t).
  • For EM-BP, parameter-learning updates supplement the density-evolution mapping, and the resulting equations are used to analyze reconstruction phase diagrams.The seeding-matrix equations extend the full-matrix case, recovered by taking L = 1.

Appendix E: Replica analysis and density evolution for seeding-measurement matrices

The appendix analyzes seeded belief propagation through replica-based dynamical equations, using them to study convergence and optimize measurement-matrix parameters. The analysis predicts reconstruction at the theoretical limit asymptotically, while finite blocks prevent exact saturation in practice.

  • Density evolution: Replica-based dynamical equations describe convergence toward the successful fixed point, where all layer MSEs and variances vanish within a chosen accuracy.The convergence time is estimated by the iterations needed to reach this fixed point.
  • Parameter optimization: The dynamical system is used to optimize α1, J1, and J2 and estimate good choices of L for Gauss-Bernoulli signals.Its numerical iteration is fast and provides the theoretical infinite-system performance.
  • Asymptotic performance: Perfect reconstruction reaches α = ρ0 asymptotically as L →∞, with corrections scaling as 1/L, for optimally chosen J1 and J2.The passage also reports rigorous work extending and proving this asymptotic claim.
  • Finite-size effects: Finite-size implementations require blocks of a few hundred variables to match theory, so practical choices of L remain limited to several dozens.Exact reconstruction is therefore possible very close to, but not quite at, α = ρ0 in practice.
  • Robustness and tuning: Performance is empirically robust to choices of J1, J2, and α1, and many parameter choices produce effective seeding matrices.A large J1/J2 ratio appears important for short convergence, although better choices may reduce convergence time, finite-size effects, and noise sensitivity.

Appendix F: Phase diagram in the variables used by [5]

Appendix F replots the phase diagram using the convention of [5]. It expresses the diagram in terms of under-sampling and over-sampling ratios.

  • Variable convention: The phase diagram uses ρDT = K/M = ρ/α for the under-sampling ratio and δDT = M/N = α for the over-sampling ratio.These variables replace the original convention while preserving the same data as Fig. 2.

Appendix G: Details on the phantom and Lena examples

The appendix documents reconstruction examples using Haar-wavelet representations of the Shepp–Logan phantom and Lena images, alongside implementation details and comparisons with ℓ1 and BP. It reports moderate s-BP runtimes and strong performance under several noise levels.

  • Image construction: The Shepp–Logan image is a 128^2 phantom using a sparse one-step Haar transform, while Lena is a 128^2 crop retaining 24% of two-step Haar coefficients.All other coefficients in the modified Lena representation are set to zero.
  • Experimental setup: The experiments compare s-BP, BP, and ℓ1 reconstruction, with s-BP coupling parameters listed for the two image examples.The ℓ1 and EM-BP experiments use the same full Gaussian random measurement matrices.
  • Runtime: At α = 0.5, Lena reconstruction takes about 500 iterations and 30 seconds; at α = 0.3, it takes around 2000 iterations and about 2 minutes on a standard laptop.The Matlab implementation is reported as faster than ℓ1-magic on the same machine.
  • Noise comparison: For measurement-noise standard deviations ∆ = 10^-3, 10^-4, and 10^-5, s-BP with L = 9, J1 = 30, J2 = 8, and α1 = 0.8 decodes very well.Under the stated comparison, ℓ1 cannot reconstruct for α < 0.75.

Appendix H: Performance of the algorithm in the presence of measurement noise

Appendix H extends the algorithm to noisy measurements and studies its mean-squared error using density-evolution dynamics. Noise handling is straightforward, and the reported results remain robust to small measurement noise.

  • Noise extension: The noisy-measurement modification of s-BP is straightforward, while the reported reconstruction results are robust to a small amount of noise.A systematic replica study of noisy measurements is identified as beyond the scope of the work.
  • Noise model: The analysis assumes homogeneous Gaussian noise with variance ∆µ = ∆ for every measurement and reports results in units of standard deviation ∆.The noise variance can also be learned through an expectation-maximization approach.
  • Analysis: Density-evolution dynamics are complemented by an iteration on ∆ to study MSE versus α for Gauss-Bernoulli priors and signals.The study compares both full and seeding measurement matrices with ℓ1 minimization.
Loading 1109.4424v4…