Technical Summary: Approximating the Fourier Transform from Non-equispaced Discrete Samples

Problem Statement

The paper addresses the numerical approximation of the continuous Fourier transform, $\hat{f}(\xi) = \int_{-\infty}^{\infty} f(x)e^{-2\pi i x \xi} dx$, using a finite number of samples. While classical approaches rely on equispaced sampling on a truncated domain (often combined with sinc interpolation or the Fast Fourier Transform), these methods face inherent limitations. Specifically, equispaced grids necessitate truncating the infinite spatial domain, which introduces truncation errors and forces a compromise between spatial resolution near the origin and the coverage of asymptotic tails. This trade-off often leads to severe phase aliasing for high-frequency components and large pre-asymptotic error constants, particularly for functions with algebraic or heavy-tailed decay.

Existing literature, such as [4], has established optimal error rates for equispaced sampling by relating the number of samples, sampling interval, and grid size. However, the authors investigate whether non-equispaced sampling strategies can outperform these classical methods, particularly by adapting the sampling distribution to the specific decay properties of the function $f$ and its Fourier transform $\hat{f}$.

Methodology

The proposed approach replaces the standard equispaced grid with a deterministic non-equispaced sampling scheme driven by a probability density function $\mu$. The core methodology involves the following steps:

1. Change of Variables: The authors utilize the inverse cumulative distribution function (CDF), $\Psi^{-1}$, of a chosen probability density $\mu$ to map equispaced points $y_j = \frac{2j-1}{2M}$ from the unit interval $[0, 1]$ to the real line $\mathbb{R}$. The sampling nodes are defined as $x_j = \Psi^{-1}(y_j)$.
2. Transformed Quadrature: This transformation compactifies the infinite domain $\mathbb{R}$ onto $[0, 1]$, eliminating the need for explicit spatial truncation. The Fourier transform is approximated by a weighted exponential sum (a Riemann sum on the transformed domain):
   $$\hat{f}(\xi) \approx \tilde{f}(\xi) = \frac{1}{M} \sum_{j=1}^{M} \frac{f(x_j)}{\mu(x_j)} e^{-2\pi i x_j \xi}.$$
   This sum can be evaluated at arbitrary frequencies using the Nonequispaced Fast Fourier Transform (NNFFT).
3. Error Analysis Framework: The error analysis combines two main theoretical tools:
   • Aliasing Identity: An exact identity for the midpoint rule (Theorem 2.1) relates the quadrature error to the decay of the Fourier coefficients of the transformed integrand $g(y) = f(\Psi^{-1}(y)) \mu(\Psi^{-1}(y))^{-1} e^{-2\pi i \Psi^{-1}(y) \xi}$.
   • Stationary Phase Estimates: The decay of these Fourier coefficients is analyzed using van der Corput's lemma (Lemma 3.1) and stationary phase arguments. The analysis distinguishes between resonant cases (where the phase derivative vanishes) and non-resonant cases.

The authors study two specific classes of weight functions to define the decay of $f$ and $\hat{f}$:

• Polynomial Decay: $v(x) = (1+x^2)^{a/2}$ and $w(\xi) = (1+\xi^2)^{b/2}$.
• (Sub-)Exponential Decay: $v(x) = e^{r|x|^a}$ and $w(\xi) = e^{s|\xi|^b}$.

For polynomial decay, the sampling density $\mu$ is chosen as a polynomial distribution (related to the Student-t distribution):
$$\mu(x) = C_{\alpha, \sigma} \left(1 + \frac{x^2}{\sigma^2}\right)^{-\alpha/2}.$$
This density introduces two tunable parameters: the tail exponent $\alpha$ and the variance $\sigma$.

Key Contributions

1. Deterministic Non-equispaced Quadrature

The paper introduces a deterministic quadrature rule for the Fourier transform that avoids random sampling (unlike random Fourier features). By mapping equispaced points via $\Psi^{-1}$, the method guarantees deterministic error bounds rather than probabilistic estimates. The transformation naturally handles the infinite domain without truncation.

2. Comprehensive Error Estimates

The authors derive error bounds for every frequency $\xi \in \mathbb{R}$ and in $L^p(\mathbb{R})$ norms for all $p \geq 1$. Unlike previous works restricted to equispaced frequencies, this framework provides global error control. The error is decomposed into a truncation error (in frequency domain) and a discretization/aliasing error.

3. Optimization of Sampling Parameters

• Polynomial Decay: For functions with polynomial decay, the authors explicitly optimize the tail exponent $\alpha$ of the sampling density to balance the error terms. They show that with an optimal $\alpha$, the method recovers the asymptotic convergence rate of equispaced sampling, $O(M^{-(a-1/2)(b-1/2)/(a+b-1)})$. Crucially, the variance parameter $\sigma$ is shown to significantly reduce pre-asymptotic error constants without altering the asymptotic rate.
• (Sub-)Exponential Decay: For functions with exponential decay, the paper demonstrates that using a polynomial density (rather than an exponential one) yields (sub-)exponential error decay rates. This is because exponential densities can introduce non-smoothness in the transformed integrand, whereas polynomial densities preserve the necessary smoothness for high-order convergence.

4. Numerical Validation

Numerical experiments confirm the theoretical findings. The results demonstrate that the proposed non-equispaced method achieves the same asymptotic convergence rates as the optimal equispaced sampling from [4] but with substantially smaller pre-asymptotic errors. The method is particularly effective for functions with sharp peaks or heavy tails, where equispaced grids struggle.

Results

• Polynomial Decay Regime: The paper establishes that the mapped quadrature scheme can achieve the same optimal asymptotic error rate as equispaced sampling. However, by tuning the variance $\sigma$, the method reduces the pre-asymptotic error constants. In the "heavy-tail" regime where optimal balancing is impossible (due to integrability constraints on $\alpha$), the method still recovers the limiting rate of equispaced sampling.
• Exponential Decay Regime: The method achieves root-exponential or sub-exponential error decay rates, comparable to classical results for analytic functions. The use of a polynomial density with a tuned variance $\sigma$ allows for the reduction of approximation errors in practice, even if the theoretical exponent is slightly conservative compared to the ideal equispaced case.
• Pre-asymptotic Performance: In all tested scenarios, the non-equispaced method outperforms the equispaced benchmark (with sinc interpolation) in terms of absolute error for a fixed number of samples $M$, primarily due to the elimination of spatial truncation and the adaptive clustering of nodes.

Significance and Claims

The paper claims that its proposed method offers a robust alternative to classical equispaced sampling for Fourier transform approximation. The primary significance lies in the elimination of spatial truncation, which structurally captures the entire real line and maintains high resolution at the origin while integrating the tails to infinity.

The authors emphasize that while the asymptotic convergence rates for polynomial decay match those of equispaced sampling (as established in [4]), the pre-asymptotic error constants are substantially reduced through the optimization of the variance parameter $\sigma$. This makes the method more efficient in practical computational regimes where the number of samples $M$ is finite.

Furthermore, the paper highlights that the choice of a polynomial density is superior to an exponential density for functions with exponential decay, as it ensures the smoothness of the transformed integrand required for high-order convergence. The work provides a deterministic framework that complements existing equispaced theories, offering a unified approach to handle both algebraic and exponential decay with tailored sampling distributions.