Sparse Sampling and the Limits of Recovery

The sampling theorem of the first article in this cluster says that a signal with bandwidth $B$ needs $2B$ samples per second, and the argument is a counting one: the signal occupies a $B$-wide region of the frequency axis, and describing that region requires a matching number of independent numbers. Compressed sensing breaks the counting argument by changing the model. If the signal is known in advance to be a sum of a few complex exponentials, or to have a sparse representation in some basis, then the number of degrees of freedom is the number of nonzero coefficients rather than the bandwidth, and a number of samples proportional to that sparsity can suffice.

The catch, and it is a substantial one, is that the guarantee is not about the signal class alone. It also constrains how the samples are taken, and the reconstruction requires solving an optimization problem rather than evaluating a formula. Both conditions are easy to violate, and when violated the reconstruction can fail silently, producing a plausible-looking signal that is not the one that was sampled.

Why fewer samples can be enough

Consider a signal of length $N$ that is a sum of $K$ complex exponentials, with $K$ small. Its DFT has $K$ nonzero coefficients out of $N$. Suppose $M < N$ samples are taken at arbitrary distinct times, not on the uniform grid. The measurements are $M$ linear combinations of the $N$ unknown samples, so the system is underdetermined and has infinitely many solutions in general.

The sparsity prior selects among them. Of all signals consistent with the $M$ measurements, the one with the fewest nonzero coefficients is the intended answer, and for a well-designed measurement set that minimum is unique and equals the original. The counting condition is roughly $M \ge C K \log(N/K)$ for a constant $C$ that depends on the measurement ensemble: the samples needed scale with the sparsity times a logarithm, not with $N$.

The logarithm is what makes the idea interesting for large $N$. For $N = 10^6$ and $K = 100$, the bound is on the order of a few thousand measurements rather than a million. The improvement is asymptotic in the sparsity but only logarithmic in the ambient dimension, so it does not help when the signal is dense; a signal with $K$ comparable to $N$ needs a comparable number of measurements, and the reconstruction problem becomes as hard as the original.

The uniform-grid samples carry a specific structure that makes them the wrong choice here. The DFT of a length-$N$ signal sampled at $M$ of the $N$ grid points is not a subset of the full DFT: it is the full DFT aliased, with each retained coefficient being a sum of $N/M$ coefficients of the original. A sparse spectrum therefore becomes a dense aliased spectrum, and the sparsity that was supposed to help is destroyed by the sampling pattern itself. This is the same folding mechanism described in Sampling Rate, Aliasing, and the Folding of the Spectrum, and it is why compressed sensing uses random or pseudo-random sample times rather than a decimated uniform grid.

Incoherence is the design requirement

The measurement matrix $\Phi$ of size $M \times N$ maps the signal to its samples. Recovery works when $\Phi$ is incoherent with the basis in which the signal is sparse. Incoherence means that the rows of $\Phi$ have small inner products with the basis vectors, so each measurement carries a little bit of every coefficient rather than a lot of a few.

For a sparse-in-frequency signal and measurements at random times, incoherence follows from the fact that a random time point samples a random phase of every exponential. The measurement is a sum of all $K$ components with phases that are effectively random, so no single component dominates any measurement. The information about each coefficient is distributed across all $M$ measurements, which is why the recovery needs all of them jointly and why losing a random subset degrades gracefully rather than catastrophically.

A uniform-grid measurement is maximally coherent with the frequency basis: each grid point aligns with the exponentials in a structured way, and the aliasing described above is the concrete consequence. The design rule is therefore not “sample faster” but “sample at times that do not share structure with the basis”. Random sampling achieves that; a periodic pattern does not, and a jittered periodic pattern achieves it only when the jitter is large enough to break the periodicity.

import numpy as np

N, K, M = 512, 8, 64
rng = np.random.default_rng(7)

# Sparse signal: K random frequencies, random phases
freqs = rng.choice(N, K, replace=False)
coef = rng.standard_normal(K) + 1j * rng.standard_normal(K)
x = np.zeros(N, dtype=complex)
x[freqs] = coef

# Measure at M random times (rows of the inverse DFT at those times)
t = np.sort(rng.choice(N, M, replace=False))
F = np.exp(2j * np.pi * np.outer(t, np.arange(N)) / N)   # M x N
y = F @ x

# l1 minimization recovers x; the check below measures the measurement
# matrix restricted to the support instead of running the solver.
sub = F[:, freqs]
s = np.linalg.svd(sub, compute_uv=False)
print(s.min() / s.max())   # well above zero means the support is distinguishable

The ratio printed is a proxy for the conditioning of the measurement restricted to the true support; when it is close to one, the measurements distinguish the coefficients cleanly, and when it falls, the problem is ill-conditioned and noise is amplified.

The restricted isometry property

The formal condition that makes recovery work is that $\Phi$ acts almost like an isometry on every set of $K$ or fewer columns. That is, for every $K$-sparse vector $v$,

$$(1 - \delta_K)|v|_2^2 \le |\Phi v|_2^2 \le (1 + \delta_K)|v|_2^2$$

with a small constant $\delta_K$. If the property holds with a sufficiently small $\delta_{2K}$, then $\ell_1$ minimization recovers the sparsest solution exactly in the noiseless case, and with error proportional to the noise level in the noisy case.

The property is a statement about all $\binom{N}{K}$ supports simultaneously, which is why verifying it directly is impractical for realistic sizes: it requires checking every submatrix. What makes the theory usable is that random matrices satisfy it with high probability. A matrix with independent Gaussian entries, properly normalized, satisfies the property for $M$ on the order of $K \log(N/K)$; so do random partial Fourier matrices, which are the ones that correspond to sampling a signal at random times and evaluating its transform.

For a random partial Fourier matrix the guarantee is not quite as strong as for a Gaussian one, because the Fourier basis vectors have constant magnitude and the resulting matrix has heavier tails in its singular values. The required number of measurements carries an extra logarithmic factor, so the practical rule of thumb is $M$ several times $K$, with a factor of four a common starting point rather than a theorem.

minimize ||a||_1  subject to  ||Phi a - y||_2 <= epsilon

The optimization above is a convex problem, and it is the entire reconstruction algorithm. It is a linear program when written in real variables, or a second-order cone program in the complex case, and it can be solved by generic convex solvers at a cost polynomial in $N$ and $M$. Greedy alternatives, such as orthogonal matching pursuit, iteratively select the column most correlated with the residual and re-fit on the selected support; they are cheaper per iteration, but they have weaker guarantees and can fail on measurement matrices that satisfy the restricted isometry condition.

Where the theory stops being useful

Several failure modes are worth naming because they are common in applications and are not visible in the reconstruction.

The sparsity assumption may be approximate rather than exact. A signal that is a sum of a few components plus a small dense tail is not sparse, and the recovery error scales with the size of the tail relative to the sparse part. The theorem handles this in the sense that the error is bounded, but the bound grows with the tail’s energy, and a tail that is a few percent of the total can dominate the error at the sparsity levels where the method is otherwise attractive.

The basis may be wrong. Sparsity is a property of a representation, not of a signal. A piecewise-constant image is sparse in a wavelet basis and dense in a Fourier basis; applying compressed sensing with the wrong basis violates the incoherence premise and the recovery produces artifacts that look like ringing around edges. Choosing the basis is a modelling decision that the algorithm cannot make.

The measurement matrix may be constrained by hardware in ways that break incoherence. A magnetic resonance scanner cannot sample arbitrary times in $k$-space; it samples along trajectories, and the design of those trajectories is a separate problem with its own constraints on gradient slew rate and acquisition time. A coded-aperture camera cannot choose its measurement matrix freely either. The theory gives a target, and the hardware determines how close the achievable matrix is to it.

Noise interacts with the conditioning in a way that is easy to underestimate. The reconstruction’s error is proportional to the noise level divided by a quantity related to the restricted isometry constant, and that constant degrades as the measurement count approaches the theoretical minimum. Operating close to the minimum is therefore a regime with a large noise amplification factor, and the practical sample count is usually chosen with margin for that reason.

Relation to classical sampling

The two frameworks answer different questions. Classical sampling assumes a known band and guarantees recovery of every signal in that band with a fixed linear reconstruction. Sparse recovery assumes a known sparsity and guarantees recovery of the sparse signals with a nonlinear reconstruction whose cost is a convex program. The first gives a formula and a hard bandwidth limit; the second gives an optimization and a probabilistic guarantee whose constant depends on the measurement design.

The two can be combined. A signal that is bandlimited and also sparse within its band can be sampled below the band’s Nyquist rate if the sampling pattern is designed for the sparsity, which is the setting of the frequency-sparse results cited above. A signal that is sparse in a basis but not bandlimited requires the measurements to be taken at a rate that captures the basis’s structure, which for a Fourier basis is the same rate the classical theorem would demand. Sparsity reshapes the sample budget; it does not repeal the requirement that the measurements actually measure something.

The reconstruction also produces a signal defined on the full grid, so the output sample rate is a choice made by the reconstruction rather than by the acquisition. That is convenient for display and for downstream processing, but it obscures the fact that the samples taken were fewer than the output points. Any uncertainty statement attached to the output must account for the reconstruction’s assumptions, not merely the number of measurements, and the assumptions are what determine whether the output is a measurement or an inference. The same distinction appears in the design of anti-alias filters, where FIR Design by Windowing and the Gibbs Ripple sets how much of the out-of-band content the acquisition admits before any algorithm sees it.