Sparse Reconstruction with Compressed Sensing: Difference between revisions

From Cornell University Computational Optimization Open Textbook - Optimization Wiki
Jump to navigation Jump to search
No edit summary
No edit summary
Line 1: Line 1:
Author: Ngoc Ly (SysEn 5800 Fall 2021)
Author: Ngoc Ly (SysEn 5800 Fall 2021)


 
==  Introduction ==
Introduction
<math>x \in \mathbb{R}^N</math>often not really sparse but approximately sparse
 
<math>x \in \mathbb{R}^N</math>often not really sparse but aproximently sparse


<math>\Phi \in \mathbb{R}^{M \times N}</math>for i<math>M \ll N</math>s a Random Gaussian or Bernoulli matrix
<math>\Phi \in \mathbb{R}^{M \times N}</math>for i<math>M \ll N</math>s a Random Gaussian or Bernoulli matrix
Line 16: Line 14:
How can we reconstruct x from  <math>y = \Phi x + e</math>
How can we reconstruct x from  <math>y = \Phi x + e</math>


The goal is to reconstruct x 0 ∈ R N <math>x \in \mathbb{R}^N</math>given <math>y</math> and <math>\Phi</math>
The goal is to reconstruct <math>x \in \mathbb{R}^N</math>given <math>y</math> and <math>\Phi</math>
 
Sensing matrix <math>\Phi</math>must satisfy RIP i.e. Random Gaussian or Bernoulli ma-
 
traxies satsifies (cite)
 
let <math>[ N ] = \{ 1, \dots , N \} </math>be an index set
 
<math>[N]</math> enumerates the colums of <math>\Phi</math> and <math>x</math>


<math>\Phi</math> is an underdetermined systems with infinite solutions. Why <math>l_2</math> norm does
Sensing matrix <math>\Phi</math>must satisfy RIP i.e. Random Gaussian or Bernoulli matrixies satisfies (cite)


not work
let <math>[ N ] = \{ 1, \dots , N \} </math>be an index set <math>[N]</math> enumerates the columns of <math>\Phi</math> and <math>x</math> <math>\Phi</math> is an under determined systems with infinite solutions. Why <math>l_2</math> norm does not work


What is compression is sunomunus with to the sparcity.
What is compression is synonymous with to the sparsity.


The problem formulation is to recover sparse data x 0 ∈ R N <math> \mathbf{x} \in \mathbb{R}^N </math>
The problem formulation is to recover sparse data <math> \mathbf{x} \in \mathbb{R}^N </math>


The support of <math>\mathbf{x}</math> is <math>supp(\mathbf{x}) = \{i \in [N] : \mathbf{x}_i \neq 0 \}</math>


we say <math>\mathbf{x}</math> is <math>k</math> sparse when <math>|supp(x)| \leq k</math>
The support of <math>\mathbf{x}</math> is <math>supp(\mathbf{x}) = \{i \in [N] : \mathbf{x}_i \neq 0 \}</math> we say <math>\mathbf{x}</math> is <math>k</math> sparse when <math>|supp(x)| \leq k</math>


We are interested in the smallest <math>supp(x)</math> , i.e. <math>min(supp(x))</math>
We are interested in the smallest <math>supp(x)</math> , i.e. <math>min(supp(x))</math>


Restricted isometry property (RIP) Introduced by Candes, Tao
Random Gausian and Bernoulli satsifies RIP


Let Φ ∈ R M × N <math>\Phi \in \mathbb{R}^{M \times N}</math> satsify RIP, Let <math>[N]</math> be an index set


Before we get into RIP lets talk about RIC
Before we get into RIP lets talk about RIC


Resiricted Isometry Constant (RIC) is the smallest δ | s s.t. s ⊆ [ N ] <math>\delta_{|s} \ s.t. \  s \subseteq [N]</math>that satsi-
Restricted Isometry Constant (RIC) is the smallest <math>\delta_{|s} \ s.t. \  s \subseteq [N]</math>that satisfies the RIP condition introduced by Candes, Tao
 
fies the RIP condition
 
For <math>s</math> is a restrciton on <math>\mathbf{x}</math> denoted by <math>x_{|s}</math> <math>x \in \mathbb{R}^N</math> to <math>s</math> k-sparse <math>\mathbf{x}</math> s.t. RIP is


satisified the s = | Γ | i.e. s ⊆ [ N ] and Φ | s ⊆ Φ where the columns of Φ | s is
Random Gaussian and Bernoulli satisfies RIP


indexed by i ∈ S
Let <math>\Phi \in \mathbb{R}^{M \times N}</math> satisfy RIP, Let <math>[N]</math> be an index set For <math>s</math> is a restriction on <math>\mathbf{x}</math> denoted by <math>x_{|s}</math> <math>x \in \mathbb{R}^N</math> to <math>s</math> k-sparse <math>\mathbf{x}</math> s.t. RIP is satisfied the s = | Γ | i.e. s ⊆ [ N ] and Φ | s ⊆ Φ where the columns of Φ | s is indexed by i ∈ S


Line 122: Line 103:
y = Φx + e = Φx s + Φx r + e = Φx s + ẽ
y = Φx + e = Φx s + Φx r + e = Φx s + ẽ


If Φ statisfies RIP for sparcity s, then the norm of error ẽ is bounded by
If Φ satisfies RIP for sparsity s, then the norm of error ẽ is bounded by


k ẽ k 2 ≤
k ẽ k 2 ≤
Line 144: Line 125:
∀ x
∀ x


Theory
== Theory ==
 
Gel’fend n-width
Gel’fend n-width


Line 164: Line 144:
We want a small µ A because it will be close ot the normal matrix, which
We want a small µ A because it will be close ot the normal matrix, which


will satisfie RIP
will satisfies RIP
 
Algorithm IHT


* Initalize <math> \Phi, \mathbf{y}, \mathbf{e} \ \mbox{with} \ \mathbf{y} = \mathbf{\Phi} \mathbf{x} | \mathbf{e} and \mathfrak{M}</math>
=== Algorithm IHT ===
* Initialize <math> \Phi, \mathbf{y}, \mathbf{e} \ \mbox{with} \ \mathbf{y} = \mathbf{\Phi} \mathbf{x} | \mathbf{e} and \mathfrak{M}</math>


* output <math>IHT(\mathbf{y}, \mathbf{\Phi}, \mathfrak{M}) </math>
* output <math>IHT(\mathbf{y}, \mathbf{\Phi}, \mathfrak{M}) </math>


* While halting criterion false do
* While halting criterion false do
* <math>x^{(n+1)} \leftarrow \P_{\mathfrak{M} \left[ x^{(n) + \Phi^* (\mathbf{y} - \mathbf{\Phi x}^{(n)})}\right]}</math>
*<math>x^{(n+1)} \leftarrow \P_{\mathfrak{M} \left[ x^{(n) + \Phi^* (\mathbf{y} - \mathbf{\Phi x}^{(n)})}\right]}</math>
* <math>n \leftarrow n + 1 </math>
*<math>n \leftarrow n + 1 </math>
end while
end while
return: <math>IHT(\mathbf{y}, \mathbf{\Phi}, \mathfrak{M}) \leftarrow \mathbf{x}^{(n)}</math>
return: <math>IHT(\mathbf{y}, \mathbf{\Phi}, \mathfrak{M}) \leftarrow \mathbf{x}^{(n)}</math>


 
== Numerical Example ==
 
Numerical Example
 
Basis Persuit


Iterative Hard Thresholding IHT
Iterative Hard Thresholding IHT
Solver
2

Revision as of 05:06, 28 November 2021

Author: Ngoc Ly (SysEn 5800 Fall 2021)

Introduction

often not really sparse but approximately sparse

for is a Random Gaussian or Bernoulli matrix

are the observed y samples

noise vector k e k 2 ≤ η

s.t.

How can we reconstruct x from

The goal is to reconstruct given and

Sensing matrix must satisfy RIP i.e. Random Gaussian or Bernoulli matrixies satisfies (cite)

let be an index set enumerates the columns of and is an under determined systems with infinite solutions. Why norm does not work

What is compression is synonymous with to the sparsity.

The problem formulation is to recover sparse data


The support of is we say is sparse when

We are interested in the smallest , i.e.


Before we get into RIP lets talk about RIC

Restricted Isometry Constant (RIC) is the smallest that satisfies the RIP condition introduced by Candes, Tao

Random Gaussian and Bernoulli satisfies RIP

Let satisfy RIP, Let be an index set For is a restriction on denoted by to k-sparse s.t. RIP is satisfied the s = | Γ | i.e. s ⊆ [ N ] and Φ | s ⊆ Φ where the columns of Φ | s is indexed by i ∈ S

x i , i f i ∈ S

( x | S ) i =

0 otherwise

RIP defined as

( 1 − δ s )k x k 22 ≤ k Φx k 22 ≤ ( 1 + δ s )k x k 22

3 Lemmas Page 267 Blumensath Davies IHT for CS

Lemma 1(Blumensath, Davis 2009 Iterative hard thresholding for compressed

sensing), For all index sets Γ and all Φ for which RIP holds with s = | Γ | that is

s = supp ( x )

1k Φ Γ T k 2 ≤

q

1 + δ | Γ | k y k 2

( 1 − δ | Γ | )k x Γ k 22 ≤ k Φ Γ T Φ Γ x Γ k 22 ≤ ( 1 + δ | Γ | )k x Γ k 22

and

k( I − Φ Γ T Φ Γ )k 2 ≤ δ | Γ | k x Γ k 2

SupposeΓ ∩ Λ = ∅

k Φ Γ T Φ Λ ) x Λ k 2 ≤ δ s k x Λ k 2

Lemma 2 (Needell Tropp, Prop 3.5 in CoSaMP: Iterative signal recovery

from incomplete and inaccurate √ samples)

If Φ satisfies RIP k Φx s k 2 ≤ 1 + δ s k x s k 2 , ∀ x s : k x s k 0 ≤ s, Then ∀ x

k Φx k 2 ≤

p

1 + δ s k x k 2 +

p

1 + δ s

k x k 1

sqrts

Lemma 3 (Needell Tropp, Prop 3.5 in CoSaMP: Iterative signal recovery

from incomplete and inaccurate samples)

Let x s be the best s-term approximation to x. Let x r = x − x s Let

y = Φx + e = Φx s + Φx r + e = Φx s + ẽ

If Φ satisfies RIP for sparsity s, then the norm of error ẽ is bounded by

k ẽ k 2 ≤

p

1 + δ s k x − x s k 2 +

p

1 + δ s

k x − x s k 1

+ k e k 2

s

∀ x

Theory

Gel’fend n-width

Errors E ( S, Φ, D )

Definition Mutual Coherence

LetA ∈ R M × N , themutualcoherenceµ A isde f inedby :

µ A =

|h a i , a j i|

k a i kk a j k

i 6 = j

We want a small µ A because it will be close ot the normal matrix, which

will satisfies RIP

Algorithm IHT

  • Initialize
  • output
  • While halting criterion false do

end while return:

Numerical Example

Iterative Hard Thresholding IHT