CoSaMP Iterative signal recovery from incomplete and inaccurate - PowerPoint PPT Presentation

CoSaMP ❦ Iterative signal recovery from incomplete and inaccurate samples Joel A. Tropp Applied & Computational Mathematics California Institute of Technology jtropp@acm.caltech.edu Joint with D. Needell (UC-Davis). Research supported in part by NSF and DARPA 1

The Sparsity Heuristic A sparse signal has fewer degrees of freedom than its nominal dimension Sparse signal Nearly sparse signal CoSaMP (DIMACS, 26 March 2009) 2

Example: Wavelet Sparsity Courtesy of J. Romberg CoSaMP (DIMACS, 26 March 2009) 3

Example: Time–Frequency Sparsity 0.01 Frequency (MHz) 0.02 0.04 0.05 0.06 0.07 40.08 80.16 120.23 160.31 200.39 Time ( µ s) Data provided by L3 Communications CoSaMP (DIMACS, 26 March 2009) 4

Quantifying Sparsity ❧ Let { ψ k : k = 1 , 2 , . . . , N } be an orthobasis for R N ❧ The coefficients of x with respect to the basis are f k = � x , ψ k � for k = 1 , 2 , . . . , N ❧ The signal is s -sparse when # { k : f k � = 0 } ≤ s ❧ Generalization: the signal is p -compressible with magnitude R if | f | ( k ) ≤ R · k − 1 /p for k = 1 , 2 , . . . , N ❧ p -compressible is slightly weaker than “in ℓ p ” for each p > 0 CoSaMP (DIMACS, 26 March 2009) 5

Approximating Compressible Signals ❧ Consider a signal p -compressible w.r.t. the standard basis | x | ( k ) ≤ R · k − 1 /p for k = 1 , 2 , 3 , . . . ❧ Approximating x by its s largest terms gives error k>s k − 2 /p � 1 / 2 �� � x − x s � 2 ≤ R · �� ∞ � 1 / 2 u − 2 /p d u ≈ R · s 1 / 2 − 1 /p ≈ R · s ❧ Compressible signals are well approximated by sparse signals ❧ Fundamental idea behind transform coding CoSaMP (DIMACS, 26 March 2009) 6

Counting Bits ❧ Consider the class of 0–1 signals in R N with exactly s ones � N � ❧ Clearly need at least log 2 bits to distinguish signals s ❧ By Stirling’s approximation, about s log( N/s ) bits ❧ When s ≪ N , signals contain much less information than the ambient dimension suggests ❧ A simple adaptive coding scheme can achieve this rate CoSaMP (DIMACS, 26 March 2009) 7

What is a Sample? ❧ A sample is the value of a linear functional applied to the signal ❧ Examples: ❧ CCD: Point intensity of an image ❧ ADC: Voltage of an electrical signal at a point in time ❧ MRI: Frequency in the 2D Fourier transform of an image ❧ CAT: Line integral of density in one direction ❧ Some of these technologies acquire samples in batches ❧ We wish to acquire signals with as few samples as possible CoSaMP (DIMACS, 26 March 2009) 8

Compressive Sampling and Signal Recovery ❧ Design linear sampling operator Φ : C N → C m ❧ Suppose x is an unknown (compressible) signal in C N ❧ Collect noisy samples u = Φ x + e ❧ Problem: Given samples u , approximate x CoSaMP (DIMACS, 26 March 2009) 9

Restricted Isometries ❧ Abstract property of sampling operator supports efficient sampling ❧ Φ has the restricted isometry property of order 2 s when (1 − c) � x � 2 2 ≤ � Φ x � 2 2 ≤ (1 + c) � x � 2 whenever � x � 0 ≤ 2 s 2 ❧ Φ preserves geometry of s -sparse signals (take x = y − z ) ❧ W.h.p., a Gaussian sampling operator has RIP( 2 s ) when m ≥ C s log( N/s ) ❧ Gaussian matrices are practically useless References: [Cand` es–Tao 2006, Rudelson–Vershynin 2006] CoSaMP (DIMACS, 26 March 2009) 10

Practical Sampling Operators ❧ Partial Fourier matrices [CRT 2006] ❧ Each row of Φ is chosen at random from rows of unitary DFT F N ❧ Random demodulator [Rice DSP 2006]     1 . . . 1 ± 1 Φ = 1 . . . 1 ± 1 F N     ... ... m × N N × N ❧ W.h.p., both have RIP( 2 s ) when m ≥ C s log α N ❧ Certain technologies can acquire these samples efficiently ❧ Fast matrix–vector multiplies! CoSaMP (DIMACS, 26 March 2009) 11

Desiderata for Recovery Algorithm ❧ Works for general sampling schemes ❧ Succeeds with minimal number of samples ❧ Tolerates noise in samples ❧ Produces approximations with optimal error bound ❧ Yields rigorous guarantees on resource requirements ❧ Exploits structured sampling matrices CoSaMP (DIMACS, 26 March 2009) 12

CoSaMP ( Φ , u , s ) Input: Sampling operator Φ , noisy sample vector u , sparsity level s Output: An s -sparse approximation a of the target signal k = 0 { Initialization } a k = 0 while halting criterion false v ← u − Φ a k { Update samples } y ← Φ ∗ v { Form signal proxy } Ω ← supp( y 2 s ) { Identification } T ← Ω ∪ supp( a k ) { Merge supports } b | T ← Φ † { Signal estimation by least squares } T u b | T c ← 0 a k +1 ← b s { Prune to obtain next approximation } k ← k + 1 end while a ← a k { Return final approximation } CoSaMP (DIMACS, 26 March 2009) 13

Cost per Iteration ❧ Update samples and form signal proxy: v ← u − Φ a k y ← Φ ∗ v and ❧ One matrix–vector multiplication each ❧ Signal approximation by least squares: b T ← Φ † T u ❧ Use conjugate gradient to apply pseudoinverse ❧ Each iteration requires two matrix–vector mulitplies ❧ Assuming RIP( 2 s ), constant number of iterations for fixed accuracy ❧ Constant number of matrix–vector multiplies per CoSaMP iteration! CoSaMP (DIMACS, 26 March 2009) 14

Performance Guarantee Theorem 1. [CoSaMP] Suppose that ❧ the sampling matrix Φ has RIP(2 s ) , ❧ the sample vector u = Φ x + e , ❧ η is a precision parameter, ❧ L bounds cost of a matrix–vector multiply with Φ or Φ ∗ . Then CoSaMP produces a 2 s -sparse approximation a such that � � η, 1 √ s � x − x s � 1 + � e � 2 � x − a � 2 ≤ C max with execution time O( L · log( � x � 2 /η )) . ❧ Need m ≥ C s log α N samples for restricted isometry hypothesis CoSaMP (DIMACS, 26 March 2009) 15

Error Bound for Compressible Signals Corollary 2. [Compressible signals] Suppose ❧ the sampling matrix Φ has RIP(2 s ) , ❧ the signal x is p -compressible with magnitude R , ❧ the sample vector u = Φ x + e , ❧ L bounds cost of a matrix–vector multiply with Φ or Φ ∗ . Then CoSaMP produces a 2 s -sparse approximation a such that Rp − 1 · s 1 / 2 − 1 /p + � e � 2 � � � x − a � 2 ≤ C with execution time O( L · p − 1 log s ) . CoSaMP (DIMACS, 26 March 2009) 16

To learn more... E-mail: ❧ jtropp@acm.caltech.edu ❧ dneedell@math.ucdavis.edu Web: http://www.acm.caltech.edu/~jtropp Relevant Papers: ❧ NTV, “CoSaMP: Iterative signal recovery from incomplete and inaccurate samples,” ACHA 2009 ❧ T and Rice DSP, “Beyond Nyquist: Efficient sampling of sparse, bandlimited signals,” submitted ❧ N and Vershynin, “Stable signal recovery from incomplete and inaccurate samples,” submitted ❧ T and Gilbert, “Signal recovery from random measurements via Orthogonal Matching Pursuit,” Trans. IT , Dec. 2007. CoSaMP (DIMACS, 26 March 2009) 17