Simple and Practical Algorithm for the Sparse Fourier Transform - - PowerPoint PPT Presentation

Simple and Practical Algorithm for the Sparse Fourier Transform

Haitham Hassanieh Piotr Indyk Dina Katabi Eric Price

MIT

2012-01-19

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 1 / 19

Outline

1

Introduction

2

Algorithm

3

Experiments

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 2 / 19

The Dicrete Fourier Transform

Discrete Fourier transform: given x ∈ Cn, find

• xi =
• xjωij

Fundamental tool

◮ Compression (audio, image, video) ◮ Signal processing ◮ Data analysis ◮ ...

FFT: O(n log n) time.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 4 / 19

Sparse Fourier Transform

Often the Fourier transform is dominated by a small number of “peaks”

◮ Precisely the reason to use for compression.

If most of mass in k locations, can we compute FFT faster?

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 5 / 19

Previous work

Boolean cube: [KM92], [GL89]. What about C? [Mansour-92]: kc logc n. Long list of other work [GGIMS02, AGS03, Iwen10, Aka10] Fastest is [Gilbert-Muthukrishnan-Strauss-05]: k log4 n.

◮ All have poor constants, many logs. ◮ Need n/k > 40,000 or ω(log3 n) to beat FFTW. ◮ Our goal: beat FFTW for smaller n/k in theory and practice. ◮ Result: n/k > 2, 000 or ω(log n) to beat FFTW. Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 6 / 19

Our result

Simple, practical algorithm with good constants. Compute the k-sparse Fourier transform in O( √ kn log3/2 n) time. Get x′ with approximation error

• x′ −

x2

∞ ≤ 1

k x − xk2

2

If x is sparse, recover it exactly. Caveats:

◮ Additional x2/nΘ(1) error. ◮ n must be a power of 2. Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 7 / 19

Structure of this section

If x is k-sparse with known support S, find xS exactly in O(k log2 n) time. In general, estimate x approximately in O( √ nk) time.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 9 / 19

Intuition

Original signal (time) Original signal (freq) Cutoff signal (time) Cutoff signal (freq) Cutoff signal (time) Cutoff signal, subsampled (freq)

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 10 / 19

n-dimensional DFT n-dimensional DFT

• f first B terms.

B-dimensional DFT

• f first B terms.

Framework

Cutoff signal (time) Cutoff subsampled signals (freq)

“Hashes” into B buckets in B log B time. Issues:

◮ “Hashing” needs a random hash function ⋆ Access x′ t = ω−atxσt, so

x′t = xσ−1t+a [GMS-05]

◮ Collisions ⋆ Have B > 4k, repeat O(log n) times and take median. [Count-Sketch,

CCF02]

◮ Leakage ◮ Finding the support. [Porat-Strauss-12], talk at 9:45am. Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 11 / 19

Leakage

Cutoff signal (time) Cutoff, subsampled signals (freq)

Let Fi = 1 i < B

• therwise

be the “boxcar” filter. (Used in [GGIMS02,GMS05]) Observe DFT(F·x, B) = subsample(DFT(F·x, n), B) = subsample( F∗ x, B). DFT F of boxcar filter is sinc, decays as 1/i. Need a better filter F!

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 12 / 19

Filters

20 40 60 80 100 0.0 0.5 1.0 1.5 2.0

Filter (time)

Bin 5 10 15 20 25

Filter (freq)

Observe subsample( F ∗ x, B) in O(B log B) time. Needs for F:

◮ supp(F) ∈ [0, B] ◮ |

F| < δ = 1/nΘ(1) except “near” 0.

◮

F ≈ 1 over [−n/2B, n/2B].

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 13 / 19

Filters

20 40 60 80 100 0.0 0.2 0.4 0.6 0.8 1.0

Filter (time)

Bin 2 4 6 8 10 12 14 16 18

Filter (freq)

Observe subsample( F ∗ x, B) in O(B log B) time. Needs for F:

◮ supp(F) ∈ [0, B] ◮ |

F| < δ = 1/nΘ(1) except “near” 0.

◮

F ≈ 1 over [−n/2B, n/2B].

Gaussians:

◮ Standard deviation σ = B/

• log n

◮ DFT has

σ = (n/B)

• log n

◮ Nontrivial leakage into O(log n) buckets. ◮ But likely trivial contribution to correct bucket. Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 13 / 19

Filters

50 100 150 200 0.0 0.2 0.4 0.6 0.8 1.0

Filter (time)

Bin 10 20 30 40 50 60 70 80

Filter (freq)

Observe subsample( F ∗ x, B) in O(B log B) time. Needs for F:

◮ supp(F) ∈ [0, B log n] ◮ |

F| < δ = 1/nΘ(1) except “near” 0.

◮

F ≈ 1 over [−n/2B, n/2B].

Gaussians:

◮ Standard deviation σ = ✘✘✘✘

✘ B/

• log n B ·
• log n

◮ DFT has

σ = ✘✘✘✘✘ ✘ (n/B)

• log n (n/B)/
• log n

◮ Nontrivial leakage into 0 buckets. ◮ But likely trivial contribution to correct bucket. Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 13 / 19

Filters

50 100 150 200 0.0 0.2 0.4 0.6 0.8 1.0

Filter (time)

Bin 10 20 30 40 50 60 70 80

Filter (freq)

Let G be Gaussian with σ = B

• log n

H be box-car filter of length n/B.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 14 / 19

Filters

50 100 150 200 0.02 0.00 0.02 0.04 0.06 0.08 0.10

Filter (time)

Bin 0.0 0.2 0.4 0.6 0.8 1.0 1.2

Filter (freq)

Let G be Gaussian with σ = B

• log n

H be box-car filter of length n/B. Use F = G ∗ H.

◮ F = G ·

H, so supp(F) ⊂ [0, B log n].

◮ |

F| < 1/nΘ(1) outside −n/B, n/B.

◮ |

F| = 1 ± 1/nΘ(1) within n/2B, n/B.

Hashes correctly to one bucket, leaks to at most 1 bucket. Replace Gaussians with “Dolph-Chebyshev window functions”: factor 2 improvement.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 14 / 19

Algorithm to estimate xS

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Algorithm to estimate xS

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Algorithm to estimate xS

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Algorithm to estimate xS

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Algorithm in general

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Algorithm in general

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to. To find S: choose all that map to the top 2k values.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Algorithm in general

For O(log n) different permutations of x, compute subsample( F ∗ x, B). Estimate each xi as median of values it maps to. To find S: choose all that map to the top 2k values. nk/B candidates to update at each iteration: total (nk B + B log n) log n = √ nk log3/2 n time.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 15 / 19

Empirical Performance: runtime

0.001 0.01 0.1 1 10 214 215 216 217 218 219 220 221 222 223 224 225 226

Run Time (sec) Signal Size (n) Run Time vs Signal Size (k=50)

sFFT 1.0 sFFT 2.0 FFTW FFTW OPT AAFFT 0.9 0.01 0.1 1 10 26 27 28 29 210 211 212

Run Time (sec) Sparsity (K) Run Time vs Signal Sparsity (n=222)

sFFT 1.0 sFFT 2.0 FFTW FFTW OPT AAFFT 0.9

Compare to FFTW, previous best sublinear algorithm (AAFFT). Offer a heuristic that improves time to O(n1/3k2/3).

◮ Filter from [Mansour ’92]. ◮ Can’t rerandomize, might miss elements.

Faster than FFTW for n/k > 2,000. Faster than AAFFT for n/k < 1,000,000.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 17 / 19

Empirical Performance: noise

1e-07 1e-06 1e-05 0.0001 0.001 0.01 0.1 1

• 20

20 40 60 80 100 Average L1 Error per Enrty SNR (dB)

Robustness vs SNR (n=222, k=50)

sFFT 1.0 sFFT 2.0 AAFFT 0.9

Just like in Count-Sketch, algorithm is noise tolerant.

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 18 / 19

Conclusions

Roughly: fastest algorithm for n/k ∈ [2 × 103, 106]. Recent improvements [HIKP12b?]

◮ O(k log n) for exactly sparse

x

◮ O(k log n

k log n) for approximation.

◮ Beats FFTW for n/k > 400 (in the exact case). Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 19 / 19

Hassanieh, Indyk, Katabi, and Price (MIT) Simple and Practical Algorithm for the Sparse Fourier Transform 2012-01-19 20 / 19