Sparse approximation by modified Prony method Gerlind Plonka and - - PowerPoint PPT Presentation

Sparse approximation by modified Prony method

Gerlind Plonka and Vlada Pototskaia

Institut f¨ ur Numerische und Angewandte Mathematik Georg-August-Universit¨ at G¨

• ttingen

Alba di Canazei, September 19, 2016

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 1 / 29

Outline

1

Sparse approximation problem for exponential sums

2

Prony’s method

3

The AAK theorem for samples of exponential sums

4

Method for sparse approximation of exponential sums

5

Numerical example

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 2 / 29

Outline

1

Sparse approximation problem for exponential sums

2

Prony’s method

3

The AAK theorem for samples of exponential sums

4

Method for sparse approximation of exponential sums

5

Numerical example

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 3 / 29

Sparse approximation of exponential sums

Consider a function of the form f(x) =

N

• j=1

ajzx

j

with |zj| < 1, where aj, zj ∈ C.

Goal:

Find a function ˜ f(x) =

n

• j=1

˜ aj˜ zx

j

with |˜ zj| < 1, such that n < N and f − ˜ f ≤ ε

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 4 / 29

Discrete sparse approximation problem

Consider a sequence of samples f := (fk)∞

k=0 given by

fk := f(k) =

N

• j=1

ajzk

j

with |zj| < 1, where aj, zj ∈ C.

Goal:

Find a sequence ˜ f := ( ˜ fk)∞

k=0 of the form

˜ fk =

n

• j=1

˜ aj˜ zk

j

with |˜ zj| < 1, such that n < N and f − ˜ fℓ2 ≤ ε

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 5 / 29

Possible applications

We consider here a structured low-rank approximation problem for model reduction. Problem: Low-rank approximation using the SVD destroys the Hankel

• structure. [Markovsky, 2008]

Applications Approximation of special functions by exponential sums, e.g. Bessel functions, or x−1/2 to avoid quadrature methods for Schr¨

• dinger
• equations. [Beylkin, Monzon, 2005], [Hackbusch, 2005]

Signal compression by sparse representation of the (discrete) Fourier transform.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 6 / 29

Our approach

(1) Given a sufficiently large number of samples fk, reconstruct zj and aj such that fk =

N

• j=1

ajzk

j

with |zj| < 1 using a Prony-like method.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 7 / 29

Our approach

(1) Given a sufficiently large number of samples fk, reconstruct zj and aj such that fk =

N

• j=1

ajzk

j

with |zj| < 1 using a Prony-like method. (2) Given the representation (1), find ˜ zj and ˜ aj such that for ˜ fk =

n

• j=1

˜ aj˜ zk

j

with |˜ zj| < 1 and n < N we have f − ˜ fℓ2 ≤ ε using the AAK Theorem [Adamjan, Arov, Krein], (1971).

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 7 / 29

Outline

1

Sparse approximation problem for exponential sums

2

Prony’s method

3

The AAK theorem for samples of exponential sums

4

Method for sparse approximation of exponential sums

5

Numerical example

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 8 / 29

Classical Prony’s Method (1795)

Assume f(k) =

N

• j=1

ajzk

j

with zj := eTj Given: N and fk := f(k) for k = 0, . . . , 2N − 1 Wanted: zj ∈ C , aj ∈ C

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 9 / 29

Classical Prony’s Method (1795)

Assume f(k) =

N

• j=1

ajzk

j

with zj := eTj Given: N and fk := f(k) for k = 0, . . . , 2N − 1 Wanted: zj ∈ C , aj ∈ C Consider the Prony polynomial P(x) :=

N

• j=1

(x − zj) =

N

• k=0

pkxk, with pN = 1.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 9 / 29

Classical Prony’s Method (1795)

Assume f(k) =

N

• j=1

ajzk

j

with zj := eTj Given: N and fk := f(k) for k = 0, . . . , 2N − 1 Wanted: zj ∈ C , aj ∈ C Consider the Prony polynomial P(x) :=

N

• j=1

(x − zj) =

N

• k=0

pkxk, with pN = 1. We have for l = 0, . . . , N − 1

N

• k=0

pkfl+k =

N

• k=0

pk

N

• j=1

ajz(l+k)

j

=

N

• j=1

ajzl

j N

• k=0

pkzk

j = 0

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 9 / 29

Classical Prony method (1795)

Since pN = 1, the last equation can be written as

N

• k=0

pkfl+k =

N−1

• k=0

pkfl+k + fl+N = 0 ⇔

N−1

• k=0

pkfl+k = −fl+N and defines a homogeneous difference equation of order N. Matrix-vector-representation:     f0 f1 · · · fN−1 f1 f2 · · · fN . . . . . . ... . . . fN−1 fN · · · f2N−2         p0 p1 . . . pN−1     = −     fN fN+1 . . . f2N−1     .

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 10 / 29

Literature

[Prony] (1795) Reconstruction of a difference equation [Schmidt] (1979) MUSIC (Multiple Signal Classification) [Roy,Kailath] (1989) ESPRIT (Estimation of signal parameters via rotational invariance techniques) [Hua,Sakar] (1990) Matrix-Pencil method [Stoica,Moses] (2000) Annihilating filters [Potts,Tasche] (2010,2011) Approximate Prony method [Kunis et al.], [Sauer] (2015) Multivariate Prony’s method Golub, Milanfar, Varah (’99); Vetterli, Marziliano, Blu (’02); Maravi´ c, Vetterli (’04); Elad, Milanfar, Golub (’04); Beylkin, Monz´

• n (’05,’10); Batenkov, Yomdin (’12,’13);

Filbir et al. (’12); Potts, Tasche (’11,’12,’13); Plonka, Wischerhoff (’13, ’16); Peter, Plonka (’13); Cuyt, Lee, Tsai (’16); Diederichs, Iske (’16) ....

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 11 / 29

Outline

1

Sparse approximation problem for exponential sums

2

Prony’s method

3

The AAK theorem for samples of exponential sums

4

Method for sparse approximation of exponential sums

5

Numerical example

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 12 / 29

AAK Theorem for Samples of Exponential Sums

Consider the sequence f := (fk)∞

k=0 given by samples

fk = f(k) =

N

• j=1

ajzk

j

with 0 < |zj| < 1 and let D := {z ∈ C : 0 < |z| < 1}. We define the infinite Hankel matrix Γf :=     f0 f1 f2 · · · f1 f2 f3 · · · f2 f3 f4 · · · . . . . . . . . . ...     = (fk+j)∞

k,j=0

with respect to f.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 13 / 29

AAK theorem for samples of exponential sums

Then Γf has the following properties: Γf defines a compact operator on ℓ2 = ℓ2(N). Γf has finite rank N. The singular values of Γf are of the form σ0 ≥ σ1 ≥ . . . ≥ σN−1 > σN = . . . = σ∞ = 0.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 14 / 29

AAK theorem for samples of exponential sums

Then Γf has the following properties: Γf defines a compact operator on ℓ2 = ℓ2(N). Γf has finite rank N. The singular values of Γf are of the form σ0 ≥ σ1 ≥ . . . ≥ σN−1 > σN = . . . = σ∞ = 0.

[Young] (1988) An Introduction to Hilbert Space [Chui, Chen] (1992) Discrete H∞ optimization [Peller] (2000) Hankel Operators and Their Applications [Beylkin,Monz´

• n]

(2005) On approximation of functions by exponential sums [Andersson et al.] (2011) Sparse approximation of functions using sums

• f exponentials and AAK theory

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 14 / 29

The AAK theorem (Adamjan, Arov, Krein, 1971)

Let f := (f(k))∞

k=0 be given as before.

Let (σn, un) be a fixed singular pair of Γf with σn / ∈ {σk}k=n and σn = 0. Then Pun(x) :=

∞

• k=0

un(k)xk has exactly n zeros ˜ z1, . . . , ˜ zn in D, repeated according to multiplicity.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 15 / 29

The AAK theorem (Adamjan, Arov, Krein, 1971)

Let f := (f(k))∞

k=0 be given as before.

Let (σn, un) be a fixed singular pair of Γf with σn / ∈ {σk}k=n and σn = 0. Then Pun(x) :=

∞

• k=0

un(k)xk has exactly n zeros ˜ z1, . . . , ˜ zn in D, repeated according to multiplicity. If the ˜ zk are pairwise different, then there are ˜ a1, . . . , ˜ an ∈ C such that for ˜ f = ( ˜ fj)∞

j=0 =

n

• k=1

˜ ak˜ zkj ∞

j=0

it holds that Γf − Γ ˜

fℓ2→ℓ2 = σn.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 15 / 29

The AAK theorem

Let Γf be of rank N and the singular values be of the form σ0 > σ1 > . . . > σN−1 > σN = . . . = σ∞ = 0. zeros of n σn Pun(x) in D ˜ f Γf − Γ ˜

f

σ0 − σ0 1 σ1 ˜ z1 ˜ fj = ˜ a˜ zj

1

σ1 2 σ2 ˜ z1, ˜ z2 ˜ fj = ˜ a1˜ zj

1 + ˜

a2˜ zj

2

σ2 3 σ3 ˜ z1, ˜ z2, ˜ z3 ˜ fj = ˜ a1˜ zj

1 + ˜

a2˜ zj

2 + ˜

a3˜ zj

3

σ3 . . . . . . . . . . . . . . . N − 1 σN−1 ˜ z1, . . . , ˜ zN−1 ˜ fj = N−1

k=1 ˜

ak˜ zj

k

σN−1

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 16 / 29

The AAK theorem

Let Γf be of rank N and the singular values be of the form σ0 > σ1 > . . . > σN−1 > σN = . . . = σ∞ = 0. zeros of n σn Pun(x) in D ˜ f Γf − Γ ˜

f

σ0 − σ0 1 σ1 ˜ z1 ˜ fj = ˜ a˜ zj

1

σ1 2 σ2 ˜ z1, ˜ z2 ˜ fj = ˜ a1˜ zj

1 + ˜

a2˜ zj

2

σ2 3 σ3 ˜ z1, ˜ z2, ˜ z3 ˜ fj = ˜ a1˜ zj

1 + ˜

a2˜ zj

2 + ˜

a3˜ zj

3

σ3 . . . . . . . . . . . . . . . N − 1 σN−1 ˜ z1, . . . , ˜ zN−1 ˜ fj = N−1

k=1 ˜

ak˜ zj

k

σN−1 Original sequence: fj = N

k=1 akzj k

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 16 / 29

Problems with application of AAK theory

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 17 / 29

Problems with application of AAK theory

Let f be given as before and let (σn, un) be a fixed singular pair of Γf such that σn / ∈ {σk}k=n and σn = σ∞. Then Pun(x) :=

∞

• k=0

un(k)xk has exactly n zeros ˜ z1, . . . , ˜ zn in D, repeated according to multiplicity. If the ˜ zk are pairwise different, then there are ˜ a1, . . . , ˜ an ∈ C such that for ˜ f = ( ˜ fj)∞

j=0 =

n

• k=1

˜ ak˜ zj

k

∞

j=0

it holds that Γf − Γ ˜

f ℓ2→ℓ2 = σn.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 17 / 29

Outline

1

Sparse approximation problem for exponential sums

2

Prony’s method

3

The AAK theorem for samples of exponential sums

4

Method for sparse approximation of exponential sums

5

Numerical example

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 18 / 29

Singular values and con-eigenvalues

For a (complex) Hankel matrix Γf we call σ ∈ C a con-eigenvalue with the corresponding con-eigenvector v ∈ ℓ2(N) if it satisfies Γfv = σv. For symmetric matrices like Γf we have We can always select a nonnegative σ. (σ, v) is a

multiplicity is 1

− − − − − − − − − ⇀ ↽ − − − − − − − − − (σ, v) is a singular pair of Γf con-eigenpair of Γf

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 19 / 29

Structure of con-eigenvectors to non-zero con-eigenvalues

Lemma: Let f be given as before, i.e. fk =

N

• j=1

ajzjk with zj ∈ D, and let σ = 0 be a fixed con-eigenvalue of Γf with the corresponding con-eigenvector u := (uk)∞

k=0.

Then u can be represented by uk =

N

• j=1

bjzjk, k = 0, 1, . . . , where bj, j = 1, . . . , N are some (complex or real) coefficients.

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 20 / 29

Dimension reduction for the con-eigenvalue problem of Γf

Γf ¯ u = σu ⇔

∞

• j=0

fj+k¯ uj = σuk, ∀k = 0, 1, 2, . . . ⇔

∞

• j=0

N

• l=1

alzk+j

l

N

• s=1

bszj

s

• = σ

N

• l=1

blzk

l

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 21 / 29

Dimension reduction for the con-eigenvalue problem of Γf

Γf ¯ u = σu ⇔

∞

• j=0

fj+k¯ uj = σuk, ∀k = 0, 1, 2, . . . ⇔

∞

• j=0

N

• l=1

alzk+j

l

N

• s=1

bszj

s

• = σ

N

• l=1

blzk

l

⇔

N

• l=1

zk

l

 al

N

• s=1

¯ bs

∞

• j=0

(zl¯ zs)j   = σ

N

• l=1

blzk

l .

⇔

N

• l=1

zk

l

• al

N

• s=1

¯ bs 1 − zl¯ zs

• =

N

• l=1

(σbl) zk

l .

⇔ al

N

• s=1

¯ bs 1 − zl¯ zs = σbl ∀ l = 1, . . . , N

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 21 / 29

Dimension reduction for the con-eigenvalue problem of Γf

The last equation can be seen as the following con-eigenvalue problem of the dimension N AZ¯ b = σb, where A :=      a1 a2 ... aN      , Z :=      

1 1−|z1|2 1 1−¯ z2z1

· · ·

1 1−¯ zNz1 1 1−¯ z1z2 1 1−|z2|2

· · ·

1 1−¯ zNz2

. . . . . . ... . . .

1 1−¯ z1zN 1 1−¯ z2zN

· · ·

1 1−|zN|2

      and b := (b1, . . . , bN)T .

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 22 / 29

Computation of the roots of con-eigenpolynomials of Γf

Let Pu(x) be the n-th con-eigenpolynomial of Γf. Then for |x| < 1 we obtain Pu(x) =

∞

• k=0

ukxk =

∞

• k=0

 

N

• j=1

bjzk

j

  xk =

N

• j=1

bj

∞

• k=0

(zjx)k =

N

• j=1

bj 1 − zjx

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 23 / 29

Norm of the Hankel Operator: Γf vs. f

Let e1 := (1, 0, 0, . . .)T . Then fℓ2 =  

∞

• j=0

|fj|2  

1/2

= Γfe1ℓ2 ≤ sup

xℓ2=1

Γfxℓ2 = Γf. Therefore we have f − ˜ fℓ2 ≤ Γf− ˜

f = σn

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 24 / 29

Algorithm for sparse approximation of exponential sums

Input: samples fk, k = 0, . . . , L, for a sufficiently large L, target approximation error ε.

• 1. Find the N nodes zj and the weights aj of the exponential

representation of f using a Prony-like method.

• 2. Compute a con-eigenvalue σn < ε of the matrix AZ and the

corresponding con-eigenvector u = un.

• 3. Compute the n zeros ˜

zj of the con-eigenpolynomial Pu(x) of Γf in D using the rational function representation.

• 4. Compute the new coefficients ˜

aj by solving min

˜ a1,...,˜ an f − ˜

f2

ℓ2 =

min

˜ a1,...,˜ an ∞

• k=0

|fk −

n

• j=1

˜ aj(˜ zj)k|2. Output: sequence ˜ fk = n

j=1 ˜

aj˜ zk

j , such that f − ˜

fℓ2 ≤ σn < ε

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 25 / 29

Outline

1

Sparse approximation problem for exponential sums

2

Prony’s method

3

The AAK theorem for samples of exponential sums

4

Method for sparse approximation of exponential sums

5

Numerical example

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 26 / 29

Numerical example

N=6 fk = 6

j=1 ajzk j

aj = 5, j = 1, . . . , 6

n = 5 n = 4 n = 3 n = 2 n = 1 z1 = 0.3500 0.3509 0.3550 0.3671 0.3985 0.4889 ˜ z1 z2 = 0.4000 0.4103 0.4365 0.4860 0.5684 ˜ z2 z3 = 0.4500 0.4802 0.5282 0.5910 ˜ z3 z4 = 0.5000 0.5456 0.5981 ˜ z4 z5 = 0.5500 0.5998 ˜ z5 z6 = 0.6000 4.5845e-10 1.6340e-07 3.1318e-05 4.3318e-03 4.8259e-01 σn

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 27 / 29

Numerical example

N=6 fk = 6

j=1 ajzk j

aj = 5, j = 1, . . . , 6

n σn f − ˜ f2

maxk |fk− ˜ fk| maxk |fk|

1 4.8259e-01 4.7095e-01 1.1013e-02 2 4.3318e-03 4.2576e-03 7.6860e-05 3 3.1318e-05 2.8624e-05 5.9415e-07 4 1.6340e-07 1.4449e-07 2.9658e-09 5 4.5845e-10 8.0184e-10 1.1560e-11

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 28 / 29

Thank You !

Gerlind Plonka, Vlada Pototskaia Prony method for sparse approximation Canazei 2016 29 / 29