# Embracing Off-the-Grid Samples

Oscar López and Özgür Yılmaz

###### Abstract

Many empirical studies suggest that samples of continuous-time signals taken
at locations randomly deviated from an equispaced grid (i.e., _off-the-grid_)
can benefit signal acquisition, e.g., undersampling and anti-aliasing.
However, explicit statements of such advantages and their respective
conditions are scarce in the literature. This paper provides some insight on
this topic when the sampling positions are known, with grid deviations
generated i.i.d. from a variety of distributions.

By solving the basis pursuit problem with an interpolation kernel we
demonstrate the capabilities of nonuniform samples for compressive sampling,
an effective paradigm for undersampling and anti-aliasing. For functions in
the Wiener algebra that admit a discrete $s$-sparse representation in some
transform domain, we show that $\mathcal{O}(s\operatorname{polylog}N)$ random
off-the-grid samples are sufficient to recover an accurate
$\frac{N}{2}$-bandlimited approximation of the signal. For sparse signals
(i.e., $s\ll N$), this sampling complexity is a great reduction in comparison
to equispaced sampling where $\mathcal{O}(N)$ measurements are needed for the
same quality of reconstruction (Nyquist-Shannon sampling theorem).

We further consider noise attenuation via oversampling (relative to a desired
bandwidth), a standard technique with limited theoretical understanding when
the sampling positions are non-equispaced. By solving a least squares problem,
we show that $\mathcal{O}(N\operatorname{polylog}N)$ i.i.d. randomly deviated
samples provide an accurate $\frac{N}{2}$-bandlimited approximation of the
signal with suppression of the noise energy by a factor
$\sim\frac{1}{\sqrt{\operatorname{polylog}N}}$.

##  1 Introduction

The Nyquist-Shannon sampling theorem is perhaps the most impactful result in
the theory of signal processing, fundamentally shaping the practice of
acquiring and processing data [7, 30] (also attributed to Kotel’nikov [80],
Ferrar [81], Cauchy [1], Ogura [42], E.T. and J.M. Whittaker [14, 37]). In
this setting, typical acquisition of a continuous-time signal involves taking
equispaced samples at a rate slightly higher than a prescribed frequency
$\omega$ Hz in order to obtain a bandlimited approximation via a quickly
decaying kernel. Such techniques provide accurate approximations of (noisy)
signals whose spectral energy is largely contained in the band
$[-\omega/2,\omega/2]$ [1, 2, 4, 22].

As a consequence, industrial signal acquisition and post-processing methods
tend to be designed to incorporate uniform sampling. However, such sampling
schemes are difficult to honor in practice due to physical constraints and
natural factors that perturb sampling locations from the uniform grid, i.e.,
_nonuniform_ or _off-the-grid_ samples. In response, nonuniform analogs of the
noise-free sampling theorem have been developed, where an average sampling
density proportional to the highest frequency $\omega$ of the signal
guarantees accurate interpolation, e.g., Landau density [32, 22, 43]. However,
non-equispaced samples are typically unwanted and regarded as a burden due to
the extra computational cost involved in regularization, i.e., interpolating
the nonuniform samples onto the desired equispaced grid.

On the other hand, many works in the literature have considered the potential
benefits of deliberate nonuniform sampling [31, 19, 20, 64, 67, 57, 53, 51,
78, 44, 6, 82, 10, 11, 71, 27, 62, 26, 48, 15, 61, 65, 60, 45, 40, 13, 33,
28]. Suppression of aliasing error, i.e., _anti-aliasing_ , is a well known
advantage of randomly perturbed samples. For example, jittered sampling is a
common technique for anti-aliasing that also provides a well distributed set
of samples [9, 63, 12, 62]. To the best of the authors’ knowledge, this
phenomenon was first noticed by Harold S. Shapiro and Richard A. Silverman
[31] (also by Federick J. Beutler [19, 20] and implicitly by Henry J. Landau
[32]) and remained unused in applications until rediscovered in Pixar
animation studios by Robert L. Cook [64]. According to our literature review,
such observations remain largely empirical or arguably uninformative for
applications. Closing this gap between theory and experiments would help the
practical design of such widely used methodologies.

To this end, in this paper we propose a practical framework that allows us to
concretely investigate the properties of randomly deviated samples for
undersampling, anti-aliasing and general noise attenuation. To elaborate (see
Section 1.1 for notation), let
$\textbf{f}:[-\frac{1}{2},\frac{1}{2})\mapsto\mathbb{C}$ be our function of
interest that belongs to some smoothness class. Our goal is to obtain a
uniform discretization $f\in\mathbb{C}^{N}$, where
$f_{k}=\textbf{f}(\frac{k-1}{N}-\frac{1}{2})$ will provide an accurate
approximation of f on $[-\frac{1}{2},\frac{1}{2})$. We are given noisy non-
equispaced samples, $b=\tilde{f}+d\in\mathbb{C}^{m}$, where
$\tilde{f}_{k}=\textbf{f}(\frac{k-1}{m}-\frac{1}{2}+\Delta_{k})$ is the
nonuniformly sampled signal and $d\in\mathbb{C}^{m}$ encompasses unwanted
additive noise. In general, we will consider functions f with support on
$[-\frac{1}{2},\frac{1}{2})$ whose periodic extension is in the Wiener algebra
$A(\Omega)$ [83], where by abuse of notation $\Omega$ denotes the interval
$[-\frac{1}{2},\frac{1}{2})$ and the torus $\mathbb{T}$.

To achieve undersampling and anti-aliasing, we assume our uniform signal
admits a sparse (or compressible) representation along the lines of
compressive sensing [36, 73, 29]. We say that $f$ is _compressible_ with
respect to a transform, $\Psi\in\mathbb{C}^{N\times n}$, if there exists some
$g\in\mathbb{C}^{n}$ such that $f=\Psi g$ and $g$ can be accurately
approximated by an $s$-sparse vector ($s\leq n$). In this scenario, our
methodology consists of constructing an interpolation kernel
$\mathcal{S}\in\mathbb{R}^{m\times N}$ that achieves
$\mathcal{S}f\approx\tilde{f}$ accurately for smooth signals, in order to
approximate $\Psi g^{\sharp}\approx f$ where

$g^{\sharp}:=\operatorname{argmin}_{h\in\mathbb{C}^{n}}\|h\|_{1}\ \
\mbox{s.t.}\ \ \|\mathcal{S}\Psi h-b\|_{2}\leq\tilde{\eta},$ (1)

and $\tilde{\eta}\geq 0$ is a parameter to be chosen appropriately. We show
that for signals in the Wiener algebra and under certain distributions
$\mathcal{D}$, $m\sim\mathcal{O}(s\operatorname{polylog}(n))$ off-the-grid
samples with i.i.d. deviations $\\{\Delta_{k}\\}_{k=1}^{m}\sim\mathcal{D}$
supply $\|f-\Psi g^{\sharp}\|_{2}$ proportional to $\|d\|_{2}$, the error of
the best $s$-sparse approximation of $g$, and the error of the best
$\frac{N}{2}$-bandlimited approximation of f in the Wiener algebra norm (see
equation 6.1 in [83]). If $s\ll N$, the average sampling rate required for our
result (step size $\sim\frac{1}{s\operatorname{polylog}(n)}$) provides a stark
contrast to standard density conditions where a rate proportional to the
highest frequency $\omega\sim N$, resulting in step size $\sim\frac{1}{N}$, is
needed for the same bandlimited approximation. The result is among the first
to formally state the anti-aliasing nature of nonuniform sampling in the
context of undersampling (see Section 3).

Removing the sparse signal model, we attenuate measurement noise (i.e.,
_denoise_) via the estimate

$f\approx
f^{\sharp}:=\operatorname{argmin}_{h\in\mathbb{C}^{N}}\|\mathcal{S}h-b\|_{2}.$
(2)

In this context, our main result states that
$m\sim\mathcal{O}(N\operatorname{polylog}N)$ i.i.d. randomly deviated samples
provide recovery error $\|f-f^{\sharp}\|_{2}$ proportional to the noise level
($\frac{\sqrt{N}}{\sqrt{m}}\|d\|_{2}$) and the error of the best
$\frac{N}{2}$-bandlimited approximation of f in the Wiener algebra norm. Thus,
by nonuniform oversampling $m\geq N$ (relative to the desired
$\frac{N}{2}$-bandwidth) we attenuate unwanted noise regardless of its
structure. While uniform oversampling is a common noise filtration technique,
our results show that general nonuniform samples also posses this denoising
property (see Section 4). The methodology and proofs of Section 4 can be
modified to account for denoising by repeating acquisitions, a common approach
related to signal averaging that increases the signal-to-noise ratio [79].

The rest of the paper is organized as follows: Section 2 provides a detailed
elaboration of our sampling scenario, signal model and methodology. Section 3
showcases our results for sampling of sparse signals and Section 4 considers
noise attenuation via oversampling. A comprehensive discussion of the results
and implications is presented in Section 5. Several experiments and
computational considerations are found in Section 6, followed by concluding
remarks in Section 7. We postpone the proofs of our statements until Section
8. Before proceeding to the next section, we find it best to introduce the
general notation that will be used throughout. However, each subsequent
section may introduce additional notation helpful in its specific context.

###  1.1 Notation

We denote complex valued functions of real variables using bold letters, e.g.,
$\textbf{f}:\mathbb{R}^{2}\to\mathbb{C}$. For any integer $n\in\mathbb{N}$,
$[n]$ denotes the set $\\{\ell\in\mathbb{N}:1\leq\ell\leq n\\}$. For
$k,\ell\in\mathbb{N}$, $b_{k}$ indicates the $k$-th entry of the vector $b$,
$D_{k\ell}$ denotes the $(k,\ell)$ entry of the matrix $D$ and $D_{k*}\
(D_{*\ell})$ denotes the $k$-th row (resp. $\ell$-th column) of the matrix. We
reserve $x,y$ to denote real variables and write the complex exponential as
$\textbf{e}(x):=e^{2\pi ix}$, where $i$ is the imaginary unit. For a vector
$f\in\mathbb{C}^{n}$ and $1\leq p<\infty$,
$\|f\|_{p}:=\big{[}\sum_{k=1}^{n}\lvert f_{k}\rvert^{p}\big{]}^{1/p}$ is the
$p$-norm, $\|f\|_{\infty}=\max_{k\in[n]}\lvert f_{k}\rvert$, and $\|f\|_{0}$
gives the total number of non-zero elements of $f$. For a matrix
$X\in\mathbb{C}^{n\times m}$, $\sigma_{k}(X)$ denotes the $k$-th largest
singular value of $X$ and $\|X\|:=\sigma_{1}(X)$ is the spectral norm.
$A(\Omega)$ is the Wiener algebra and $H^{k}(\Omega)$ is the Sobolev space
$W^{k,2}(\Omega)$ (with domains $\Omega$), $S^{n-1}$ is the unit sphere in
$\mathbb{C}^{n}$, and the adjoint of an operator $\mathcal{A}$ is denoted by
$\mathcal{A}^{*}$.

##  2 Assumptions and Methodology

In this section we develop the signal model, deviation model and interpolation
kernel, in Sections 2.1, 2.2 and 2.3 respectively. This will allow us to
proceed to Sections 3 and 4 where the main results are elaborated. We note
that the deviation model (Section 2.2) and sparse signal model at the end of
Section 2.1 only apply to the compressive sensing results in Section 3.
However, the _sampling on the torus_ assumption in Section 2.2 does apply to
the results in Section 4 as well.

###  2.1 Signal Model

For the results in Section 3 and 4, let $\Omega=[-\frac{1}{2},\frac{1}{2})$
and let $\textbf{f}:\Omega\to\mathbb{C}$ be the function of interest to be
sampled. We assume $\textbf{f}\in A(\Omega)$ with Fourier expansion

$\textbf{f}(x)=\sum_{\ell=-\infty}^{\infty}c_{\ell}\textbf{e}(\ell x),$ (3)

on $\Omega$. Note that our regularity assumption implies that

$\sum_{\ell=-\infty}^{\infty}\lvert c_{\ell}\rvert<\infty,$

which will be crucial for our error bound. Further, $H^{k}(\Omega)\subset
A(\Omega)$ for $k\geq 1$ so that our context applies to many signals of
interest.

Henceforth, let $N\in\mathbb{N}$ be odd. We denote the discretized regular
data vector by $f\in\mathbb{C}^{N}$, which is obtained by sampling f on the
uniform grid $\\{t_{1},\cdots,t_{N}\\}\subset\Omega$, with
$t_{k}:=\frac{k-1}{N}-\frac{1}{2}$, (which is a collection of equispaced
points) so that $f_{k}=\textbf{f}(t_{k})$. Here, $f$ will be our discrete
signal of interest to recover via few nonuniform samples. Similar results can
be obtained in the case $N$ is even, our current assumption is adopted to
simplify the exposition.

The observed nonuniform data is encapsulated in the vector
$\tilde{f}\in\mathbb{C}^{m}$ with underlying unstructured grid
$\tilde{\tau}=\\{\tilde{t}_{1},\cdots,\tilde{t}_{m}\\}\subset\Omega$ where
$\tilde{t}_{k}:=\frac{k-1}{m}-\frac{1}{2}+\Delta_{k}$ is now a collection of
generally non-equispaced points. The entries of the perturbation vector
$\Delta\in\mathbb{R}^{m}$ define the pointwise deviations of $\tilde{\tau}$
from the equispaced grid $\\{\frac{k-1}{m}-\frac{1}{2}\\}_{k=1}^{m}$, where
$\tilde{f}_{k}=\textbf{f}(\tilde{t}_{k})$. We assume that the nonuniform
points remain in our sampling space, i.e., $\tilde{\tau}\subset\Omega$ so that
(3) remains valid for these measurements. See the _sampling on the torus_
premise in Section 2.2 for a detailed discussion of this crucial restriction.

In our sampling scenario, noisy nonuniform samples are given as

$b=\tilde{f}+d\in\mathbb{C}^{m},$

where the noise model, $d$ with $\|d\|_{2}\leq\eta$, does not incorporate off-
the-grid corruption. We assume that we know $\tilde{\tau}$ and $\eta$.

_Sparse signal model:_ For the results of Section 3 only, we impose a
compressibility condition on $f\in\mathbb{C}^{N}$. To this end, let
$\Psi\in\mathbb{C}^{N\times N}$ be a basis with $0<\sigma_{N}(\Psi)=:\alpha$
and $\sigma_{1}(\Psi)=:\beta$. We assume there exists some
$g\in\mathbb{C}^{N}$ such that $f=\Psi g$, where $g$ can be accurately
approximated by an $s\leq N$ sparse vector. To be precise, for $s\in[N]$ we
define the error of best $s$-sparse approximation of $g$ as

$\epsilon_{s}(g):=\min_{\|h\|_{0}\leq s}\|h-g\|_{1},$

and assume $s$ has been chosen so that $\epsilon_{s}(g)$ is within a
prescribed error tolerance determined by the practitioner. In Section 8.1, we
will relax the condition that $\Psi$ be a basis by allowing full column rank
matrices. While such transforms are not typical in compressive sensing, we
argue that they may be of practical interest since our results will show that
if $\Psi$ can be selected as tall matrix then the sampling complexity will
solely depend on its number of columns (i.e., the smallest dimension).
Further, we require $m\leq N$ (in Section 3 only, needed in the proof of
Theorem 7, see Lemma 3.6 from [49]).

The transform $\Psi$ will have to be incoherent with respect to the 1D
centered discrete Fourier basis $\mathcal{F}\in\mathbb{C}^{N\times N}$ (see
Section 2.3 for definition of $\mathcal{F}$). We define the DFT-incoherence
parameter as

$\gamma=\max_{\ell\in[N]}\sum_{k=1}^{N}\lvert\langle\mathcal{F}_{*k},\Psi_{*\ell}\rangle\rvert,$

which provides a uniform bound on the $\ell_{1}$-norm of the discrete Fourier
coefficients of the columns of $\Psi$. This parameter will play a role in the
sampling complexity of our result in Section 3. We discuss $\gamma$ in detail
in Section 5.3, including examples for several transforms common in
compressive sensing.

###  2.2 Deviation Model

Section 3 will apply to random deviations $\Delta\in\mathbb{R}^{m}$ with
$m\leq N$ and whose entries are i.i.d. with any distribution, $\mathcal{D}$,
that obeys the following: for $\delta\sim\mathcal{D}$, there exists some
$\theta\geq 0$ such that for all integers $0<\lvert
j\rvert\leq\frac{2(N-1)}{m}$ we have

$\frac{2N}{m}\lvert\mathbb{E}\textbf{e}(jm\delta)\rvert\leq\theta.$

This will be known as our _deviation model_. In our results, distributions
with smaller $\theta$ parameter will require less samples and provide reduced
error bounds. We postpone further discussion of the deviation model until
Section 5.2, where we will also provide examples of deviations that fit our
criteria. We note that the case $m>N$ is discussed in Section 4, which no
longer requires this deviation model.

_Sampling on the torus:_ For the results in Section 3 and 4, it is important
to notice that the deviations can allow for nonuniform samples outside of
$\Omega:=[-\frac{1}{2},\frac{1}{2})$. This violates our assumption in Section
2.1, as (3) will no longer hold for $\tilde{\tau}$. To remedy this, we modify
our sampling scheme to be on the torus. In other words, if
$\textbf{f}\lvert_{\Omega}(x)$ is given as

$\textbf{f}\lvert_{\Omega}(x)=\begin{cases}\textbf{f}(x)&\text{if}\ \
x\in[-\frac{1}{2},\frac{1}{2})\\\ 0&\text{if}\ \
x\notin[-\frac{1}{2},\frac{1}{2}),\end{cases}$

then we define $\tilde{\textbf{f}}(x)$ as the periodic extension of
$\textbf{f}\lvert_{\Omega}(x)$ to the whole line

$\tilde{\textbf{f}}(x)=\sum_{\ell=-\infty}^{\infty}\textbf{f}\lvert_{\Omega}(x+\ell).$

We now apply samples generated from our deviations to $\tilde{\textbf{f}}(x)$.
Indeed, for any $\tilde{t}_{k}$ generated outside of $\Omega$ will have
$\tilde{\textbf{f}}(\tilde{t}_{k})=\textbf{f}(t^{*})$ for some
$t^{*}\in\Omega$. In this way we may proceed with our sampling scenario, which
will provide nonuniform samples for f in $\Omega$ (in particular the expansion
(3) will remain valid for these samples).

###  2.3 Dirichlet Kernel

In Section 3 and 4, we model our nonuniform samples via an interpolation
kernel $\mathcal{S}\in\mathbb{R}^{m\times N}$ that achieves
$\mathcal{S}f\approx\tilde{f}$ accurately. We consider the Dirichlet kernel
defined by
$\mathcal{S}=\mathcal{N}\mathcal{F}^{*}:\mathbb{C}^{N}\to\mathbb{C}^{m}$,
where $\mathcal{F}\in\mathbb{C}^{N\times N}$ is a 1D centered discrete Fourier
transform (DFT) and $\mathcal{N}\in\mathbb{C}^{m\times N}$ is a 1D centered
nonuniform discrete Fourier transform (NDFT, see [38, 46]) with normalized
rows and non-harmonic frequencies chosen according to $\tilde{\tau}$. In other
words, let $\tilde{N}=\frac{N-1}{2}$, then the $(k,\ell)\in[m]\times[N]$ entry
of $\mathcal{N}$ is given as

$\mathcal{N}_{k\ell}=\frac{1}{\sqrt{N}}\textbf{e}(-\tilde{t}_{k}(\ell-\tilde{N}-1)).$

This NDFT is referred to as a nonuniform discrete Fourier transform of type 2
in [46]. Thus, the action of $\mathcal{S}$ on $f\in\mathbb{C}^{N}$ can be
given as follows: we first apply the centered inverse DFT to our discrete
uniform data

$\displaystyle\check{f}_{u}:=(\mathcal{F}^{*}f)_{u}=\sum_{p=1}^{N}f_{p}\mathcal{F}^{*}_{up}:=\frac{1}{\sqrt{N}}\sum_{p=1}^{N}f_{p}\textbf{e}(t_{p}(u-\tilde{N}-1)),\
\ \forall u\in[N],$ (4)

followed by the NDFT in terms of $\tilde{\tau}$:

$\displaystyle(\mathcal{S}f)_{k}:=(\mathcal{N}\check{f})_{k}=\sum_{u=1}^{N}\check{f}_{u}\mathcal{N}_{ku}:=\frac{1}{\sqrt{N}}\sum_{u=1}^{N}\check{f}_{u}\textbf{e}(-\tilde{t}_{k}(u-\tilde{N}-1)),\
\ \forall k\in[m].$ (5)

Equivalently,

$\displaystyle(\mathcal{S}f)_{k}=\frac{1}{N}\sum_{p=1}^{N}f_{p}\textbf{K}(\tilde{t}_{k}-t_{p})\
\ \ \ \forall k\in[m],$ (6)

where $\textbf{K}(\theta)=\frac{\sin{(N\pi\theta})}{\sin{(\pi\theta)}}$ is the
Dirichlet kernel. This equality is well known and holds by applying the
geometric series formula upon expansion (notice that $\mathcal{S}$ is real
valued). This kernel is commonly used for trigonometric interpolation and is
therefore accurate when acting on signals that can be well approximated by
trigonometric polynomials of finite order, as we show in the following
theorem.

######  Theorem 1.

Let $\mathcal{S},f$ and $\tilde{f}$ be defined as above with
$\tilde{\tau}\subset\Omega$. For each $k\in[m]$, if $\tilde{t}_{k}=t_{p}$ for
some $p\in[N]$ then

$\left(\tilde{f}-\mathcal{S}f\right)_{k}=0$ (7)

and otherwise

$\left(\tilde{f}-\mathcal{S}f\right)_{k}=\sum_{\lvert\ell\rvert>\tilde{N}}c_{\ell}\left(\emph{{e}}(\ell\tilde{t}_{k})-(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}\emph{{e}}(r(\ell)\tilde{t}_{k})\right),$
(8)

where $r(\ell)=\operatorname{rem}(\ell+\tilde{N},N)-\tilde{N}$ with
$\operatorname{rem}(\ell,N)$ giving the remainder after division of $\ell$ by
$N$. As a consequence, for any $1\leq p<\infty$

$\displaystyle\|\tilde{f}-\mathcal{S}f\|_{p}\leq
2m^{\frac{1}{p}}\sum_{\lvert\ell\rvert>\tilde{N}}\lvert c_{\ell}\rvert,$ (9)

and

$\displaystyle\|\tilde{f}-\mathcal{S}f\|_{\infty}\leq
2\sum_{\lvert\ell\rvert>\tilde{N}}\lvert c_{\ell}\rvert.$ (10)

The proof of this theorem is postponed until Section 8.3. Therefore, the error
due to $\mathcal{S}$ is proportional to the $1$-norm (or Wiener algebra norm)
of the Fourier coefficients of f that correspond to frequencies larger than
$\tilde{N}=\frac{N-1}{2}$. In particular notice that if $c_{\ell}=0$ for all
$\ell>\tilde{N}$ we obtain perfect interpolation, as expected from standard
results in signal processing (i.e., bandlimited signals consisting of
trigonometric polynomials with finite degree $\leq\tilde{N}$). Despite the
wide usage of trigonometric interpolation in applications [74, 16, 34], such a
result that gives a sharp error bound does not seem to exist in the
literature.

Notice that Theorem 1 only holds for $\tilde{\tau}\subset\Omega$ as restricted
in Section 2.1. However, the results continues to hold for $\tilde{\tau}$
unrestricted if we sample on the torus as imposed in Section 2.2. Therefore,
the error bound will always hold under our setup.

##  3 Compressive Nonuniform Sampling

With the definitions and assumptions introduced in Section 2, our methodology
in this chapter will consist of modeling our nonuniform measurements via
$\mathcal{S}$ and approximating the $s$ largest coefficients of $f$ in $\Psi$
(in the representation $f=\Psi g$).

The following is a simplified statement, assuming that $\Psi$ is an
orthonormal basis. We focus on this cleaner result for ease of exposition,
presented as a corollary of the main result in Section 8.1. The full statement
considers more general and practical $\Psi$ that also allow for reduced sample
complexity.

######  Theorem 2.

Under our signal model with Fourier expansion (3), let
$\Psi\in\mathbb{C}^{N\times N}$ be an orthonormal basis with DFT-incoherence
parameter $\gamma$. Define the interpolation kernel $\mathcal{S}$ as in
Section 2.3 with the entries of $\Delta$ i.i.d. from any distribution
satisfying our deviation model from Section 2.2 with
$\theta<\frac{1}{\sqrt{2}}$. Let

$g^{\sharp}=\operatorname{argmin}_{h\in\mathbb{C}^{N}}\|h\|_{1}\ \
\mbox{s.t.}\ \ \|\mathcal{S}\Psi
h-b\|_{2}\leq\eta+2\sqrt{m}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert,$ (11)

where $m$ is the number of nonuniform samples. If

$m\geq\frac{C_{1}\gamma^{2}s\log^{4}(N)\left(1+\theta+\tau\right)}{\tau^{2}}$
(12)

where $\tau:=\frac{1}{4\sqrt{2}}-\frac{\theta}{4}>0$ and $C_{1}$ is an
absolute constant, then

$\|f-\Psi
g^{\sharp}\|_{2}\leq\frac{C_{2}\epsilon_{s}(g)}{\sqrt{s}}+\frac{C_{3}\sqrt{N}\eta}{\sqrt{m}}+2C_{3}\sqrt{N}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ (13)

with $C_{2},C_{3}>0$ depending only on $\theta$ and probability exceeding

$1-\exp\left(-\frac{m\tau}{8s\gamma^{2}}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+1+2\theta}\right)\right)\right).$
(14)

Therefore, with $m\sim s\log^{4}(N)$ randomly perturbed samples we can recover
$f$ with error (13) proportional to the sparse model mismatch
$\epsilon_{s}(g)$, the noise level $\eta$, and the error of the best
$\frac{N-1}{2}$-bandlimited approximation of f in the Wiener algebra norm
(i.e., $\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert c_{\ell}\rvert$). As a
consequence, we can approximate $\textbf{f}(x)$ for all $x\in\Omega$ as stated
in the following corollary.

######  Corollary 1.

Let $\emph{{h}}:\Omega\to\mathbb{C}^{N}$ be the vector valued function defined
entry-wise for $\ell\in[N]$ as

$\emph{{h}}(x)_{\ell}:=\frac{1}{\sqrt{N}}\textbf{e}(-x(\ell-\tilde{N}-1)),$

and define the function $\emph{{f}}^{\sharp}:\Omega\to\mathbb{C}$ via

$\emph{{f}}^{\sharp}(x)=\langle\emph{{h}}(x),\mathcal{F}^{*}\Psi
g^{\sharp}\rangle,$ (15)

where $g^{\sharp}$ is given by (11).

Then, under the assumptions of Theorem 2,

$\lvert\emph{{f}}(x)-\emph{{f}}^{\sharp}(x)\rvert\leq\frac{C_{2}\epsilon_{s}(g)}{\sqrt{s}}+\frac{C_{3}\sqrt{N}\eta}{\sqrt{m}}+2(C_{3}\sqrt{N}+1)\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ (16)

holds for all $x\in\Omega=[-\frac{1}{2},\frac{1}{2})$ with probability
exceeding (14).

The proof of this corollary is presented in Section 8.3. In the case
$\epsilon_{s}(g)=\eta=0$, the result intuitively says that we can recover a
$\frac{N-1}{2}$-bandlimited approximation of f with
$\mathcal{O}(s\operatorname{polylog}(N))$ random off-the-grid samples. In the
case of equispaced samples, $\mathcal{O}(N)$ measurements are needed for the
same quality of reconstruction by the Nyquist-Shannon sampling theorem (or by
Theorem 1 directly). Thus, for compressible signals with $s\ll N$, random
nonuniform samples provide a significant reduction in sampling complexity
(undersampling) and simultaneously allow recovery of frequency components
exceeding the sampling density (anti-aliasing). See Section 5 for further
discussion.

The term $2\sqrt{m}\sum_{\lvert\ell\rvert>\tilde{N}}\lvert c_{\ell}\rvert$ in
the constraint of (11) is due to the error of our kernel $\mathcal{S}$ given
in (9). The introduction of this interpolation error may render our
methodology cumbersome, since a practitioner may not be supplied with this
term or an appropriate bound. We offer the following proposition to illustrate
that our approach is still applicable when this quantity is not known.

######  Proposition 3.

Let $\Psi\in\mathbb{C}^{N\times N}$ be an orthonormal basis with DFT-
incoherence parameter $\gamma$. Define $\mathcal{S}$ as in (6) with each
sampling point $\\{\tilde{t}_{k}\\}_{k\in[m]}$ chosen uniformly at random from
the set

$\Bigg{\\{}\frac{\ell-1}{N}-\frac{1}{2}\Bigg{\\}}_{\ell=1}^{N}$

and

$g^{\sharp}=\operatorname{argmin}_{h\in\mathbb{C}^{N}}\|h\|_{1}\ \
\mbox{s.t.}\ \ \|\mathcal{S}\Psi h-b\|_{2}\leq\eta.$ (17)

If (12) holds with $\theta=0$, then $\emph{{f}}^{\sharp}$ defined as in (15)
satisfies

$\lvert\emph{{f}}(x)-\emph{{f}}^{\sharp}(x)\rvert\leq\frac{C_{4}\epsilon_{s}(g)}{\sqrt{s}}+\frac{C_{5}\sqrt{N}\eta}{\sqrt{m}}+2\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ (18)

for all $x\in\Omega=[-\frac{1}{2},\frac{1}{2})$ with absolute constants
$C_{4},C_{5}>0$ and probability exceeding (14).

This approach removes the need to estimate
$\sum_{\lvert\ell\rvert>\tilde{N}}\lvert c_{\ell}\rvert$ for implementation,
due to Theorem 1 which shows that samples chosen from this dense uniform grid
no have interpolation error (7). The sampling model of this statement does not
obey our deviation model, so a brief proof using our main result does not seem
possible. For the sake of brevity, we will not prove this proposition, but
note that it can be established by a method of proof very similar to that of
Corollary 1 and Theorem 2. In fact, the sampling model of Proposition 3 is
considered by [34] and the proofs therein can be set within our context to
establish the claim (see Section 5.1 for further discussion on this related
work).

Notice that general denoising is not guaranteed via our undersampling scenario
($m\leq N$), due to the term $\frac{\sqrt{N}\eta}{\sqrt{m}}$ in (13), (16) and
(18). In other words, one cannot expect to reduce the measurement noise $\eta$
since $\frac{\sqrt{N}}{\sqrt{m}}\geq 1$ appearing in our error bound implies
an amplification of the input noise level. This observation is demonstrated
numerically in Section 6.3. In general a practitioner must oversample (i.e.,
$N<m$) to attenuate the effects of noise. However, the main result states that
nonuniform samples specifically attenuate aliasing noise.

##  4 Denoising via Nonuniform Oversampling

In this section, we show that reduction in the noise level introduced during
acquisition is possible given nonuniform samples whose average density exceeds
the Nyquist rate (relative to a desired bandwidth). While the implications of
this section are not surprising in the context of classical sampling theory,
to the best of our knowledge such guarantees do not exist in the literature
when the sampling points are nonuniform.

By removing the sparse signal model (Section 2.1), deviation model (Section
2.2), and requiring $m\geq N$, we now use the numerically cheaper program of
least squares. To reiterate, $\textbf{f}\in A(\Omega)$ with Fourier expansion
$\sum_{\ell=-\infty}^{\infty}c_{\ell}\textbf{e}(\ell x)$ is our continuous
signal of interest. With $N$ odd, $f\in\mathbb{C}^{N}$ is the discrete signal
to be approximated, where $f_{k}=\textbf{f}(t_{k})$ for
$t_{k}:=\frac{k-1}{N}-\frac{1}{2}$. $\tilde{f}\in\mathbb{C}^{m}$ encapsulates
the nonuniformly sampled data where $\tilde{f}_{k}=\textbf{f}(\tilde{t}_{k})$
for $\tilde{t}_{k}:=\frac{k-1}{m}-\frac{1}{2}+\Delta_{k}$. Noisy nonuniform
samples are given as

$b=\tilde{f}+d\in\mathbb{C}^{m},$

where the additive noise model, $d$, does not incorporate off-the-grid
corruption.

A similar method of proof to obtain Theorem 2 provides the following denoising
result for a more general set of deviations.

######  Theorem 4.

Let the entries of $\Delta$ be i.i.d. from any distribution and define

$f^{\sharp}=\operatorname{argmin}_{g\in\mathbb{C}^{N}}\|\mathcal{S}g-b\|_{2}.$
(19)

If for some $0<\tau<\frac{1}{2}$

$m\geq\frac{36N\log(2N)\left(1+\tau\right)}{\tau^{2}},$ (20)

then

$\displaystyle\|f-f^{\sharp}\|_{2}\leq\frac{\tau\|d\|_{2}}{3\sqrt{(1-\tau-2\tau^{2})\log(2N)}}+\frac{2\sqrt{N}}{\sqrt{1-2\tau}}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ (21)

with probability exceeding

$1-\exp\left(-\frac{m\tau}{4N}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+1}\right)\right)\right).$

The proof can be found in Section 8.2. In this scenario, we oversample
relative to the $\frac{N-1}{2}$-bandlimited output by generating a set samples
with average density exceeding the Nyquist rate (average step size
$\frac{1}{N}$). With $m\sim N\log(N)$, bound (21) tells us that we can
diminish the measurement noise level $\|d\|_{2}$ by a factor
$\sim\frac{\tau}{\sqrt{\log(N)}}$. Notice that the parameter $\tau$ may be
varied for increased attenuation at the cost of denser sampling. In
particular, the deviations $\Delta$ need not be from our deviation model in
Section 2.2 and instead the result applies to perturbations generated by any
distribution. This includes the degenerate distribution (deterministic), so
the claim also holds in the case of equispaced samples. Furthermore, we no
longer require the sparsity assumption and the result applies to all functions
in the Wiener algebra. Finally, the recovery method (19) consists of standard
least squares which can be solved cheaply relative to the basis pursuit
problem (11) from the previous section.

We may proceed analogously to Corollary 1 and show that the recovered discrete
signal provides a continuous approximation to $\textbf{f}(x)$ for all
$x\in\Omega$. The error is proportional to the error of the best
$\frac{N-1}{2}$-bandlimited approximation in the Wiener algebra norm while
attenuating the introduced measurement noise. In the result, the structure of
the deviated samples is quite general and accounts for many practical cases.

While related results exist in the equispaced case (see for example Section 4
of [4]), Theorem 4 is the first such statement in a general non-equispaced
case. The result therefore provides insight into a widely applied techniques
for the removal of unwanted noise, without making any assumptions on the noise
structure.

_Repeated acquisitions:_ We end this section by noting that our setup can be
modified to obtain a similar conclusion for noise reduction by repeating
signal acquisition, a methodology related to signal averaging [79]. The
discussion is informal and we do not provide proofs to our statements for
brevity. However, the results can be obtained by analogous arguments used to
establish Theorem 4 and Theorem 2.

In this context, practitioners repeat acquisition several times to obtain a
sequence of independent and corrupted observations
$b^{1},b^{2},\cdots,b^{P}\in\mathbb{C}^{m}$. Each $b^{\ell}$ is of the form

$b^{\ell}=\tilde{f}^{\ell}+d^{\ell}\in\mathbb{C}^{m},$

where $\\{\tilde{f}^{\ell}\\}_{\ell=1}^{P}$ and $\\{d^{\ell}\\}_{\ell=1}^{P}$
are $P$ independently generated copies of $\tilde{f}$ and $d$ (respectively)
from our previous set up. In other words, each of the acquisitions samples f
at independently generated grids and introduces noise (we do not assume the
noise terms are distinct or independent).

This situation also corresponds to oversampling, but rather than probing on a
dense grid we repeat our acquisition on a fixed grid (not necessarily dense)
with perturbations. We proceed by forming $\mathcal{S}\in\mathbb{R}^{mP\times
N}$ as a concatenation of Dirichlet kernels
$\mathcal{S}^{1},\mathcal{S}^{2},\cdots,\mathcal{S}^{P}\in\mathbb{R}^{m\times
N}$ constructed with the nonuniform grids of their respective acquisition
(i.e., so that $\mathcal{S}^{\ell}f\approx\tilde{f}^{\ell}$). Similarly,
construct $b\in\mathbb{C}^{mP}$ by concatenating $b^{1},b^{2},\cdots,b^{P}$.
For simplicity, assume $\|d^{\ell}\|_{2}\leq\eta$ for all $\ell\in[P]$.

In the case $m\geq 2N$ with i.i.d. deviations and $0<\tau<1/2$, a proof
analogous to that of Theorem 4 gives that $P\geq
18\log(N)\frac{1+\tau}{\tau^{2}}$ repeated acquisitions provide $f^{\sharp}$
from (19) satisfying

$\displaystyle\|f-f^{\sharp}\|_{2}\leq\frac{\tau\eta}{3\sqrt{2(1-\tau-2\tau^{2})}}+\frac{2\sqrt{N}}{\sqrt{1-2\tau}}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ (22)

with high probability. The result implies recovery of the
$\frac{N-1}{2}$-bandlimited approximation of f with reduction of the noise
level by a factor of $\frac{\tau}{3\sqrt{2(1-\tau-2\tau^{2})}}$ which is less
than 1 if $\tau<\frac{1}{2.02}$. Once again, the parameter $\tau$ allows the
practitioner to control the noise attenuation level at the cost of more
acquisitions.

The case $m<2N$ can also be handled, but requires us to adopt our deviation
model from Section 2.2 with parameter $\theta$. In the simplest case, we can
consider distributions with $\theta=0$ (see Section 5.2). We obtain that
$P\geq 36\frac{N}{m}\log(N)\frac{1+\tau}{\tau^{2}}$ repeated acquisitions
provide a similar error bound but with the noise term in (22) replaced by
$\frac{\tau\sqrt{N}\eta}{3\sqrt{m(1-\tau-2\tau^{2})}}$. In particular, this
case of sampling on a less dense grid still allows for denoising if roughly
$\tau<\frac{\sqrt{m}}{\sqrt{N}}$.

This setup is related to signal averaging for noise removal. In this practice,
the corrupted observations $b^{\ell}$ are averaged in order to amplify
components due to the desired signal f and diminish noise effects $d^{\ell}$
via stochastic cancellation. This technique requires one to presume random
noise and that each acquisition is taken on the same grid points. Our
methodology is more practical since we do not impose any noise structure and
do not require exact repetition in sampling. However, in the case of random
noise, we expect signal averaging to improve the error bound (22) due to the
added cancellation of noise terms. We do not investigate these claims here but
note that such a result can be obtained by our analysis, thereby allowing
signal averaging on distinct grids.

##  5 Discussion

This section elaborates on several aspects of the results. Section 5.1
discusses relevant work in the literature. Section 5.2 provides examples of
distributions that satisfy our deviation model and intuition of its meaning.
Section 5.3 explores the $\gamma$ parameter with examples of transformations
$\Psi$ that produce a satisfiable sampling complexity.

###  5.1 Related Work

Several studies in the compressive sensing literature are similar to our
results in Section 3 ([35, 34]). Our main proof under the sparse signal model
adopts the approach of [49, 35, 34]. However, we derive recovery guarantees
for non-orthonormal systems (when $\theta\neq 0$) while focusing the scope of
the paper within the context of classical sampling theory (introducing error
according to bandlimited approximation). The work in [34] considers sampling
of sparse trigonometric polynomials and overlaps with our application in the
case $\Psi=\mathcal{F}$. Our results generalize this work (in the basis
pursuit case) to allow for other signal models and sparsifying transforms.
Furthermore, [34] assumes that the samples are chosen uniformly at random from
a continuous interval or a discrete set of $N$ equispaced points. In contrast,
our results pertain to general deviations from an equispaced grid with average
sampling density $\sim s\operatorname{polylog}(N)$ and allow for many other
distributions of the perturbations.

###  5.2 Deviation Model

In this section, we present several examples of distributions that are
suitable for our results in Section 3. Notice that our deviation model
utilizes the characteristic function of a given distribution, evaluated at a
finite set of points. This allows one to easily consider many distributions
for our purpose by consulting the relevant and exhaustive literature of
characteristic functions (see for example [21]).

  * •

1) Uniform continuous: $\mathcal{D}=\mathcal{U}[-\frac{1}{2m},\frac{1}{2m}]$
gives $\theta=0$. To generalize this example, we may take
$\mathcal{D}=\mathcal{U}[\mu-\frac{p}{2m},\mu+\frac{p}{2m}]$, for any
$\mu\in\mathbb{R}$ and $p\in\mathbb{N}/\\{0\\}$ to obtain $\theta=0$. Notice
that with $p=m$, we obtain i.i.d. samples chosen uniformly from the whole
interval $\Omega$ (as in [34]).

  * •

2) Uniform discrete:
$\mathcal{D}=\mathcal{U}\\{-\frac{1}{2m}+\frac{k}{m\bar{n}}\\}_{k=0}^{\bar{n}-1}$
with $\bar{n}:=\lceil\frac{2(N-1)}{m}\rceil+1$ gives $\theta=0$. To generalize
this example, we may take
$\mathcal{D}=\mathcal{U}\\{\mu-\frac{p}{2m}+\frac{pk}{m\bar{n}}\\}_{k=0}^{\bar{n}-1}$,
for any $\mu\in\mathbb{R}$ and $p\in\mathbb{N}/\\{0\\}$ and integer
$\bar{n}>\frac{2(N-1)p}{m}$. Here $\theta=0$.

  * •

3) Normal: $\mathcal{D}=\mathcal{N}(\mu,\bar{\sigma}^{2})$, for any mean
$\mu\in\mathbb{R}$ and variance $\bar{\sigma}^{2}>0$. Here
$\theta=\frac{2N}{m}e^{-2(\bar{\sigma}\pi m)^{2}}$. In particular, for fixed
$\bar{\sigma}$, $m$ may be chosen large enough to satisfy the conditions of
Theorem 5 and vice versa.

  * •

4) Laplace: $\mathcal{D}=\mathcal{L}(\mu,b)$, for any location
$\mu\in\mathbb{R}$ and scale $b>0$ gives $\theta=\frac{2N}{m(1+(2\pi
bm)^{2})}$.

  * •

5) Exponential: $\mathcal{D}=\mbox{Exp}(\lambda)$, for any rate $\lambda>0$
gives $\theta=\frac{2N}{m\sqrt{1+4\pi^{2}m^{2}\lambda^{-2}}}$.

In particular, notice that examples 1) and 2) include cases of jittered
sampling [9, 63, 12, 62]. Indeed, with $p=1$ these examples partition $\Omega$
into $m$ regions of equal size and these distributions will choose a point
randomly from each region (in a continuous or discrete sense). The jittered
sampling list can be expanded by considering other distributions to generate
samples within each region (e.g., Poisson distribution). In general we will
have $\theta>0$, which implies deteriorated quality of reconstruction
according to Theorem 2. Arguably, our deviation model introduces a notion of
_optimal_ jitter when the chosen distribution achieves $\theta=0$, ideal in
our results. This observation may be of interest in the active literature of
jittered sampling techniques [62].

Intuitively, $\theta$ is measuring how biased a given distribution is in
generating deviations. If $\delta\sim\mathcal{D}$,
$\lvert\mathbb{E}\textbf{e}(jm\delta)\rvert\approx 0$ means that the
distribution is nearly centered and impartial. On the other hand,
$\lvert\mathbb{E}\textbf{e}(jm\delta)\rvert\approx 1$ gives the opposite
interpretation where the deviations will be generated favoring a certain
direction in an almost deterministic sense. Our result is not applicable to
such biased distributions, since in Theorem 2 as
$\theta\rightarrow\frac{1}{\sqrt{2}}$ we have $\tau\rightarrow 0$ and
$C_{2},C_{3}\rightarrow\infty$ (see Section 8).

###  5.3 Signal Model

In this section we discuss the DFT-incoherence parameter $\gamma$, introduced
in Section 2.1 as

$\gamma=\max_{\ell\in[n]}\sum_{k=1}^{N}\lvert\langle\mathcal{F}_{*k},\Psi_{*\ell}\rangle\rvert,$

where we now let $\Psi\in\mathbb{C}^{N\times n}$ be a full column rank matrix.
The parameter $\gamma$ a uniform upper bound on the $1$-norm of the discrete
Fourier coefficients of the columns of $\Psi$. Since the decay of the Fourier
coefficients of a function is related to its smoothness, intuitively $\gamma$
can be seen as a measure of the smoothness of the columns of $\Psi$.
Implicitly, this also measures the smoothness of f since its uniform
discretization admits a representation via this transformation $f=\Psi g$.
Therefore, the role of $\gamma$ on the sampling complexity is clear, smaller
$\gamma$ implies that our signal of interest is smooth and therefore requires
less samples. This observation is intuitive, since non-smooth functions will
require additional samples to capture discontinuities in accordance with Gibbs
phenomenon. This argument is validated numerically in Section 6.1, where we
compare reconstruction via an infinitely differentiable ensemble (FFT) and a
discontinuous wavelet (Daubechies 2).

We now consider several common choices for $\Psi$ and discuss the respective
$\gamma$ parameter:

  * •

1) $\Psi=\mathcal{F}$ (the DFT), then $\gamma=1$ which is optimal. However,
most appropriate and common is the choice $\Psi=\mathcal{F}^{*}$ which can be
shown to exhibit $\gamma\sim\mathcal{O}(1)$ by a simple calculation.

  * •

2) When $\Psi=\mathcal{H}^{*}$ is the inverse 1D Haar wavelet transform, we
have $\gamma\sim\mathcal{O}(\log(N))$. In [24] it is shown that the inner
products between rows of $\mathcal{F}$ and rows of $\mathcal{H}$ decay
according to an inverse power law of the frequency (see Lemma 1 therein). A
similar proof shows that
$\lvert\langle\mathcal{F}_{*k},\mathcal{H}^{*}_{*\ell}\rangle\rvert\sim\frac{1}{\lvert
k\rvert}$, which gives the desired upper bound for $\gamma$ via an integral
comparison. Notice that these basis vectors have jump discontinuities, and yet
we still obtain an acceptable DFT-incoherence parameter for nonuniform
undersampling.

  * •

3) $\Psi=\mathcal{I}$ (the $N\times N$ identity) gives $\gamma=\sqrt{N}$. This
is the worst case scenario for normalized transforms since

$\max_{x\in
S^{N-1}}\sum_{k=1}^{N}\lvert\langle\mathcal{F}_{*k},x\rangle\rvert=\max_{x\in
S^{N-1}}\sum_{k=1}^{N}\lvert\langle\mathcal{F}_{*k},\mathcal{F}^{*}x\rangle\rvert=\max_{x\in
S^{N-1}}\sum_{k=1}^{N}\lvert x_{k}\rvert\leq\sqrt{N}.$

In general, our smooth signals of interest are not fit for this sparsity
model.

  * •

4) Let $p\geq 1$ be an integer. Matrices $\Psi$ whose columns are uniform
discretizations of $p$-differentiable functions, with $p-1$ periodic and
continuous derivatives and $p$-th derivative that is piecewise continuous. In
this case $\gamma\sim\mathcal{O}(\log(N))$ if $p=1$ and
$\gamma\sim\mathcal{O}(1)$ if $p\geq 2$. For sake of brevity we do not provide
this calculation, but refer the reader to Section 2.8 in [59] for an informal
argument.

Example 4) in particular is informative due to its generality and ability to
somewhat formalize the intuition behind $\gamma$ previously discussed. This
example implies the applicability of our result to a general class of smooth
functions that agree nicely with our signal model defined in Section 2.1
(functions in $A(\Omega)$).

##  6 Numerical Experiments

In this section we present numerical experiments to elaborate several aspects
of our methodology and results. Specifically, we consider the effects of the
DFT-incoherence and $\theta$ parameter in Section 6.1 and 6.2 respectively.
Section 6.3 investigates the noise attenuation properties of nonuniform
samples. We first introduce several terms and models to describe the setup of
the experiments. Throughout we let $N=2015$ be the uniformly discretized
signal size.

Program (1) is solved by SPGL1 [17], a Pareto curve approach that uses duality
theory to solve the basis pursuit problem via a sequence of numerically
cheaper LASSO subproblems. Each basis pursuit problem is solved by limiting
the number of SPGl1 iterations to 200.

We implement the Dirichlet kernel using (6) directly to construct a Spot
operator [18], a matrix free linear operator toolbox that is compatible with
SPGL1. We warn the reader that in this section we have not dedicated much
effort to optimize the numerical complexity of the interpolation kernel. For a
faster implementation, we recommend instead applying the DFT/NDFT
representation $\mathcal{S}=\mathcal{N}\mathcal{F}^{*}$ (see Section 2.3)
using NFFT 3 software from [38] or its parallel counterpart [52].

Given output $f^{\sharp}=\Psi g^{\sharp}$ of any of our programs with true
solution $f$, we consider the relative reconstruction error norm as a measure
of output quality, given as

$\mbox{Relative Reconstruction
Error}=\frac{\|f^{\sharp}-f\|_{2}}{\|f\|_{2}}.\\\ $

_Grid perturbations:_ To construct the nonuniform grid $\tilde{\tau}$, we
introduce an irregularity parameter $\rho\geq 0$. We define our perturbations
by sampling from a uniform distribution, so that each $\Delta_{k}$ is drawn
uniformly at random from $[-\frac{\rho}{m},\frac{\rho}{m}]$ for all $k\in[m]$
independently. $\tilde{\tau}$ is generated independently for each signal
reconstruction experiment.

_Complex exponential signal model:_ We consider bandlimited complex
exponentials with random harmonic frequencies. With bandwidth
$\omega=\frac{N-1}{2}=1007$, and sparsity level $s=50$ we generate
$\vec{\omega}\in\mathbb{Z}^{s}$ by choosing $s$ frequencies uniformly at
random from $\\{-\omega,-\omega+1,\cdots,\omega\\}$ and let

$\textbf{f}(x)=\sum_{k=1}^{s}\textbf{e}\left(\vec{\omega}_{k}x\right).$

We use the DFT as a sparsifying transform $\Psi=\mathcal{F}$ so that
$g=\Psi^{*}f=\Psi^{-1}f$ is a $50$-sparse vector. This transform is
implemented as a Spot operator which utilizes MATLAB’s fft function. The
frequency vector, $\vec{\omega}$, is generated randomly for each independent
set of experiments. Note that in this case we have optimal DFT-incoherence
parameter $\gamma=1$ (see Section 5.3).

_Gaussian signal model:_ We consider a non-bandlimited signal consisting of
sums of Gaussian functions. This signal model is defined as

$\textbf{f}(x)=-e^{-100x^{2}}+e^{-100(x-.104)^{2}}-e^{-100(x+.217)^{2}}.$

For this dataset, we use the Daubechies 2 wavelet as a sparsifying transform
$\Psi$. This operator is implemented as a Spot operator which utilizes the
Rice Wavelet Toolbox [66]. This provides $g=\Psi^{*}f=\Psi^{-1}f$ that can be
well approximated by a $50$-sparse vector. In other words, all entries of $g$
are non-zero but $\epsilon_{50}(g)<.088\approx\frac{\|f\|_{2}}{250}$. In this
case we have $\gamma\approx 40.78$, which was computed numerically.

###  6.1 Effect of DFT-incoherence

This section is dedicated to exploring the effect of the DFT-incoherence
parameter in signal reconstruction. We consider the complex exponential and
Gaussian signal models described above, notice that both signals are
effectively $50$-sparse. Recall that in the complex exponential model we have
$\Psi=\mathcal{F}$ (the DFT) with optimal DFT-incoherence parameter
$\gamma=1$. In the Gaussian model $\Psi$ is the Daubechies 2 wavelet with
$\gamma\approx 40.78$.

Here we set $\rho=\frac{1}{2}$ to generate the deviations (so that $\theta=0$)
and vary the average step size of the nonuniform samples. We do so by letting
$m$ vary through the set
$\\{\lfloor\frac{N}{1.5}\rfloor,\lfloor\frac{N}{2}\rfloor,\lfloor\frac{N}{2.5}\rfloor,\cdots,\lfloor\frac{N}{9.5}\rfloor\\}$.
For each fixed value of $m$, the relative reconstruction error is obtained by
averaging the result of 10 independent reconstruction experiments. The results
are shown in Figure 1, where we plot the average step size vs average relative
reconstruction error.

Figure 1: Plot of average relative reconstruction error vs average step size
for both signal models. In the complex exponential model ($\Psi=$ DFT) we have
$\gamma=1$ and in the Gaussian signal model we have $\gamma\approx 40.78$
(Daubechies 2 wavelet). Notice that the complex exponential model allows for
reconstruction from larger step sizes in comparison to the Gaussian signal
model.

These experiments demonstrate the negative effect of larger DFT-incoherence
parameters in signal reconstruction. Indeed, in Figure 1 we see that the
complex exponential model with $\gamma=1$ allows for accurate reconstruction
from larger step sizes. This is to be expected from Section 3, where the
results imply that the Daubechies 2 wavelet will require more samples with
diminished probability of successful reconstruction according to its parameter
$\gamma\approx 40.78$.

###  6.2 Effect of Deviation Model Parameter $\theta$

In this section we generate the deviations in such a way that vary the
parameter $\theta$, in order to explore its effect on signal reconstruction.
We only consider the complex exponential signal model for this purpose and fix
$m=\lfloor\frac{N}{7}\rfloor=287$.

We vary $\theta$ by generating deviations with irregularity parameter $\rho$
varying over $\\{.06,.07,.08,\cdots,.5\\}$. For each fixed $\rho$ value we
compute the average relative reconstruction error of 50 independent
experiments. Notice that for each $k\in[m]$ and any $j$

$\mathbb{E}\textbf{e}\left(jm\Delta_{k}\right)=\frac{\sin\left(2\pi
j\rho\right)}{2\pi j\rho}.$

Given $\rho$, we use this observation to compute the respective $\theta$ value
by considering all $0<\lvert j\rvert\leq\lfloor\frac{2(N-1)}{m}\rfloor=14$.
Figure 2 plots $\theta$ vs average relative reconstruction error.

Figure 2: Plot of average relative reconstruction error vs corresponding
$\theta$ parameter. The plot emphasizes via red outlines the $\theta$ values
that satisfy the conditions of Theorem 2 (i.e., $\theta<\frac{1}{\sqrt{2}}$).
Although our results only hold for three $\theta$ values ($0,.29,.58$), the
experiments demonstrate that accurate recovery is possible otherwise.

In the plot, our main result (Theorem 2) is only strictly applicable in three
cases (outlined in red, $\theta=0,.29,.58$). However, Figure 2 demonstrates
that decent signal reconstruction may be achieved when the condition
$\theta<\frac{1}{\sqrt{2}}$ does not hold. Therefore, the applicability of the
methodology goes beyond the restrictions of the theorem.

###  6.3 Noise Attenuation

This section explores the robustness of the methodology when presented with
measurement noise, in both the undersampled and oversampled case relative to
the target bandwidth $\frac{N-1}{2}$ (Sections 3 and 4 respectively). We only
solve the basis pursuit problem (1), and avoid the least squares problem (19)
for brevity. However, we note that both programs produce similar results and
conclusions in the oversampled case. We only consider the bandlimited complex
exponential signal model for this purpose. We generate additive random noise
$d\in\mathbb{R}^{m}$ from a uniform distribution. Each entry of
$d\in\mathbb{R}^{m}$ is i.i.d. from $[-\frac{\chi}{1000},\frac{\chi}{1000}]$
where $\chi=\frac{1}{\sqrt{m}}\|f\|_{1}$, chosen to maintain $\|d\|_{2}$
relatively constant as $m$ varies.

We set $\rho=\frac{1}{2}$ to generate the deviations (so that $\theta=0$) and
vary the average step size of the nonuniform samples. We do so by letting $m$
vary through the set
$\\{\lfloor\frac{N}{.5}\rfloor,\lfloor\frac{N}{.75}\rfloor,N,\cdots,\lfloor\frac{N}{6.75}\rfloor,\lfloor\frac{N}{7}\rfloor\\}$,
notice that only the first two cases correspond to oversampling. For each
fixed value of $m$, the relative reconstruction error is obtained by averaging
the result of 50 independent experiments. The results are shown in Figure 3,
where we plot the average step size vs average relative reconstruction error
and average relative noise level $\|d\|_{2}/\|f\|_{2}$.

Figure 3: Plot of average relative reconstruction error
($\|f-f^{\sharp}\|_{2}/\|f\|_{2}$) vs average step size (blue curve) and
average input relative measurement error ($\|d\|_{2}/\|f\|_{2}$) vs average
step size (red curve). Notice that for the first ten step size values noise
attenuation is achieved, i.e., the reconstruction error is lower than the
input noise level.

The first two cases
($m=\lfloor\frac{N}{.5}\rfloor,\lfloor\frac{N}{.75}\rfloor$) correspond to
oversampling and illustrate the results from Section 4, where attenuation of
the input noise level is achieved. Surprisingly, these experiments demonstrate
that nonuniform undersampling also allows for denoising. This is seen in
Figure 3, where the values
$m=\lfloor\frac{N}{1.25}\rfloor,\lfloor\frac{N}{1.5}\rfloor,\cdots,\lfloor\frac{N}{2.75}\rfloor$
output an average relative reconstruction error smaller than the input
measurement noise level. Thus, when nonuniform samples are not heavily
undersampled, the negative effects of random noise can be reduced.

##  7 Conclusions

This paper provides a concrete framework to study the benefits of random
nonuniform samples for signal acquisition (in comparison to equispaced
sampling), with explicit statements that are informative for practitioners.
Related observations are extensive but largely empirical in the sampling
theory literature. Therefore, this work supplies novel theoretical insights on
this widely discussed phenomenon. In the context of compressive sensing, we
extend the applicability of this novel acquisition paradigm by demonstrating
how it naturally intersects with standard sampling techniques. We hope that
these observations will prompt a broader usage of compressive sensing in real
world applications that rely of classical sampling theory.

There are several avenues for future research. First, the overall methodology
requires the practitioner to know the nonuniform sampling locations
$\tilde{\tau}$ accurately. While this is typical for signal reconstruction
techniques that involve non-equispaced samples, it would be of practical
interest to extend the methodology is such a way that allows for robustness to
inaccurate sampling locations and even self-calibration. Furthermore, the
implementation requires an estimate of the best bandlimited approximation
error in the Wiener algebra norm $\sum_{\lvert\ell\rvert>(N-1)/2}\lvert
c_{\ell}\rvert$ which is generally unknown a priori. Corollary 3 offers a
specific example where this quantity is not needed, but a general result that
removes this prerequisite would be most practical (e.g., by establishing a
variant of the quotient property [73])

As mentioned in Section 6, this work has not dedicated much effort to a
numerically efficient implementation of the Dirichlet kernel $\mathcal{S}$.
This is crucial for large-scale applications, where a direct implementation of
the Dirichlet kernel via its Fourier or Dirichlet representation (see [58])
may be too inefficient for practical purposes. As future work, it would be
useful to consider other interpolation kernels with greater numerical
efficiency (e.g., a low order Lagrange interpolation operator).

## Acknowledgments

This work was in part financially supported by the Natural Sciences and
Engineering Research Council of Canada (NSERC) Collaborative Research and
Development Grant DNOISE II (375142-08). This research was carried out as part
of the SINBAD II project with support from the following organizations: BG
Group, BGP, CGG, Chevron, ConocoPhillips, DownUnder GeoSolutions, Hess
Corporation, Petrobras, PGS, Sub Salt Solutions, WesternGeco, and Woodside. Ö.
Yılmaz also acknowledges an NSERC Discovery Grant (22R82411) and an NSERC
Accelerator Award (22R68054).

##  8 Proofs

We now provide proofs to all of our claims. In Section 8.1 we prove Theorem 2
via a more general result. Theorem 4 is proven in Section 8.2. Section 8.3
establishes the Dirichlet kernel error bounds from Theorem 1 and Corollary 1.

###  8.1 Proof of Theorem 2

In this section, we will prove a more general result than Theorem 2 assuming
that $\Psi$ is a full column-rank matrix. Theorem 2 will follow from Theorem 5
by taking $\alpha=\beta=1$.

######  Theorem 5.

Let $n\leq N$ and $\Psi\in\mathbb{C}^{N\times n}$ be a full column rank matrix
with DFT-incoherence parameter $\gamma$ and extreme singular values

$\sigma_{n}(\Psi):=\alpha>\sqrt{1-\frac{1}{\sqrt{2}}},\ \ \ \
\sigma_{1}(\Psi):=\beta\leq\sqrt{1+\frac{1}{\sqrt{2}}}.$

Let the entries of $\Delta$ be i.i.d. from any distribution satisfying our
deviation model with

$\theta<\max\bigg{\\{}\frac{1}{\beta^{2}}\left(1+\frac{1}{\sqrt{2}}\right)-1,\
1+\frac{1}{\alpha^{2}}\left(\frac{1}{\sqrt{2}}-1\right)\bigg{\\}}.$

Define

$g^{\sharp}=\operatorname{argmin}_{h\in\mathbb{C}^{n}}\|h\|_{1}\ \
\mbox{s.t.}\ \ \|\mathcal{S}\Psi
h-b\|_{2}\leq\eta+2\sqrt{m}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ (23)

and let
$\tau:=\frac{1}{4\sqrt{2}}-\frac{\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}}{4}>0$.

If

$\sqrt{m}\geq\tilde{c}_{1}\gamma\sqrt{s}\log^{2}(n)\frac{\sqrt{\beta^{2}+\beta^{2}\theta+\tau}}{\tau}$
(24)

where $\tilde{c}_{1}$ is an absolute constant, then for
$\tilde{c}_{2},\tilde{c}_{3}>0$ depending only on $\theta,\alpha$ and $\beta$,

$\|f-\Psi
g^{\sharp}\|_{2}\leq\frac{\beta\tilde{c}_{2}\epsilon_{s}(g)}{\sqrt{s}}+\frac{\beta\tilde{c}_{3}\sqrt{N}\eta}{\sqrt{m}}+2\beta\tilde{c}_{3}\sqrt{N}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$

with probability exceeding

$1-\exp\left(-\frac{m\tau}{8s\gamma^{2}}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+\beta^{2}(1+2\theta)}\right)\right)\right).$

This theorem generalizes Theorem 2 to more general transformations $\Psi$ for
sparse representation. This is more practical since the columns of $\Psi$ need
not be orthogonal, instead linear independence suffices (with knowledge of the
singular values $\alpha,\beta$). In particular notice that (24) requires
$m\sim s\log^{4}(n)$ as opposed to $m\sim s\log^{4}(N)$ in (12). Since $n\leq
N$, this general result allows for a potential reduction in sample complexity
if the practitioner may construct $\Psi$ in such an efficient manner while
still allowing a sparse and accurate representation of $f$.

For the proof we will use the restricted isometry property (RIP), a common
analytical tool used to establish sparse vector recovery results via basis
pursuit. We will obtain our main result by establishing this property for our
interpolation kernel $\mathcal{S}$ with respect to the transform $\Psi$.
Specifically, we wish to show that $\mathcal{S}\Psi$ satisfies the following
[73, 77]:

_Definition:_ Suppose $\mathcal{A}\in\mathbb{C}^{m\times n}$ is a measurement
matrix and $1\leq s\leq n$ is an integer. The restricted isometry constant
(RIC) of order $s$ of $\mathcal{A}$ is the smallest number $\delta_{s}$ such
that for all $s$-sparse vectors $v\in\mathbb{C}^{n}$ (i.e., $\|v\|_{0}\leq
s$),

$(1-\delta_{s})\|v\|_{2}^{2}\leq\|\mathcal{A}v\|_{2}^{2}\leq(1+\delta_{s})\|v\|_{2}^{2}.$
(25)

Informally, the RIP is said to hold if the RIC is small enough for sparse
vector recovery. Many results exist in the literature that give such
conditions on the RIC. To the best of the author’s knowledge, [77] provides
the best result in terms of the RIC, which was extended to the complex valued
case in [72]:

######  Theorem 6 (Theorem 6 in [72]).

Let $\mathcal{A}\in\mathbb{C}^{m\times n}$, $g\in\mathbb{C}^{n}$ and

$g^{\sharp}:=\operatorname{argmin}_{h\in\mathbb{C}^{n}}\|h\|_{1}\ \
\mbox{s.t.}\ \ \|\mathcal{A}h-b\|_{2}\leq\eta$

where $b=\mathcal{A}(g)+d$ and $\|d\|_{2}\leq\eta$. If the RIC $\delta_{2s}$
of $\mathcal{A}$ satisfies

$\delta_{2s}<\frac{1}{\sqrt{2}},$

then

$\|g-g^{\sharp}\|_{2}\leq\frac{\tilde{c}_{2}\epsilon_{s}(g)}{\sqrt{s}}+\tilde{c}_{3}\eta,$

where $\tilde{c}_{2},\tilde{c}_{3}$ only depend on $\delta_{2s}$.

In our context, we will see that $\delta_{2s}$ only depends on $\theta,\alpha$
and $\beta$ (and therefore $\tilde{c}_{2}$ and $\tilde{c}_{3}$ as well). This
result provides recovery error bounds proportional to the noise level and
sparse model mismatch. We will thus compute the RIC of $\mathcal{S}\Psi$ and
then apply Theorem 6 to obtain our main result.

####  8.1.1 Restricted Isometry Property of $\mathcal{S}\Psi$

We will bound the RIC of
$\mathcal{A}:=\frac{\sqrt{N}}{\sqrt{m}}\mathcal{S}\Psi\in\mathbb{C}^{m\times
n}$. For an index set $T\subset[n]$, let
$\mathcal{A}_{T}\in\mathbb{C}^{m\times\lvert T\rvert}$ be the sub-matrix of
$\mathcal{A}$ that results by only considering the columns of $\mathcal{A}$
specified by $T$. We consider the term

$X_{m}:=\sup_{\lvert T\rvert\leq
s}\|\mathcal{A}_{T}^{*}\mathcal{A}_{T}-\mathbb{E}\mathcal{A}_{T}^{*}\mathcal{A}_{T}\|.$

Notice that in contrast to many works in the literature that adopt the same
approach ([34, 49]), we have
$\mathbb{E}\mathcal{A}_{T}^{*}\mathcal{A}_{T}\neq\mathcal{I}_{s\times s}$ in
general (equality only holds in the case $\theta=0$ and $\alpha=\beta=1$). We
must therefore deal with the term
$\mathbb{E}\mathcal{A}_{T}^{*}\mathcal{A}_{T}$ explicitly and modify our
approach respectively.

In order to bound the RIC, our goal is to show

$X_{m}\leq\tau.$ (26)

for some $\tau\geq 0$. This will in turn show that the order $s$ RIC of
$\mathcal{A}$ is bounded as

$\delta_{s}\leq\tau+\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\},$

where the additional term
$\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}$ is due
to our deviation model and assumptions on $\Psi$. This is shown in the
following argument:

The variational characterization of the spectral norm for a Hermitian matrix
gives

$\displaystyle X_{m}$ $\displaystyle=\sup_{\lvert T\rvert\leq s}\max_{w\in
S^{s-1}}\lvert\langle\left(\mathcal{A}_{T}^{*}\mathcal{A}_{T}-\mathbb{E}\mathcal{A}_{T}^{*}\mathcal{A}_{T}\right)w,w\rangle\rvert$
$\displaystyle=\sup_{w\in S^{n-1},\|w\|_{0}\leq
s}\lvert\langle\left(\mathcal{A}^{*}\mathcal{A}-\mathbb{E}\mathcal{A}^{*}\mathcal{A}\right)w,w\rangle\rvert$
$\displaystyle=\sup_{w\in S^{n-1},\|w\|_{0}\leq
s}\lvert\|\mathcal{A}w\|_{2}^{2}-\mathbb{E}\|\mathcal{A}w\|_{2}^{2}\rvert.$

Therefore, $X_{m}\leq\tau$ is equivalent to

$\mathbb{E}\|\mathcal{A}w\|_{2}^{2}-\tau\|w\|_{2}^{2}\leq\|\mathcal{A}w\|_{2}^{2}\leq\mathbb{E}\|\mathcal{A}w\|_{2}^{2}+\tau\|w\|_{2}^{2},$

and to establish our claim it suffices to show that

$\mathbb{E}\|\mathcal{A}w\|_{2}^{2}\leq\beta^{2}(1+\theta)\|w\|_{2}^{2}\ \ \ \
\mbox{and}\ \ \ \
\mathbb{E}\|\mathcal{A}w\|_{2}^{2}\geq\alpha^{2}(1-\theta)\|w\|_{2}^{2}.$ (27)

To this end, let $w\in\mathbb{C}^{n}$ and normalize
$\tilde{\mathcal{N}}:=\sqrt{N}\mathcal{N}$ so that for $k\in[m],\ell\in[N]$

$\tilde{\mathcal{N}}_{k\ell}:=\textbf{e}(-\tilde{t}_{k}(\ell-\tilde{n}-1)).$

Throughout, let $\tilde{\Delta}\in\mathbb{R}$ be an independent copy of the
entries of $\Delta\in\mathbb{R}^{m}$. Then with $v:=\mathcal{F}^{*}\Psi w$,

$\displaystyle\mathbb{E}\|\mathcal{A}w\|_{2}^{2}=\frac{1}{m}\mathbb{E}\|\tilde{\mathcal{N}}\mathcal{F}^{*}\Psi
w\|_{2}^{2}:=\frac{1}{m}\mathbb{E}\|\tilde{\mathcal{N}}v\|_{2}^{2}$
$\displaystyle=\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\lvert\langle\tilde{\mathcal{N}}_{k*},v\rangle\rvert^{2}=\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\Biggl{\lvert}\sum_{\ell=1}^{N}\textbf{e}(\tilde{t}_{k}(\ell-\tilde{n}-1))v_{\ell}\Biggr{\rvert}^{2}$
$\displaystyle=\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\left(\sum_{\ell=1}^{N}\sum_{\tilde{\ell}=1}^{N}\textbf{e}(\tilde{t}_{k}(\ell-\tilde{\ell}))v_{\ell}\bar{v}_{\tilde{\ell}}\right)$
$\displaystyle=\sum_{\ell=1}^{N}\sum_{\tilde{\ell}=1}^{N}v_{\ell}\bar{v}_{\tilde{\ell}}\left(\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\textbf{e}(\tilde{t}_{k}(\ell-\tilde{\ell}))\right)$
$\displaystyle=\sum_{\ell=1}^{N}\lvert
v_{\ell}\rvert^{2}+\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\sum_{\ell-\tilde{\ell}=jm}v_{\ell}\bar{v}_{\tilde{\ell}}\mathbb{E}\textbf{e}(jm(\tilde{\Delta}-1/2))$
$\displaystyle+\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\sum_{\ell-\tilde{\ell}=-jm}v_{\ell}\bar{v}_{\tilde{\ell}}\mathbb{E}\textbf{e}(-jm(\tilde{\Delta}-1/2)).$

The last equality can be obtained as follows,

$\displaystyle\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\textbf{e}(\tilde{t}_{k}(\ell-\tilde{\ell}))=\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\textbf{e}\left(\left(\frac{k-1}{m}-\frac{1}{2}+\Delta_{k}\right)(\ell-\tilde{\ell})\right)$
$\displaystyle=\frac{1}{m}\sum_{k=1}^{m}\textbf{e}\left(\left(\frac{k-1}{m}-\frac{1}{2}\right)(\ell-\tilde{\ell})\right)\mathbb{E}\textbf{e}(\Delta_{k}(\ell-\tilde{\ell}))$
$\displaystyle=\frac{1}{m}\sum_{k=1}^{m}\textbf{e}\left(\left(\frac{k-1}{m}-\frac{1}{2}\right)(\ell-\tilde{\ell})\right)\mathbb{E}\textbf{e}(\tilde{\Delta}(\ell-\tilde{\ell}))$
$\displaystyle=\mathbb{E}\textbf{e}\left((\tilde{\Delta}-1/2)(\ell-\tilde{\ell})\right)\sum_{k=1}^{m}\frac{1}{m}\textbf{e}\left(\frac{k-1}{m}(\ell-\tilde{\ell})\right)$
$\displaystyle=\left\\{\begin{array}[]{ll}1&\mbox{if }\ell=\tilde{\ell}\\\
\mathbb{E}\textbf{e}\left(jm(\tilde{\Delta}-1/2)\right)&\mbox{if
}\ell-\tilde{\ell}=jm,j\in\mathbb{Z}/\\{0\\}\\\
0&\mbox{otherwise.}\end{array}\right.$

The third equality uses the fact that
$\mathbb{E}\textbf{e}(\Delta_{k}(\ell-\tilde{\ell}))=\mathbb{E}\textbf{e}(\tilde{\Delta}(\ell-\tilde{\ell}))$
for all $k\in[m]$ in order to properly factor out this constant from the sum
in the fourth equality. The last equality is due to the geometric series
formula.

Returning to our original calculation, we bound the last term using our
deviation model assumptions

$\displaystyle\Biggl{\lvert}\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\sum_{\ell-\tilde{\ell}=-jm}v_{\ell}\bar{v}_{\tilde{\ell}}\mathbb{E}\textbf{e}(-jm(\tilde{\Delta}-1/2))\Biggr{\rvert}$
$\displaystyle=\Biggl{\lvert}\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\sum_{\ell\in
Q_{j}}v_{\ell}\bar{v}_{\ell+jm}\mathbb{E}\textbf{e}(-jm(\tilde{\Delta}-1/2))\Biggr{\rvert}$
$\displaystyle\leq\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\sum_{\ell\in Q_{j}}\lvert
v_{\ell}\rvert\lvert
v_{\ell+jm}\rvert\lvert\mathbb{E}\textbf{e}(-jm(\tilde{\Delta}-1/2))\rvert$
$\displaystyle\leq\frac{\theta
m}{2N}\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\sum_{\ell\in Q_{j}}\lvert
v_{\ell}\rvert\lvert v_{\ell+jm}\rvert\leq\frac{\theta
m}{2N}\sum_{j=1}^{\lfloor(N-1)/m\rfloor}\lvert\langle v,v\rangle\rvert$
$\displaystyle=\frac{\theta
m\|v\|_{2}^{2}}{2N}\sum_{j=1}^{\lfloor(N-1)/m\rfloor}1=\frac{\theta
m\|v\|_{2}^{2}}{2N}\lfloor\frac{N-1}{m}\rfloor\leq\frac{\theta\|v\|_{2}^{2}}{2}.$

$Q_{j}\subset[N]$ is the index set of allowed $\ell$ indices according to $j$,
i.e., that satisfy $\ell\in[N]$ and $\ell+jm\in[N]$. The remaining sum can be
bounded similarly. Combine these inequalities with the singular values of
$\Psi$ to obtain

$\mathbb{E}\|\mathcal{A}w\|_{2}^{2}\leq\|v\|_{2}^{2}+\frac{2\theta\|v\|_{2}^{2}}{2}:=\|\Psi
w\|_{2}^{2}\left(1+\theta\right)\leq\beta^{2}\|w\|_{2}^{2}\left(1+\theta\right),$

and

$\mathbb{E}\|\mathcal{A}w\|_{2}^{2}\geq\alpha^{2}\|w\|_{2}^{2}\left(1-\theta\right).$

We will apply this inequality and similar orthogonality properties in what
follows (and in Section 8.2), and ask the reader to keep this in mind.

To establish (26), we will apply a concentration inequality (Theorem 5.2 in
[34], established in [75]). The following lemma will be useful for this
purpose.

######  Lemma 1.

Let $v\in\mathbb{C}^{n}$ be an $s$-sparse vector and define
$\tilde{\mathcal{N}}\in\mathbb{C}^{m\times N}$ as above with the entries of
$\Delta$ i.i.d. from any distribution satisfying our deviation model for some
$\theta\geq 0$ (see Section 2.2). Let $\gamma$ be the DFT-incoherence
parameter of $\Psi\in\mathbb{C}^{N\times n}$ (see Section 2.1). Then for all
$k\in[m]$

$\displaystyle\lvert\langle\mathcal{A}_{k*},v\rangle\rvert\leq\frac{\gamma\sqrt{s}\|v\|_{2}}{\sqrt{m}},$

and

$\mathbb{E}\sum_{k=1}^{m}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{4}\leq\frac{s\beta^{2}\gamma^{2}(1+2\theta)\|v\|_{2}^{4}}{m}$

The proof of this lemma can be found in Section 8.1.3.

We are now in a position to bound the RIC. The proof will follow the arguments
in [34, 49]. In particular, the main difference here is that the rows of
$\mathcal{A}$ are not identically distributed and
$\mathbb{E}\mathcal{A}^{*}\mathcal{A}\neq\mathcal{I}_{n\times n}$ in general
(equality only holds in the case $\theta=0$ and $\alpha=\beta=1$). However,
the proof goes through similarly, where we need only to apply the previous
lemma.

######  Theorem 7.

Let the entries of $\Delta$ be i.i.d. from any distribution satisfying our
deviation model for some $\theta\geq 0$ and $\gamma$ be the DFT-incoherence
parameter of $\Psi$. If

$\frac{\sqrt{m}}{\sqrt{\log(m)}}\geq
2c\gamma\sqrt{s}\log(s)\sqrt{\log(n)}\frac{\sqrt{1+\theta+\tau}}{\tau},$

where $c$ is an absolute constant, then the order $s$ restricted isometry
constant of
$\mathcal{A}:=\frac{\sqrt{N}}{\sqrt{m}}\mathcal{S}\Psi\in\mathbb{C}^{m\times
n}$ satisfies

$\delta_{s}\leq
2\tau+\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}$

with probability exceeding

$1-\exp\left(-\frac{m\tau}{4s\gamma^{2}}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+\beta^{2}(1+2\theta)}\right)\right)\right).$

###### Proof:.

In what follows, for $k\in[m]$ let $\mathcal{A}_{T}^{k}\in\mathbb{C}^{\lvert
T\rvert}$ denote the $k$-th row of
$\mathcal{A}_{T}\in\mathbb{C}^{m\times\lvert T\rvert}$. Define

$X_{m}:=\sup_{\lvert T\rvert\leq
s}\|\mathcal{A}_{T}^{*}\mathcal{A}_{T}-\mathbb{E}\mathcal{A}_{T}^{*}\mathcal{A}_{T}\|=\sup_{\lvert
T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\left(\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}-\mathbb{E}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\right)\Biggr{\|}.$

As noted in the beginning of the section, our RIC bound will be established if
we show $X_{m}\leq 2\tau$ with high probability. We begin by bounding
$\mathbb{E}_{\mathcal{D}}X_{m}$, and then show that $X_{m}$ does not deviate
much from its expectation with high probability.

Symmetrizing as Lemma 6.3 in [47], we have

$\mathbb{E}_{\mathcal{D}}X_{m}\leq
2\mathbb{E}_{\mathcal{D}}\mathbb{E}_{\epsilon}\sup_{\lvert T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\epsilon_{k}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\Biggr{\|},$

where the $\epsilon_{k}$ are independent Rademacher random variables. We will
now apply Lemma 3.6 in [49], which requires $m\leq n$ and an upper bound of
$\|\mathcal{A}_{k*}\|_{\infty}$ for all $k\in[m]$ accordingly. For any
$k\in[m]$, we obtain

$\displaystyle\|\mathcal{A}_{k*}\|_{\infty}=\max_{\ell\in[n]}\lvert\mathcal{A}_{k\ell}\rvert=\max_{\ell\in[n]}\frac{1}{\sqrt{m}}\lvert\langle\tilde{\mathcal{N}}_{k*},(\mathcal{F}^{*}\Psi)_{*\ell}\rangle\rvert=\max_{\ell\in[n]}\frac{1}{\sqrt{m}}\lvert\sum_{p=1}^{N}\tilde{\mathcal{N}}_{kp}(\mathcal{F}^{*}\Psi)_{p\ell}\rvert$
$\displaystyle\leq\frac{1}{\sqrt{m}}\max_{\ell\in[n]}\sum_{p=1}^{N}\lvert\tilde{\mathcal{N}}_{kp}\rvert\lvert(\mathcal{F}^{*}\Psi)_{p\ell}\rvert=\frac{1}{\sqrt{m}}\max_{\ell\in[n]}\sum_{p=1}^{N}\lvert(\mathcal{F}^{*}\Psi)_{p\ell}\rvert:=\frac{\gamma}{\sqrt{m}}.$

Applying Lemma 3.6 in [49] with this observation gives that

$2\mathbb{E}_{\epsilon}\sup_{\lvert T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\epsilon_{k}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\Biggr{\|}\leq
Q(s,m,n)\sup_{\lvert T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\Biggr{\|}^{1/2}.$

We have defined

$Q(s,m,n):=\frac{2c\gamma\sqrt{s}\log(s)\sqrt{\log(n)\log(m)}}{\sqrt{m}},$

where $c$ is an absolute constant given in Lemma 3.6 from [49]. We may
therefore bound

$\displaystyle\mathbb{E}_{\mathcal{D}}X_{m}\leq
Q(s,m,n)\mathbb{E}_{\mathcal{D}}\sup_{\lvert T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\Biggr{\|}^{1/2}$
$\displaystyle\leq Q(s,m,n)\left(\mathbb{E}_{\mathcal{D}}\sup_{\lvert
T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\Biggr{\|}\right)^{1/2}$
$\displaystyle\leq Q(s,m,n)\left(\mathbb{E}_{\mathcal{D}}\sup_{\lvert
T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\left(\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}-\mathbb{E}_{\mathcal{D}}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\right)\Biggr{\|}+\sup_{\lvert
T\rvert\leq
s}\Biggl{\|}\mathbb{E}_{\mathcal{D}}\mathcal{A}_{T}^{*}\mathcal{A}_{T}\Biggr{\|}\right)^{1/2}$
$\displaystyle\leq Q(s,m,n)\left(\mathbb{E}_{\mathcal{D}}\sup_{\lvert
T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\left(\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}-\mathbb{E}_{\mathcal{D}}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\right)\Biggr{\|}+\beta^{2}+\beta^{2}\theta\right)^{1/2}$
$\displaystyle:=Q(s,m,n)\left(\mathbb{E}_{\mathcal{D}}X_{m}+\beta^{2}+\beta^{2}\theta\right)^{1/2}.$

The last inequality holds as in the proof of (27), since

$\displaystyle\sup_{\lvert T\rvert\leq
s}\|\mathbb{E}_{\mathcal{D}}\mathcal{A}_{T}^{*}\mathcal{A}_{T}\|\leq\|\mathbb{E}_{\mathcal{D}}\mathcal{A}^{*}\mathcal{A}\|=\sup_{v\in
S^{n-1}}\mathbb{E}_{\mathcal{D}}\langle\mathcal{A}^{*}\mathcal{A}v,v\rangle$
$\displaystyle=\sup_{v\in
S^{n-1}}\mathbb{E}_{\mathcal{D}}\|\mathcal{A}v\|_{2}^{2}\leq\beta^{2}+\beta^{2}\theta.$

In conclusion,

$\frac{\mathbb{E}_{\mathcal{D}}X_{m}}{\left(\mathbb{E}_{\mathcal{D}}X_{m}+\beta^{2}+\beta^{2}\theta\right)^{1/2}}\leq\frac{2c\gamma\sqrt{s}\log(s)\sqrt{\log(n)\log(m)}}{\sqrt{m}}$

and we may achieve

$\mathbb{E}_{\mathcal{D}}X_{m}\leq\tau$

if

$\frac{2c\gamma\sqrt{s}\log(s)\sqrt{\log(n)\log(m)}}{\sqrt{m}}\leq\frac{\tau}{\sqrt{\tau+\beta^{2}+\beta^{2}\theta}}.$
(28)

We now apply a concentration inequality to show that $X_{m}$ is close to its
expected value with high probability. We use the following result (Theorem 5.2
in [34], established in [75]).

######  Theorem 8.

Let $Y_{1},...,Y_{m}$ be a sequence of independent random variables with
values in some Polish space $\mathcal{X}$. Let $\mathcal{F}$ be a countable
collection of real-valued measurable and bounded functions $f$ on
$\mathcal{X}$ with $\|f\|_{\infty}\leq B$ for all $f\in\mathcal{F}$. Let $Z$
be the random variable

$Z=\sup_{f\in F}\sum_{k=1}^{m}f(Y_{k}).$

Assume $\mathbb{E}f(Y_{k})=0$ for all $k\in[m]$ and all $f\in\mathcal{F}$.
Define
$\sigma^{2}=\sup_{f\in\mathcal{F}}\mathbb{E}\sum_{k=1}^{m}f(Y_{\ell})^{2}$.
Then for all $t\geq 0$

$\mathbb{P}\left(Z\geq\mathbb{E}Z+t\right)\leq\exp\left(-\frac{t}{4B}\log\left(1+2\log\left(1+\frac{Bt}{2B\mathbb{E}Z+\sigma^{2}}\right)\right)\right).$

To apply the theorem, let

$\mathcal{D}_{s}:=\\{v\in S^{n-1}:\|v\|_{0}\leq s\\},$

and

$\mathcal{X}:=\\{Z\in\mathbb{C}^{n\times n}:\langle Zv,v\rangle\in\mathbb{R}\
\ \forall v\in\mathcal{D}_{s}\ \ \mbox{and}\ \
\sup_{v\in\mathcal{D}_{s}}\langle Zv,v\rangle\leq\frac{\gamma^{2}s}{m}\\},$

which is a closed subset of $\mathbb{C}^{n\times n}$ (a Polish space), and
therefore itself a Polish space (see for example [25]). For
$v\in\mathcal{D}_{s}$, we define the function
$f_{v}:\mathcal{X}\mapsto\mathbb{R}$ as

$f_{v}(Z):=\langle Zv,v\rangle,$

which is real valued by definition of $\mathcal{X}$. Then notice that

$\displaystyle X_{m}=\sup_{\lvert T\rvert\leq
s}\Biggl{\|}\sum_{k=1}^{m}\left(\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}-\mathbb{E}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\right)\Biggr{\|}$
$\displaystyle=$ $\displaystyle\sup_{\lvert T\rvert\leq s}\sup_{v\in
S^{s-1}}\sum_{k=1}^{m}\langle\left(\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}-\mathbb{E}\mathcal{A}_{T}^{k*}\mathcal{A}_{T}^{k}\right)v,v\rangle$
$\displaystyle=$
$\displaystyle\sup_{v\in\mathcal{D}_{s}}\sum_{k=1}^{m}\langle\left(\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}-\mathbb{E}\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}\right)v,v\rangle=\sup_{v\in\mathcal{D}_{s}}\sum_{k=1}^{m}f_{v}\left(\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}-\mathbb{E}\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}\right)$
$\displaystyle=$
$\displaystyle\sup_{v\in\mathcal{D}_{s}^{*}}\sum_{k=1}^{m}f_{v}\left(\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}-\mathbb{E}\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}\right),$

where $\mathcal{D}_{s}^{*}$ is a dense countable subset of $\mathcal{D}_{s}$.
Furthermore, by the first part of Lemma 1, for all $k\in[m]$ and
$v\in\mathcal{D}_{s}$

$f_{v}(\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}-\mathbb{E}\mathcal{A}_{k*}^{*}\mathcal{A}_{k*})=\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{2}-\mathbb{E}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{2}\leq\frac{\gamma^{2}s}{m},$

so that
$\\{\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}-\mathbb{E}\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}\\}_{k=1}^{m}\subset\mathcal{X}$.

Otherwise, we meet all the assumptions of the Theorem with
$B:=\frac{\gamma^{2}s}{m}$, and need only to compute $\sigma^{2}$. To this end
we apply the second part of Lemma 1, and obtain for any $v\in\mathcal{D}_{s}$

$\displaystyle\mathbb{E}\sum_{k=1}^{m}f_{v}(\mathcal{A}_{k*}^{*}\mathcal{A}_{k*}-\mathbb{E}\mathcal{A}_{k*}^{*}\mathcal{A}_{k*})^{2}=\mathbb{E}\sum_{k=1}^{m}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{4}-\left(\mathbb{E}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{2}\right)^{2}$
$\displaystyle\leq\mathbb{E}\sum_{k=1}^{m}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{4}\leq\frac{s\beta^{2}\gamma^{2}(1+2\theta)}{m}:=\sigma^{2}.$

To finish, apply Theorem 8 with $t=\tau>0$ and assume

$\frac{\sqrt{m}}{\sqrt{\log(m)}}\geq
2c\gamma\sqrt{s}\log(s)\sqrt{\log(n)}\frac{\sqrt{\beta^{2}+\beta^{2}\theta+\tau}}{\tau},$

which gives $\mathbb{E}X_{m}\leq\tau$ according to (28). Then

$X_{m}\leq\mathbb{E}X_{m}+t\leq 2\tau$

with probability exceeding

$1-\exp\left(-\frac{m\tau}{4s\gamma^{2}}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+\beta^{2}(1+2\theta)}\right)\right)\right).$

Thus, by our remarks in the beginning of the section,
$\mathcal{A}:=\frac{\sqrt{N}}{\sqrt{m}}\mathcal{S}\Psi$ satifies the RIP
inequality with constant
$2\tau+\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}$
and the derived probability.

∎

####  8.1.2 Proof of Theorem 5

We finish by applying the RIC bound established in the previous section to
well known results in compressive sensing and obtain our recovery error bound.

###### Proof of Theorem 5:.

Let $f=\Psi g$ and $e=\tilde{f}-\mathcal{S}f$. Then for any
$\bar{g}\in\mathbb{C}^{n}$

$\displaystyle\|\mathcal{S}\Psi\bar{g}-\tilde{f}-d\|_{2}=\|\mathcal{S}\Psi\bar{g}-(\mathcal{S}f+e)-d\|_{2}$
$\displaystyle=$ $\displaystyle\|\mathcal{S}\Psi\bar{g}-\mathcal{S}\Psi
g-e-d\|_{2}.$

Thus, $g^{\sharp}$ equivalently solves

$\displaystyle\min_{\bar{g}\in\mathbb{C}^{n}}\|\bar{g}\|_{1}\hskip
5.0pt\mbox{subject to}\hskip 5.0pt\|\mathcal{S}\Psi\bar{g}-\mathcal{S}\Psi
g-e-d\|_{2}\leq\eta+2\sqrt{m}\sum_{\lvert\ell\rvert\geq\frac{N-1}{2}}\lvert
c_{\ell}\rvert$ $\displaystyle=$
$\displaystyle\min_{\bar{g}\in\mathbb{C}^{n}}\|\bar{g}\|_{1}\hskip
5.0pt\mbox{subject to}\hskip
5.0pt\Biggl{\|}\mathcal{A}(\bar{g}-g)-\frac{\sqrt{N}}{\sqrt{m}}(e+d)\Biggr{\|}_{2}\leq\frac{\sqrt{N}}{\sqrt{m}}\left(\eta+2\sqrt{m}\sum_{\lvert\ell\rvert\geq\frac{N-1}{2}}\lvert
c_{\ell}\rvert\right),$

where as before $\mathcal{A}:=\frac{\sqrt{N}}{\sqrt{m}}\mathcal{S}\Psi$.
Notice that $g$ is feasible since

$\displaystyle\frac{\sqrt{N}}{\sqrt{m}}\|e+d\|_{2}\leq\frac{\sqrt{N}}{\sqrt{m}}\|e\|_{2}+\frac{\sqrt{N}}{\sqrt{m}}\|d\|_{2}\leq\frac{\sqrt{N}}{\sqrt{m}}\left(2\sqrt{m}\sum_{\lvert\ell\rvert\geq\frac{N-1}{2}}\lvert
c_{\ell}\rvert\right)+\frac{\sqrt{N}\eta}{\sqrt{m}},$

where the last equality follows since $\|e\|_{2}\leq
2\sqrt{m}\sum_{\lvert\ell\rvert\geq\frac{N-1}{2}}\lvert c_{\ell}\rvert$ by
Theorem 1. Apply Theorem 7 with

$\tau:=\frac{1}{4\sqrt{2}}-\frac{\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}}{4}>0,$

which is positive by our assumptions on $\theta,\alpha$ and $\beta$. Then if

$\frac{\sqrt{m}}{\sqrt{\log(m)}}\geq
2c\gamma\sqrt{2s}\log(2s)\sqrt{\log(n)}\frac{\sqrt{\beta^{2}+\beta^{2}\theta+\tau}}{\tau},$
(29)

$\mathcal{A}$ satisfies the $2s$-restricted isometry property with constant
satisfying

$\displaystyle\delta_{2s}\leq
2\tau+\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}$
$\displaystyle:=$
$\displaystyle\frac{1}{2\sqrt{2}}+\frac{\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}}{2}<\frac{1}{\sqrt{2}}$

and probability exceeding

$1-\exp\left(-\frac{m\tau}{8s\gamma^{2}}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+\beta^{2}(1+2\theta)}\right)\right)\right).$

Since $\delta_{2s}<\frac{1}{\sqrt{2}}$, Theorem 6 and the singular values of
$\Psi$ give

$\displaystyle\|f-\Psi
g^{\sharp}\|_{2}=\|\Psi(g-g^{\sharp})\|_{2}\leq\beta\|g-g^{\sharp}\|_{2}$
$\displaystyle\leq\frac{\beta\tilde{c}_{2}\epsilon_{s}(g)}{\sqrt{s}}+\beta\tilde{c}_{3}\frac{\sqrt{N}}{\sqrt{m}}\left(\eta+2\sqrt{m}\sum_{\lvert\ell\rvert\geq\frac{N-1}{2}}\lvert
c_{\ell}\rvert\right),$

where $\tilde{c}_{2},\tilde{c}_{3}>0$ depend only on
$\frac{1}{2\sqrt{2}}+\frac{\max\\{\beta^{2}\theta+\beta^{2}-1,\alpha^{2}\theta+1-\alpha^{2}\\}}{2}$.

In particular, since $s\leq n$ and by assumption $m\leq n$, if

$\sqrt{m}\geq
2c\gamma\sqrt{2s}\log^{2}(2n)\frac{\sqrt{\beta^{2}+\beta^{2}\theta+\tau}}{\tau},$
(30)

then (29) holds. In the statement of the theorem we have absorbed all the
absolute constants in this expression into the definition of $\tilde{c}_{1}$.
∎

####  8.1.3 Proof of Lemma 1

We end this section by establishing Lemma 1. This lemma was used in the RIC
bound to apply the respective concentration inequality.

###### Proof of Lemma 1:.

We upper bound $\lvert\langle\mathcal{A}_{k*},v\rangle\rvert$ as follows,

$\displaystyle\lvert\langle\mathcal{A}_{k*},v\rangle\rvert=\Biggl{\lvert}\langle\frac{1}{\sqrt{m}}\tilde{\mathcal{N}}_{k*},\mathcal{F}^{*}\Psi
v\rangle\Biggr{\rvert}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{m}}\Biggl{\lvert}\sum_{\ell=1}^{N}\tilde{\mathcal{N}}_{k\ell}(\mathcal{F}^{*}\Psi
v)_{\ell}\Biggr{\rvert}\leq\frac{1}{\sqrt{m}}\sum_{\ell=1}^{N}\lvert\tilde{\mathcal{N}}_{k\ell}\rvert\lvert(\mathcal{F}^{*}\Psi
v)_{\ell}\rvert$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{m}}\sum_{\ell=1}^{N}\lvert(\mathcal{F}^{*}\Psi
v)_{\ell}\rvert=\frac{1}{\sqrt{m}}\sum_{\ell=1}^{N}\Biggl{\lvert}\sum_{p=1}^{n}\left(\mathcal{F}^{*}\Psi\right)_{\ell
p}v_{p}\Biggr{\rvert}=\frac{1}{\sqrt{m}}\sum_{\ell=1}^{N}\Biggl{\lvert}\sum_{p=1}^{n}\langle\mathcal{F}_{*\ell},\Psi_{*p}\rangle
v_{p}\Biggr{\rvert}$ $\displaystyle\leq$
$\displaystyle\frac{1}{\sqrt{m}}\sum_{p=1}^{n}\sum_{\ell=1}^{N}\lvert\langle\mathcal{F}_{*\ell},\Psi_{*p}\rangle\rvert\lvert
v_{p}\rvert\leq\frac{1}{\sqrt{m}}\left(\max_{q\in[n]}\sum_{\ell=1}^{N}\lvert\langle\mathcal{F}_{*\ell},\Psi_{*q}\rangle\rvert\right)\sum_{p=1}^{n}\lvert
v_{p}\rvert$ $\displaystyle:=$
$\displaystyle\frac{\gamma\|v\|_{1}}{\sqrt{m}}\leq\frac{\gamma\sqrt{s}\|v\|_{2}}{\sqrt{m}}.$

Note that we have shown $\|\mathcal{F}^{*}\Psi
v\|_{1}\leq\gamma\sqrt{s}\|v\|_{2}$ (this will be used in the proof of the
remaining claim).

We now bound the sum of the fourth moments, using orthogonality and our
deviation model as in the proof of (27). Throughout, let
$\tilde{\Delta}\in\mathbb{R}$ be an independent copy of the entries of
$\Delta\in\mathbb{R}^{m}$. Expanding and taking expected value we obtain

$\displaystyle\mathbb{E}\sum_{k=1}^{m}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{4}=\frac{1}{m^{2}}\mathbb{E}\sum_{k=1}^{m}\lvert\langle\tilde{\mathcal{N}}_{k*},\mathcal{F}^{*}\Psi
v\rangle\rvert^{4}=\frac{1}{m^{2}}\mathbb{E}\sum_{k=1}^{m}\Biggl{\lvert}\sum_{\ell=1}^{N}\tilde{\mathcal{N}}_{k\ell}(\mathcal{F}^{*}\Psi
v)_{\ell}\Biggr{\rvert}^{4}$
$\displaystyle=\frac{1}{m^{2}}\mathbb{E}\sum_{k=1}^{m}\left(\sum_{\ell_{1},\ell_{2},\ell_{3},\ell_{4}=1}^{N}\tilde{\mathcal{N}}_{k\ell_{1}}\overline{\tilde{\mathcal{N}}}_{k\ell_{2}}\tilde{\mathcal{N}}_{k\ell_{3}}\overline{\tilde{\mathcal{N}}}_{k\ell_{4}}(\mathcal{F}^{*}\Psi
v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{4}}\right)$
$\displaystyle:=\sum_{\ell_{1},\ell_{2},\ell_{3},\ell_{4}=1}^{N}\frac{1}{m^{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{4}}\left(\mathbb{E}\sum_{k=1}^{m}\textbf{e}(\tilde{t}_{k}(\ell_{1}-\ell_{2}+\ell_{3}-\ell_{4}))\right)$
$\displaystyle=\sum_{\ell_{1}-\ell_{2}=\ell_{4}-\ell_{3}}\frac{1}{m}(\mathcal{F}^{*}\Psi
v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{4}}$
$\displaystyle+\sum_{j=1}^{\lfloor
2(N-1)/m\rfloor}\sum_{\ell_{1}-\ell_{2}=\ell_{4}-\ell_{3}+jm}\frac{1}{m}(\mathcal{F}^{*}\Psi
v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{4}}\mathbb{E}\textbf{e}(jm(\tilde{\Delta}-1/2))$
$\displaystyle+\sum_{j=1}^{\lfloor
2(N-1)/m\rfloor}\sum_{\ell_{1}-\ell_{2}=\ell_{4}-\ell_{3}-jm}\frac{1}{m}(\mathcal{F}^{*}\Psi
v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{4}}\mathbb{E}\textbf{e}(-jm(\tilde{\Delta}-1/2)).$

The last equality holds as in the proof of (27), where

$\displaystyle\mathbb{E}\frac{1}{m}\sum_{k=1}^{m}\textbf{e}(\tilde{t}_{k}(\ell_{1}-\ell_{2}+\ell_{3}-\ell_{4}))$
$\displaystyle=\left\\{\begin{array}[]{ll}1&\mbox{if
}\ell_{1}-\ell_{2}+\ell_{3}-\ell_{4}=0\\\
\mathbb{E}\textbf{e}\left(jm(\tilde{\Delta}-1/2)\right)&\mbox{if
}\ell_{1}-\ell_{2}+\ell_{3}-\ell_{4}=jm,j\in\mathbb{Z}/\\{0\\}\\\
0&\mbox{otherwise.}\end{array}\right.$

We bound each of these last three terms. The first term can be bounded as

$\displaystyle\Biggl{\lvert}\sum_{\ell_{1}-\ell_{2}=\ell_{4}-\ell_{3}}\frac{1}{m}(\mathcal{F}^{*}\Psi
v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{4}}\Biggr{\rvert}$
$\displaystyle=\Biggl{\lvert}\frac{1}{m}\sum_{\ell_{3},\ell_{4}=1}^{N}(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi v)}_{\ell_{4}}\sum_{\ell_{2}\in
Q_{\ell_{3},\ell_{4}}}(\mathcal{F}^{*}\Psi
v)_{\ell_{4}+\ell_{2}-\ell_{3}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{2}}\Biggr{\rvert}$
$\displaystyle\leq\frac{1}{m}\sum_{\ell_{3},\ell_{4}=1}^{N}\lvert(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\rvert\lvert(\mathcal{F}^{*}\Psi
v)_{\ell_{4}}\rvert\Biggl{\lvert}\sum_{\ell_{2}\in
Q_{\ell_{3},\ell_{4}}}(\mathcal{F}^{*}\Psi
v)_{\ell_{4}+\ell_{2}-\ell_{3}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{2}}\Biggr{\rvert}$
$\displaystyle\leq\frac{1}{m}\sum_{\ell_{3},\ell_{4}=1}^{N}\lvert(\mathcal{F}^{*}\Psi
v)_{\ell_{3}}\rvert\lvert(\mathcal{F}^{*}\Psi
v)_{\ell_{4}}\rvert\lvert\langle\mathcal{F}^{*}\Psi v,\mathcal{F}^{*}\Psi
v\rangle\rvert$ $\displaystyle=\frac{1}{m}\|\mathcal{F}^{*}\Psi
v\|_{2}^{2}\|\mathcal{F}^{*}\Psi v\|_{1}^{2}=\frac{1}{m}\|\Psi
v\|_{2}^{2}\|\mathcal{F}^{*}\Psi v\|_{1}^{2}$
$\displaystyle\leq\frac{1}{m}\beta^{2}\|v\|_{2}^{4}\gamma^{2}s.$

In the first equality, $Q_{\ell_{3},\ell_{4}}$ is the set of allowed index
values for $\ell_{2}$ according to $\ell_{3}$ and $\ell_{4}$, i.e., so that
$\ell_{2}\in[N]$ and $\ell_{4}+\ell_{2}-\ell_{3}\in[N]$. The last bound
consists of bounding $\|\mathcal{F}^{*}\Psi
v\|_{1}\leq\gamma\sqrt{s}\|v\|_{2}$ as in the computations of
$\lvert\langle\mathcal{A}_{k*},v\rangle\rvert$, and $\|\Psi
v\|_{2}\leq\beta\|v\|_{2}$ from the singular values of $\Psi$.

The remaining terms can be bounded similarly. Using our deviation model,
notice that

$\sum_{j=1}^{\lfloor\frac{2(N-1)}{m}\rfloor}\lvert\mathbb{E}\textbf{e}(\pm
jm(\tilde{\Delta}-1/2))\rvert\leq\sum_{j=1}^{\lfloor
2(N-1)/m\rfloor}\frac{\theta m}{2N}\leq\theta,$

and therefore we can similarly obtain

$\displaystyle\Biggl{\lvert}\sum_{j=1}^{\lfloor\frac{2(N-1)}{m}\rfloor}\sum_{\ell_{1}-\ell_{2}=\ell_{4}-\ell_{3}\pm
jm}\frac{1}{m}(\mathcal{F}^{*}\Psi v)_{\ell_{1}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{2}}(\mathcal{F}^{*}\Psi v)_{\ell_{3}}\overline{(\mathcal{F}^{*}\Psi
v)}_{\ell_{4}}\mathbb{E}\textbf{e}(\pm jm(\tilde{\Delta}-1/2))\Biggr{\rvert}$
$\displaystyle\leq\frac{\theta}{m}\beta^{2}\|v\|_{2}^{4}\gamma^{2}s.$

This gives,

$\mathbb{E}\sum_{k=1}^{m}\lvert\langle\mathcal{A}_{k*},v\rangle\rvert^{4}\leq\frac{1+2\theta}{m}\beta^{2}\|v\|_{2}^{4}\gamma^{2}s,$

as desired.

∎

###  8.2 Proof of Theorem 4

To establish the claim, we aim to show that

$\inf_{v\in S^{N-1}}\|\mathcal{S}v\|_{2}\geq\frac{1}{\delta},$ (31)

holds with high probability. By optimality of $f^{\sharp}$, this will give

$\displaystyle\|f-f^{\sharp}\|_{2}\leq\delta\|\mathcal{S}(f-f^{\sharp})\|_{2}\leq\delta\|\mathcal{S}f-b\|_{2}+\delta\|b-\mathcal{S}f^{\sharp}\|_{2}$
$\displaystyle\leq 2\delta\|\mathcal{S}f-b\|_{2}\leq
2\delta\|d\|_{2}+2\delta\sqrt{m}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert,$

where the last inequality is due to our noise model and trigonometric
interpolation error (Theorem 1).

To this end, we begin by noting that our sampling operator is isometric in the
sense that for any $v\in\mathbb{C}^{N}$

$\mathbb{E}\|\mathcal{S}v\|_{2}^{2}=\mathbb{E}\|\mathcal{N}v\|_{2}^{2}=\frac{m}{N}\|v\|_{2}^{2}.$
(32)

For brevity, we omit the proof of this statement but remark that it can be
shown as in the calculation of (27) by noting that
$\lvert\ell-\tilde{\ell}\rvert\leq N-1<m$ for all $\ell,\tilde{\ell}\in[N]$
and by assumption (20) on $m$.

We normalize by letting
$\tilde{\mathcal{S}}=\frac{1}{\sqrt{m}}\tilde{\mathcal{N}}\mathcal{F}^{*}:=\frac{\sqrt{N}}{\sqrt{m}}\mathcal{N}\mathcal{F}^{*}$
and proceed to the main component of our result.

######  Theorem 9.

Let the entries of $\Delta$ be i.i.d. with any distribution. If

$\sqrt{m}\geq 6\sqrt{N\log(N)}\frac{\sqrt{1+\tau}}{\tau},$

then

$\|\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}-\mathbb{E}\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}\|\leq
2\tau$

with probability exceeding

$1-\exp\left(-\frac{m\tau}{4N}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+1}\right)\right)\right).$

###### Proof:.

Notice that

$\displaystyle
X_{m}:=\|\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}-\mathbb{E}\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}\|:=\Biggl{\|}\mathcal{F}\frac{1}{m}\left(\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}-\mathbb{E}\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}\right)\mathcal{F}^{*}\Biggr{\|}$
$\displaystyle=\frac{1}{m}\|\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}-\mathbb{E}\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}\|=\frac{1}{m}\Biggl{\|}\sum_{k=1}^{m}\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right)\Biggr{\|}.$

We begin by bounding $\mathbb{E}_{\mathcal{D}}X_{m}$. Symmetrizing as Lemma
6.3 in [47], we have

$\mathbb{E}_{\mathcal{D}}X_{m}\leq\frac{2}{m}\mathbb{E}_{\mathcal{D}}\mathbb{E}_{\epsilon}\Biggl{\|}\sum_{k=1}^{m}\epsilon_{k}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\Biggr{\|},$

where the $\epsilon_{k}$ are independent Rademacher random variables. Notice
that for $k\in[m]$

$\displaystyle\|\tilde{\mathcal{N}}_{k*}\|_{2}=\sqrt{N},$

so that applying an explicit form of Rudelson’s Lemma (Theorem 1 in [68], with
$p=1$) with this observation gives that

$\displaystyle\mathbb{E}_{\epsilon}\Biggl{\|}\sum_{k=1}^{m}\epsilon_{k}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\Biggr{\|}\leq\left(\sqrt{2\log(2N)}+\sqrt{\frac{\pi}{2}}\right)\|\tilde{\mathcal{N}}_{k*}\|_{2}\Biggl{\|}\sum_{k=1}^{m}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\Biggr{\|}^{1/2}$
$\displaystyle\leq
3\sqrt{N\log(2N)}\Biggl{\|}\sum_{k=1}^{m}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\Biggr{\|}^{1/2},$

where the last inequality holds since $\frac{\pi}{2}\leq\log(2N)$ for $N\geq
2$. We may therefore bound

$\displaystyle\mathbb{E}_{\mathcal{D}}X_{m}\leq\frac{6\sqrt{N\log(2N)}}{m}\mathbb{E}_{\mathcal{D}}\Biggl{\|}\sum_{k=1}^{m}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\Biggr{\|}^{1/2}$
$\displaystyle\leq\frac{6\sqrt{N\log(2N)}}{\sqrt{m}}\left(\mathbb{E}_{\mathcal{D}}\Biggl{\|}\sum_{k=1}^{m}\frac{1}{m}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\Biggr{\|}\right)^{1/2}$
$\displaystyle\leq\frac{6\sqrt{N\log(2N)}}{\sqrt{m}}\left(\mathbb{E}_{\mathcal{D}}\Biggl{\|}\sum_{k=1}^{m}\frac{1}{m}\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}_{\mathcal{D}}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right)\Biggr{\|}+\frac{1}{m}\|\mathbb{E}_{\mathcal{D}}\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}\|\right)^{1/2}$
$\displaystyle=\frac{6\sqrt{N\log(2N)}}{\sqrt{m}}\left(\mathbb{E}_{\mathcal{D}}X_{m}+1\right)^{1/2}.$

The last equality holds by definition of $X_{m}$ and isometry (32), since for
Hermitian operators

$\|\mathbb{E}_{\mathcal{D}}\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}\|=\sup_{v\in
S^{N-1}}\mathbb{E}_{\mathcal{D}}\langle\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}v,v\rangle=\sup_{v\in
S^{N-1}}\mathbb{E}_{\mathcal{D}}\|\tilde{\mathcal{N}}v\|_{2}^{2}=m.$

Concluding,

$\frac{\mathbb{E}_{\mathcal{D}}X_{m}}{\sqrt{\mathbb{E}_{\mathcal{D}}X_{m}+1}}\leq\frac{6\sqrt{N\log(2N)}}{\sqrt{m}}$

and we achieve $\mathbb{E}_{\mathcal{D}}X_{m}\leq\tau$ if

$\frac{6\sqrt{N\log(2N)}}{\sqrt{m}}\leq\frac{\tau}{\sqrt{\tau+1}}.$ (33)

We apply Theorem 8 from the previous section to show that $X_{m}$ is close to
its expected value with high probability. To this end, define

$\mathcal{X}:=\\{Z\in\mathbb{C}^{N\times N}:\langle Zv,v\rangle\in\mathbb{R}\
\ \forall v\in S^{N-1}\ \ \mbox{and}\ \ \sup_{v\in S^{N-1}}\langle
Zv,v\rangle\leq N\\},$

which is a closed subset of $\mathbb{C}^{N\times N}$ (a Polish space), and
therefore itself a Polish space (see [25]). For $v\in S^{N-1}$, define the
function $f_{v}:\mathcal{X}\mapsto\mathbb{R}$ as

$f_{v}(Z):=\frac{1}{m}\langle Zv,v\rangle,$

which gives a collection of real valued functions by definition of
$\mathcal{X}$. Then notice that

$\displaystyle
X_{m}=\frac{1}{m}\Biggl{\|}\sum_{k=1}^{m}\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right)\Biggr{\|}$
$\displaystyle=$ $\displaystyle\frac{1}{m}\sup_{v\in
S^{N-1}}\sum_{k=1}^{m}\langle\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right)v,v\rangle$
$\displaystyle=$ $\displaystyle\sup_{v\in
S^{N-1}}\sum_{k=1}^{m}f_{v}\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right)$
$\displaystyle=$ $\displaystyle\sup_{v\in
S^{*}}\sum_{k=1}^{m}f_{v}\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right),$

where $S^{*}$ is a dense countable subset of $S^{N-1}$ and the second equality
holds since
$\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}-\mathbb{E}\tilde{\mathcal{N}}^{*}\tilde{\mathcal{N}}$
is a Hermitian matrix. For all $k\in[m]$ we have

$\displaystyle
f_{v}(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}):=\frac{1}{m}\langle\left(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\frac{1}{m}\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\right)v,v\rangle$
$\displaystyle=\frac{1}{m}\lvert\langle\tilde{\mathcal{N}}_{k*},v\rangle\rvert^{2}-\frac{1}{m}\mathbb{E}\lvert\langle\tilde{\mathcal{N}}_{k*},v\rangle\rvert^{2}\leq\frac{\|\mathcal{N}_{k*}\|_{2}^{2}\|v\|_{2}^{2}}{m}=\frac{N}{m},$

so that
$\\{\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}\\}_{k=1}^{m}\subset\mathcal{X}$
and we may take $B=\frac{N}{m}$. We now compute $\sigma^{2}$.

$\displaystyle\mathbb{E}\sum_{k=1}^{m}f_{v}(\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*}-\mathbb{E}\tilde{\mathcal{N}}_{k*}^{*}\tilde{\mathcal{N}}_{k*})^{2}=\frac{1}{m^{2}}\mathbb{E}\sum_{k=1}^{m}\lvert\langle\tilde{\mathcal{N}}_{k*},v\rangle\rvert^{4}-\left(\mathbb{E}\lvert\langle\tilde{\mathcal{N}}_{k*},v\rangle\rvert^{2}\right)^{2}$
$\displaystyle\leq\frac{1}{m^{2}}\mathbb{E}\sum_{k=1}^{m}\lvert\langle\tilde{\mathcal{N}}_{k*},v\rangle\rvert^{4}:=\frac{1}{m^{2}}\mathbb{E}\sum_{k=1}^{m}\Biggl{\lvert}\sum_{p=1}^{N}\textbf{e}\left(\tilde{t}_{k}(p-\tilde{N}-1)\right)\bar{v}_{p}\Biggr{\rvert}^{4}$
$\displaystyle=\frac{1}{m^{2}}\sum_{p_{1},p_{2},p_{3},p_{4}=1}^{N}\bar{v}_{p_{1}}v_{p_{2}}\bar{v}_{p_{3}}v_{p_{4}}\left(\mathbb{E}\sum_{k=1}^{m}\textbf{e}\left(\tilde{t}_{k}(p_{1}-p_{2}+p_{3}-p_{4})\right)\right)$
$\displaystyle=\frac{1}{m}\sum_{p_{1}-p_{2}+p_{3}-p_{4}=0}\bar{v}_{p_{1}}v_{p_{2}}\bar{v}_{p_{3}}v_{p_{4}}$
$\displaystyle=\frac{1}{m}\sum_{p_{1},p_{2}=1}^{N}\sum_{p_{3}\in
Q_{p_{1},p_{2}}}\bar{v}_{p_{1}}v_{p_{2}}\bar{v}_{p_{3}}v_{p_{1}-p_{2}+p_{3}}$
$\displaystyle\leq\frac{1}{m}\sum_{p_{1},p_{2}=1}^{N}\lvert
v_{p_{1}}v_{p_{2}}\rvert\lvert\langle
v,v\rangle\rvert=\frac{\|v\|_{1}^{2}\|v\|_{2}^{2}}{m}\leq\frac{N}{m}:=\sigma^{2}.$

In the fourth equality $Q_{p_{1},p_{2}}$ denotes the set of allowed indices
for $p_{3}$ depending on $p_{1}$ and $p_{2}$, i.e., so that
$p_{3},p_{1}-p_{2}+p_{3}\in[N]$. The third equality is due to orthogonality of
the complex exponentials, similar to our isometry calculation (32) where our
assumption (20) gives that $m\geq 2N>\lvert p_{1}-p_{2}+p_{3}-p_{4}\rvert$ for
all $p_{1},p_{2},p_{3},p_{4}\in[N]$.

To finish, pick $t=\tau>0$ in Theorem 8 and assume

$\sqrt{m}\geq 6\sqrt{N\log(2N)}\frac{\sqrt{1+\tau}}{\tau},$

which gives $\mathbb{E}X_{m}\leq\tau$ according to (33). Then

$X_{m}\leq\mathbb{E}X_{m}+t\leq 2\tau$

with probability exceeding

$1-\exp\left(-\frac{m\tau}{4N}\log\left(1+2\log\left(1+\frac{\tau}{2\tau+1}\right)\right)\right).$

∎

With our remarks in the beginning of the section, we can now easily establish
the main result.

###### Proof of Theorem 4:.

Under our assumptions, apply Theorem 9 to obtain that

$\|\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}-\mathbb{E}\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}\|\leq
2\tau$

with the prescribed probability. Therefore, since
$\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}-\mathbb{E}\tilde{\mathcal{S}}^{*}\tilde{\mathcal{S}}$
is a Hermitian matrix and by isometry (32), for all $v\in S^{N-1}$

$\lvert\|\tilde{\mathcal{S}}v\|_{2}^{2}-\mathbb{E}\|\tilde{\mathcal{S}}v\|_{2}^{2}\rvert=\lvert\|\tilde{\mathcal{S}}v\|_{2}^{2}-\|v\|_{2}^{2}\rvert\leq
2\tau.$

In particular, for our un-normalized ensemble, we have

$\inf_{v\in
S^{N-1}}\|\mathcal{S}v\|_{2}\geq\frac{\sqrt{1-2\tau}\sqrt{m}}{\sqrt{N}},$

where $1-2\tau>0$ by assumption. This establishes (31) with
$\delta=\frac{\sqrt{N}}{\sqrt{1-2\tau}\sqrt{m}}$. The remainder of the proof
follows from our outline in the beginning of the section. ∎

###  8.3 Interpolation Error of Dirichlet Kernel: Proof

In this section we provide the error term of our interpolation operator when
applied to our signal model (Theorem 1) and also the error bound given in
Corollary 1.

###### Proof of Theorem 1:.

We begin by showing (7), i.e., if $\tilde{t}_{k}=t_{\tilde{p}}$ for some
$\tilde{p}\in[N]$ (our “nonuniform” sample lies on the dense interpolation
grid) then the error is zero. This is easy to see by orthogonality of the
complex exponentials, combining (4), (5) (recall that
$\tilde{N}=\frac{N-1}{2}$) we have

$\displaystyle(\mathcal{S}f)_{k}=\langle
f,\mathcal{S}_{k*}\rangle=\sum_{p=1}^{N}f_{p}\mathcal{S}_{kp}:=\frac{1}{N}\sum_{p=1}^{N}f_{p}\left(\sum_{u=-\tilde{N}}^{\tilde{N}}\textbf{e}(ut_{p})\textbf{e}(-u\tilde{t}_{k})\right)$
$\displaystyle:=\frac{1}{N}\sum_{p=1}^{N}f_{p}\left(\sum_{u=-\tilde{N}}^{\tilde{N}}\textbf{e}(ut_{p})\textbf{e}(-ut_{\tilde{p}})\right)=f_{\tilde{p}}=\mathbf{f}(t_{\tilde{p}})=\mathbf{f}(\tilde{t}_{k})=\tilde{f}_{k}.$

We now deal with the general case (8). Recall the Fourier expansion of our
underlying function

$\mathbf{f}(x)=\sum_{\ell=-\infty}^{\infty}c_{\ell}\mathbf{e}(\ell x).$

Again, using (4), (5) and the Fourier expansion at $\textbf{f}(t_{p})=f_{p}$
we obtain

$\displaystyle(\mathcal{S}f)_{k}=\langle
f,\mathcal{S}_{k*}\rangle=\sum_{p=1}^{N}f_{p}\mathcal{S}_{kp}$
$\displaystyle:=\frac{1}{N}\sum_{p=1}^{N}\left(\sum_{\ell=-\infty}^{\infty}c_{\ell}\mathbf{e}(\ell
t_{p})\right)\left(\sum_{u=-\tilde{N}}^{\tilde{N}}\textbf{e}(ut_{p})\textbf{e}(-u\tilde{t}_{k})\right).$

At this point, we wish to switch the order of summation and sum over all
$p\in[N]$. We must assume the corresponding summands are non-zero. To this
end, we continue assuming $f_{p},\mathcal{S}_{kp}\neq 0$ for all $p\in[N]$. We
will deal with these cases separately afterward. In particular we will remove
this assumption for the $f_{p}$’s and show that $\mathcal{S}_{kp}\neq 0$ under
our assumption $\tilde{\tau}\subset\Omega$.

Proceeding, we may now sum over all $p\in[N]$ to obtain

$\displaystyle(\mathcal{S}f)_{k}$
$\displaystyle=\frac{1}{N}\sum_{u=-\tilde{N}}^{\tilde{N}}\sum_{\ell=-\infty}^{\infty}c_{\ell}\textbf{e}(-u\tilde{t}_{k})\sum_{p=1}^{N}\textbf{e}((u+\ell)t_{p})$
$\displaystyle=\sum_{u=-\tilde{N}}^{\tilde{N}}\sum_{j=-\infty}^{\infty}(-1)^{jN}c_{jN+u}\textbf{e}(u\tilde{t}_{k})=\sum_{j=-\infty}^{\infty}(-1)^{\lfloor\frac{j+\tilde{N}}{N}\rfloor}c_{j}\textbf{e}(r(j)\tilde{t}_{k}).$

The second equality is obtained by orthogonality of the exponential basis
functions, $\sum_{p=1}^{N}\textbf{e}((u+\ell)t_{p})=0$ when $\ell+u\notin
N\mathbb{Z}$ and otherwise equal to $N(-1)^{jN}$ for some $j\in\mathbb{Z}$
(where $u+\ell=jN$). The last equality results from a reordering of the series
where the mapping $r$ is defined as in the statement of Theorem 1. To
illustrate, we consider $j\geq 0$ (for simplicity) and first notice that
$(-1)^{jN}=(-1)^{j}$ since $N$ is assumed to be odd in Section 2.1.
Aesthetically expanding the previous sum gives

$\displaystyle\sum_{u=-\tilde{N}}^{\tilde{N}}\sum_{j=0}^{\infty}(-1)^{j}c_{jN+u}\textbf{e}(u\tilde{t}_{k})$
$\displaystyle=\hskip
15.0pt\textbf{e}(-\tilde{N}\tilde{t}_{k})\left(c_{-\tilde{N}}-c_{N-\tilde{N}}+c_{2N-\tilde{N}}-\cdots\right)$
$\displaystyle+\textbf{e}((-\tilde{N}+1)\tilde{t}_{k})\left(c_{-\tilde{N}+1}-c_{N-\tilde{N}+1}+c_{2N-\tilde{N}+1}-\cdots\right)$
$\displaystyle\cdots$ $\displaystyle+\hskip 44.0ptc_{0}-c_{N}+c_{2N}-\cdots$
$\displaystyle\cdots$ $\displaystyle+\hskip
23.0pt\textbf{e}(\tilde{N}\tilde{t}_{k})\left(c_{\tilde{N}}-c_{N+\tilde{N}}+c_{2N+\tilde{N}}-\cdots\right).$

Notice that in the first row we have indices $N-\tilde{N}=\tilde{N}+1$ and
$2N-\tilde{N}=N+\tilde{N}+1$. Therefore, if start at the term $c_{-\tilde{N}}$
and traverse this array of Fourier coefficients we will obtain the ordered
sequence
$\\{(-1)^{\lfloor\frac{j+\tilde{N}}{N}\rfloor}c_{j}\\}_{j=-\tilde{N}}^{\infty}$
(with no repetitions). The coefficients in row $q\in[N]$ correspond to
frequency value $-\tilde{N}+q-1$ and have indices of the form
$c_{pN-\tilde{N}+q-1}$ for some $p\in\mathbb{N}$. Let us check that for a
given index the mapping $r$ gives the correct frequency value, i.e.,
$r(pN-\tilde{N}+q-1)=-\tilde{N}+q-1$ for all $q\in[N]$:

$\displaystyle
r(pN-\tilde{N}+q-1):=\mbox{rem}(pN-\tilde{N}+q-1+\tilde{N},N)-\tilde{N}$
$\displaystyle=\mbox{rem}(pN+q-1,N)-\tilde{N}=q-1-\tilde{N}.$

We can therefore reorder the series as desired and incorporate the sum over
$j<0$ via the same logic to establish the equality. Since for
$\ell\in\\{-\tilde{N},-\tilde{N}+1,\cdots,\tilde{N}\\}$ we have $r(\ell)=\ell$
and $(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}=1$, we finally obtain

$\displaystyle\textbf{f}(\tilde{t}_{k})-(\mathcal{S}f)_{k}=\sum_{\lvert\ell\rvert>\tilde{N}}c_{\ell}\left(\textbf{e}(\ell\tilde{t}_{k})-(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}\textbf{e}(r(\ell)\tilde{t}_{k})\right).$

The definition of the $p$-norms along with the triangle inequality give the
remaining claim. In particular,

$\displaystyle\|\tilde{f}-\mathcal{S}f\|_{p}=\left(\sum_{k=1}^{m}\lvert\textbf{f}(\tilde{t}_{k})-(\mathcal{S}f)_{k}\rvert^{p}\right)^{1/p}$
$\displaystyle=\left(\sum_{k=1}^{m}\Biggl{\lvert}\sum_{\lvert\ell\rvert>\tilde{N}}c_{\ell}\left(\textbf{e}(\ell\tilde{t}_{k})-(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}\textbf{e}(r(\ell)\tilde{t}_{k})\right)\Biggr{\rvert}^{p}\right)^{1/p}\leq\left(\sum_{k=1}^{m}\left(\sum_{\lvert\ell\rvert>\tilde{N}}2\lvert
c_{\ell}\rvert\right)^{p}\right)^{1/p}$
$\displaystyle=\left(m\left(\sum_{\lvert\ell\rvert>\tilde{N}}2\lvert
c_{\ell}\rvert\right)^{p}\right)^{1/p}=2m^{1/p}\sum_{\lvert\ell\rvert>\tilde{N}}\lvert
c_{\ell}\rvert.$

This finishes the proof in the case $f_{p},\mathcal{S}_{kp}\neq 0$ for all
$p\in[N]$. To remove this condition for the $f_{p}$’s, we may find a real
number $\mu$ such that the function

$\mathbf{h}(x):=\mathbf{f}(x)+\mu=\sum_{\ell\in(-\infty,\infty)\cap\mathbb{Z}/\\{0\\}}c_{\ell}\mathbf{e}(\ell
x)+c_{0}+\mu$

is non-zero when $x\in\\{t_{p}\\}_{p=1}^{N}$. In particular notice that if we
define $h=f+\mu\in\mathbb{C}^{N}$, then $h_{p}\neq 0$ for all $p\in[N]$.
Therefore, only assuming now that $\mathcal{S}_{kp}\neq 0$ for $p\in[N]$, the
previous argument can be applied to conclude

$\displaystyle\textbf{h}(\tilde{t}_{k})-(\mathcal{S}h)_{k}=\sum_{\lvert\ell\rvert>\tilde{N}}c_{\ell}\left(\textbf{e}(\ell\tilde{t}_{k})-(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}\textbf{e}(r(\ell)\tilde{t}_{k})\right).$

However, if $1_{N}\in\mathbb{C}^{N}$ denotes the all ones vector and
$e_{\tilde{N}+1}\in\mathbb{C}^{N}$ is the $\tilde{N}+1$-th standard basis
vector, notice that

$\displaystyle(\mathcal{S}h)_{k}=\langle\mathcal{S}_{k*},h\rangle=\langle\mathcal{S}_{k*},f\rangle+\mu\langle\mathcal{S}_{k*},1_{N}\rangle=\langle\mathcal{S}_{k*},f\rangle+\mu\langle\mathcal{N}_{k*},\mathcal{F}^{*}1_{N}\rangle$
$\displaystyle=\langle\mathcal{S}_{k*},f\rangle+\mu\sqrt{N}\langle\mathcal{N}_{k*},e_{\tilde{N}+1}\rangle=\langle\mathcal{S}_{k*},f\rangle+\mu=(\mathcal{S}f)_{k}+\mu.$

The fourth equality holds by orthogonality of $\mathcal{F}^{*}$ and since
$\mathcal{F}_{(\tilde{N}+1)*}^{*}=\frac{1}{\sqrt{N}}1_{N}$. Therefore

$\textbf{h}(\tilde{t}_{k})-(\mathcal{S}h)_{k}=\textbf{f}(\tilde{t}_{k})+\mu-\left((\mathcal{S}f)_{k}+\mu\right)=\textbf{f}(\tilde{t}_{k})-(\mathcal{S}f)_{k},$

and the claim holds in this case as well.

The assumption $\mathcal{S}_{kp}\neq 0$ will always hold under our assuption
$\tilde{\tau}\subset\Omega$, i.e.,
$\tilde{t}_{k}\in[-\frac{1}{2},\frac{1}{2})$ for all $k\in[m]$. We show this
case by deriving conditions under which this occurs. As noted before, we have

$\displaystyle\mathcal{S}_{kp}:=\sum_{u=-\tilde{N}}^{\tilde{N}}\textbf{e}(u(t_{p}-\tilde{t}_{k}))=\sum_{u=0}^{N-1}\textbf{e}(u(t_{p}-\tilde{t}_{k}))\textbf{e}(-\tilde{N}(t_{p}-\tilde{t}_{k}))$
$\displaystyle=\textbf{e}(-\tilde{N}(t_{p}-\tilde{t}_{k}))\frac{1-\textbf{e}(N(t_{p}-\tilde{t}_{k}))}{1-\textbf{e}(t_{p}-\tilde{t}_{k})}$

and we see that $\mathcal{S}_{kp}=0$ iff
$N(t_{p}-\tilde{t}_{k})\in\mathbb{Z}/\\{0\\}$ and
$t_{p}-\tilde{t}_{k}\notin\mathbb{Z}$. However, notice that

$N(t_{p}-\tilde{t}_{k})=N\left(\frac{p-1}{N}-\frac{k-1}{m}-\Delta_{k}\right)=p-1-\frac{N(k-1)}{m}-N\Delta_{k},$

so that $N(t_{p}-\tilde{t}_{k})\in\mathbb{Z}/\\{0\\}$ iff
$\frac{N(k-1)}{m}+N\Delta_{k}=N\tilde{t}_{k}+\frac{N}{2}\in\mathbb{Z}/\\{p-1\\}$.
This condition equivalently requires $\tilde{t}_{k}=\frac{j}{N}-\frac{1}{2}$
for some $j\in\mathbb{Z}/\\{p-1\\}$. Since this must hold for all $p\in[N]$,
we have finally have that

$N(t_{p}-\tilde{t}_{k})\in\mathbb{Z}/\\{0\\}\ \ \ \mbox{iff}\ \ \
\tilde{t}_{k}=\frac{j}{N}-\frac{1}{2}\ \ \mbox{for some}\ \
j\in\mathbb{Z}/\\{0,1,\cdots,N-1\\}.$

We see that such a condition would imply that
$\tilde{t}_{k}\notin\Omega:=[-\frac{1}{2},\frac{1}{2})$, which violates our
assumption $\tilde{\tau}\subset\Omega$. This finishes the proof.

∎

We end this section with the proof of Corollary 1.

###### Proof of Corollary 1:.

The proof will consist of applying Theorem 2 (under identical assumptions) and
Theorem 1.

By Theorem 2, we have that

$\|f-\Psi
g^{\sharp}\|_{2}\leq\frac{c_{2}\epsilon_{s}(g)}{\sqrt{s}}+\frac{c_{3}\sqrt{N}\eta}{\sqrt{m}}+2c_{3}\sqrt{N}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert$

with probability exceeding (14). By Theorem 1, for $x\in\Omega$

$\textbf{f}(x)-\langle\textbf{h}(x),\mathcal{F}^{*}f\rangle=\sum_{\lvert\ell\rvert>\tilde{N}}c_{\ell}\left(\textbf{e}(\ell
x)-(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}\textbf{e}(r(\ell)x)\right).$

Therefore

$\displaystyle\lvert\textbf{f}(x)-\textbf{f}^{\sharp}(x)\rvert:=\lvert\textbf{f}(x)-\langle\textbf{h}(x),\mathcal{F}^{*}\Psi
g^{\sharp}\rangle\rvert$
$\displaystyle\leq\lvert\textbf{f}(x)-\langle\textbf{h}(x),\mathcal{F}^{*}f\rangle\rvert+\lvert\langle\textbf{h}(x),\mathcal{F}^{*}f\rangle-\langle\textbf{h}(x),\mathcal{F}^{*}\Psi
g^{\sharp}\rangle\rvert$
$\displaystyle\leq\Biggl{\lvert}\sum_{\lvert\ell\rvert>\tilde{N}}c_{\ell}\left(\textbf{e}(\ell
x)-(-1)^{\lfloor\frac{\ell+\tilde{N}}{N}\rfloor}\textbf{e}(r(\ell)x)\right)\Biggr{\rvert}+\|\textbf{h}(x)\|_{2}\|\mathcal{F}^{*}(f-\Psi
g^{\sharp})\|_{2}$ $\displaystyle\leq 2\sum_{\lvert\ell\rvert>\tilde{N}}\lvert
c_{\ell}\rvert+\frac{c_{2}\epsilon_{s}(g)}{\sqrt{s}}+\frac{c_{3}\sqrt{N}\eta}{\sqrt{m}}+2c_{3}\sqrt{N}\sum_{\lvert\ell\rvert>\frac{N-1}{2}}\lvert
c_{\ell}\rvert.$

The last inequality holds since $\|\textbf{h}(x)\|_{2}=1$ (here $x$ is
considered fixed and $\textbf{h}(x)\in\mathbb{C}^{N}$). This finishes the
proof.

∎

## References

  * [1] A.J. Jerri [1977], The Shannon Sampling Theorem-its various extensions and applications: a tutorial review, _Proc. IEEE_ , 65(11), 1565-1596.
  * [2] A.I. Zayed [1993], Advances in Shannon’s Sampling Theory. _CRC Press_.
  * [3] A. Griffin and P. Tsakalides [2008], Compressed sensing of audio signals using multiple sensors. _16th European Signal Processing Conference_.
  * [4] A.V. Oppenheim and R.W. Schafer [2009], Discrete-Time Signal Processing, _Prentice Hall Press_ , 3rd edition.
  * [5] A.Y. Aravkin, R. Kumar, H. Mansour, B. Recht and F.J. Herrmann [2014], Fast methods for denoising matrix completion formulations, with applications to robust seismic data interpolation. _SIAM Journal on Scientific Computing_ , 36, S237-S266.
  * [6] A. Zandieh, A. Zareian, M. Azghani and F. Marvasti [2014], Reconstruction of Sub-Nyquist Random Sampling for Sparse and Multi-Band Signals. _arXiv_.
  * [7] C.E. Shannon [1949], Communication in the Presence of Noise, _Proc. IRE_ , 37(1), 10-21.
  * [8] C. Li, C.C. Mosher and S.T. Kaplan [2012], Interpolated compressive sensing for seismic data reconstruction. _82nd annual international meeting, SEG, Expanded Abstracts_.
  * [9] D.R. Bellhouse [1981], Area Estimation by Point-counting Techniques. _Biometrics_ , 37(2), 303-312.
  * [10] D.P. Mitchell [1987], Generating Antialiased Images at Low Sampling Densities. _ACM SIGGRAPH Computer Graphics_ , 21(4), 65-72.
  * [11] D.P. Mitchell [1990], The Antialiasing Problem in Ray Tracing. _SIGGRAPH 90_.
  * [12] D.P. Dobkin, D. Eppstein and D.P. Mitchell [1996], Computing the Discrepancy with Applications to SuperSampling Patterns. _ACM Transactions on Graphics_ , 15(4), 354-376.
  * [13] D.M. Bechir and B. Ridha [2009], Non-uniform Sampling Schemes for RF Bandpass Sampling Receiver. _International Conference on Signal Processing Systems_.
  * [14] E.T. Whittaker [1915], On the Functions Which are Represented by the Expansion of Interpolating Theory, _Proc. Roy. Soc. Edinburgh_ , 35, 181-194.
  * [15] E. Shlomot and Y.Y. Zeevi [1989], A Nonuniform Sampling and Representation Scheme for Images Which Are Not bandlimited. _The Sixteenth Conference of Electrical and Electronics Engineers in Israel_.
  * [16] E. Margolis and Y.C. Eldar [2008], Nonuniform Sampling of Periodic Bandlimited Signals. _IEEE Transactions on Signal Processing_ , 56(7), 2728 - 2745.
  * [17] Ewout van den Berg and Michael P. Friedlander [2008], Probing the pareto frontier for basis pursuit solutions. _SIAM Journal on Scientific Computing_ , 31(2), 890-912.
  * [18] E. van den Berg and M.P. Friedlander [2014], Spot – A Linear-Operator Toolbox. [Online]. https://www.cs.ubc.ca/labs/scl/spot/index.html. [Accessed: 13- June- 2019].
  * [19] F.J. Beutler [1966], Error-Free Recovery of Signals from Irregularly Spaced Samples. _Society for Industrial and Applied Mathematics_ , 8(3), 328-335.
  * [20] F. Beutler [1970], Alias-free Randomly Timed Sampling of Stochastic Processes. _IEEE Transactions on Information Theory_ , 16(2), 147-152.
  * [21] F. Oberhettinger [1973], Fourier Transforms of Distributions and Their Inverses: A Collection of Tables, _Academic Press_.
  * [22] F. Marvasti [2001], Nonuniform sampling: theory and practice. _Springer_.
  * [23] F.J. Herrmann, M.P. Friedlander and Ö. Yılmaz [2012], Fighting the Curse of dimensionality: compressive sensing in exploration seismology. _IEEE Signal Processing Magazine_ , 29(3), 88 - 100.
  * [24] F. Krahmer and R. Ward [2014], Stable and Robust Sampling Strategies for Compressive Imaging. _IEEE Transactions on Image Processing_ , 23(2), 612-622.
  * [25] G. Beer [1991], A Polish Topology for the Closed Subsets of a Polish Space. _Proceedings of the American Mathematical Society_ , 113(4), 1123 - 1133.
  * [26] G.L. Bretthorst [2001], Nonuniform Sampling: Bandwidth and Aliasing. _AIP Conference Proceedings_ , 567(1).
  * [27] G. Hennenfent and F.J. Herrmann [2006], Seismic Denoising with Nonuniformly Sampled Curvelets. _Computing in Science and Engineering_ , 8(3), 16-25.
  * [28] G. Hennenfent and F.J. Herrmann [2008], Simply Denoise: Wavefield Reconstruction via Jittered Undersampling. _Geophysics_ , 73(3), V19-V28.
  * [29] G.E. Pfander [2015], Sampling Theory, A Renaissance, _Birkhäuser Basel_.
  * [30] H. Nyquist [1928], Certain Topics in Telegraph Transmission Theory, _AIEE Trans._ , 47, 617-644.
  * [31] H.S. Shapiro and R.A. Silverman [1960], Alias-free Sampling of Random Noise. _Journal of Society for Industrial and Applied Mathematics_ , 8(2), 225-248.
  * [32] H. Landau [1967], Necessary Density Condition for Sampling and Interpolation of Certain Entire Functions, _Acta Math_ , 117, 37-52.
  * [33] H. Lee and Z. Bien [2005], Sub-Nyquist Nonuniform Sampling and Perfect Reconstruction of Speech Signals. _TENCON 2005 - 2005 IEEE Region 10 Conference_.
  * [34] H. Rauhut [2008], Stability Results for Random Sampling of Sparse Trigonometric Polynomials. _IEEE Transactions on Information Theory_ , 54(12), 5661 - 5670.
  * [35] H. Rauhut [2011], Compressive Sensing and Structured Random Matrices. _Radon Series Comp. Appl._ , 1-94.
  * [36] H. Boche, R. Calderbank, G. Kutyniok, J. Vybíral [2013], Compressed sensing and its applications. _Birkhäuser_.
  * [37] J.M. Whittaker [1929], The Fourier Theory of the Cardinal Functions, _Proc. Math. Soc. Edinburgh_ , 1, 169-176.
  * [38] J. Keiner, S. Kunis and D. Potts [2008], Using NFFT 3 – a software library for various non-equispaced fast Fourier transforms. _ACM Trans. Math. Softw._ , 36, 19:1–19:30.
  * [39] J.B. Harley and J.M.F. Moura [2013], Sparse recovery of the multimodal and dispersive characteristics of Lamb waves. _The Journal of the Acoustic Society of America_ , 133(5), 2732-2745.
  * [40] J. Koh, W. Lee, T.K. Sarkar and M. Salazar-Palma [2013], Calculation of Far-Field Radiation Pattern Using Nonuniformly Spaced Antennas by a Least Square Method. _IEEE Transactions on Antennas and Propagation_ , 62(4), 1572-1578.
  * [41] J.A. Tropp [2015], Convex Recovery of a Structured Signal from Independent Random Linear Measurements. _Pfander G. (eds) Sampling Theory, a Renaissance. Applied and Numerical Harmonic Analysis_ , _Birkhäuser_.
  * [42] K. Ogura [1920], On a Certain Transcendental Function in the Theory of Interpolation. _Tôhoku Math. J._ , 17, 64-72.
  * [43] K. Grochenig and H. Razafinjatovo [1996], On Landau’s Necessary Density Conditions for Sampling and Interpolation of Band-Limited Functions, _Journal of the London Mathematical Society_ , 54(3), 557-565.
  * [44] K. Czyż [2004], Nonuniformly Sampled Active Noise Control System. _IFAC Proceedings Volumes_ , 37(20), 351-355.
  * [45] K. Han, Y. Wei and X. Ma [2016], An Efficient Non-uniform Filtering Method for Level-crossing Sampling. _IEEE International Conference on Digital Signal Processing_.
  * [46] L. Greengard and J. Lee [2004], Accelerating the Nonuniform fast Fourier transform. _Applied and Computational Harmonic Analysis_ , 35, 111-129.
  * [47] M. Ledoux and M. Talagrand [1991], Probability in Banach spaces. _Springer_.
  * [48] M. Gastpar and Y. Bresler [2000], On The Necessary Density for Spectrum-Blind Nonuniform Sampling Subject to Quantization. _Proc. IEEE Int. Conf. Acoustics, Speech, Signal Process._ , 1, 348-351.
  * [49] M. Rudelson and R. Vershynin [2007], On sparse reconstruction from Fourier and Gaussian measurements. _Communications on Pure and Applied Mathematics_ , 61(8), 1025-1045.
  * [50] M.A. Herman and T. Strohmer [2009]. High-resolution radar via compressed sensing. _IEEE Transactions on Signal Processing_ , 57(6), 2275-2284.
  * [51] M.W. Maciejewski, H.Z. Qui, M. Mobli and J.C. Hoch [2009], Nonuniform Sampling and Spectral Aliasing. _Journal of Magnetic Resonance_ , 199(1), 88-93.
  * [52] M. Pippig and D. Potts [2013], Parallel three-dimensional non-equispaced fast Fourier transforms and their applications to particle simulation. _SIAM Journal of Scientific Computing_ , 35(4), C411-C437.
  * [53] M. Jia, C. Wang, K. Ting Chen and T. Baba [2013], An Non-uniform Sampling Strategy for Physiological Signals Component Analysis. _Digest of Technical Papers - IEEE International Conference on Consumer Electronics_ , 526-529.
  * [54] M. Kabanava and H. Rauhut [2015]. Cosparsity in compressed sensing. _Boche H., Calderbank R., Kutyniok G., Vybíral J. (eds) Compressed Sensing and its Applications. Applied and Numerical Harmonic Analysis_. _Birkhäuser_.
  * [55] M. Kabanava, H. Rauhut and U. Terstiege [2015]. Analysis of low rank matrix recovery via Mendelson’s small ball method. _2015 International Conference on Sampling Theory and Applications (SampTA)_.
  * [56] M. Kabanava, R. Kueng, H. Rauhut and U. Terstiege [2016]. Stable low-rank matrix recovery via null space properties. _Information and Inference: A Journal of the IMA_ , 5(4), 405-441.
  * [57] M. Hajar, M. El Badaoui, A. Raad and F. Bonnardot [2019], Discrete Random Sampling: Theory and Practice in Machine Monitoring. _Mechanical Systems and Signal Processing_ , 123, 386-402.
  * [58] O. López, R. Kumar, Ö. Yılmaz and F.J. Herrmann [2016], Off-the-Grid Low-Rank Matrix Recovery and Seismic Data Reconstruction. _IEEE Journal of Selected Topics in Signal Processing_ , 10(4).
  * [59] O. López [2019], Embracing Nonuniform Samples (T), _University of British Columbia_. Retrieved from https://open.library.ubc.ca/collections/ubctheses/24/items/1.0380720
  * [60] P.W. Cary [1997], 3D Stacking of Irregularly Sampled Data by Wavefield Reconstruction. _SEG Technical Program Expanded Abstracts_.
  * [61] P.S. Penev and L.G. Iordanov [2001], Optimal Estimation of Subband Speech from Nonuniform Non-recurrent Signal-driven Sparse Samples. _IEEE International Conference on Acoustics, Speech and Signal Processing Proceedings_.
  * [62] P. Christensen, A. Kensler and C. Kilpatrick [2018], Progressive Multi-Jittered Sample Sequences. _Eurographics Symposium on Rendering_ , 37(4).
  * [63] R. L. Cook, T. Porter and L. Carpenter [1984], Distributed Ray Tracing. _ACM Siggraph 84 Conference Proceedings_ , 18(4), 165-174.
  * [64] R.L. Cook [1986], Stochastic Sampling in Computer Graphics. _ACM Transactions on Graphics_ , 6(1).
  * [65] R.D. Wisecup [1998], Unambiguous Signal Recovery Above the Nyquist Using Random‐sample‐interval Imaging. _Geophysics_ , 63(2), 331-789.
  * [66] R. Baraniuk, H. Choi, F. Fernandes, B. Hendricks, R. Neelamani, V. Ribeiro, J. Romberg, R. Gopinath, H. Guo, M. Lang, J.E. Odegard and D. Wei [2001], Rice Wavelet Toolbox. [Online]. https://www.ece.rice.edu/dsp/software/rwt.shtml. [Accessed: 13- June- 2019].
  * [67] R. Venkataramani and Y. Bresler [2001], Optimal Sub-Nyquist Nonuniform Sampling and Reconstruction for Multiband Signals. _IEEE Transactions on Signal Processing_ , 48(10), 2301 - 2313.
  * [68] R.I. Oliveira [2010], “Sums of Random Hermitian Matrices and an Inequality by Rudelson, _Elect. Comm. in Probability_ , 15, 203–212.
  * [69] R. Kumar, O. López, E. Esser and F.J. Herrmann [2015], Matrix completion on unstructured grids : 2-D seismic data regularization and interpolation. _EAGE Annual Conference Proceedings_.
  * [70] R. Kueng and P. Jung [2017]. Robust Nonnegative Sparse Recovery and the Nullspace Property of 0/1 Measurements. _IEEE Transactions on Information Theory_ , 64(2), 689-703.
  * [71] S. Maymon and A.V. Oppenheim [2011], Sinc Interpolation of Nonuniform Samples. _IEEE Transactions on Signal Processing_ , 59(10), 4745 - 4758.
  * [72] S. Foucart [2017], Flavors of Compressive Sensing. _Approximation Theory XV: San Antonio 2016_ , 201, 61 - 104.
  * [73] S. Foucart and H. Rauhut [2013], A Mathematical Introduction to Compressive Sensing. _Birkhäuser_.
  * [74] T. Strohmer [2000], Numerical Analysis of the Non-Uniform Sampling Problem. _Journal of Computational and Applied Mathematics_ , 122, 297-316.
  * [75] T. Klein and E. Rio [2005], Concentration around the mean for maxima of empirical processes, _Annals of Probability_ , 33(3), 1060–1077.
  * [76] T.T.Y. Lin and F.J, Herrmann [2009], Designing simultaneous acquisitions with compressive sensing. _EAGE Annual Conference Proceedings_.
  * [77] T.T. Cai and A. Zhang [2014], Sparse Representation of a Polytope and Recovery of Sparse Signals and Low-Rank Matrices. _IEEE Transactions on Information Theory_ , 60(1), 122-132.
  * [78] T. Wu, S. Dey, and Mike Shuo-Wei Chen [2016], A Nonuniform Sampling ADC Architecture With Reconfigurable Digital Anti-Aliasing Filter. _IEEE Journal of Selected Topics in Signal Processing_ , 63(10), 1639 - 1651.
  * [79] U. Hassan and Muhammad S. Anwar [2010], Reducing Noise by Repetition: Introduction to Signal Averaging, _European Journal of Physics_ , 31, 453-465.
  * [80] V.A. Kotel’nikov [1933], On the Transmission Capacity of “ether” and wire in electrocommunications, _(material for the first all-union conference on questions of communication) Izd. Red. Upr. Svyazi RKKA_.
  * [81] W.L. Ferrar [1928], On the Consistency of Cardinal Function Interpolation, _Proc. Roy. Soc. Edinburgh_ , 47, 230-242.
  * [82] Y.Y. Zeevi and E. Shlomot [1993], Nonuniform Sampling and Antialiasing in Image Representation. _IEEE Transactions on Signal Processing_ , 41(3), 1223-1236.
  * [83] [2004], An Introduction to Harmonic Analysis (Third ed.). _Cambridge University Press_.