# Recovery of surfaces and functions in high dimensions: sampling theory and
links to neural networks333 Parts of this work were presented at SIAM
conference on applied algebraic geometry (SIAM AG 2019), at the IEEE
International Symposium on Biomedical Imaging (ISBI 2020) and at the IEEE
International Conference on Acoustics, Speech, and Signal Processing (ICASSP
2020).

Qing Zou111Applied Mathematics and Computational Sciences, University of Iowa,
Iowa City, IA 52242<EMAIL_ADDRESS>and Mathews Jacob222Department of
Electrical and Computer Engineering, University of Iowa, Iowa City, IA 52242.
<EMAIL_ADDRESS>
Abstract. The learning of generative models for multi-dimensional signals is a
key problem in several machine learning tasks, with wide-ranging applications.
Similarly, the learning of functions of multidimensional datasets is a key
problem in several applications. Deep learning methods provide powerful tools
to address the above problems, but the theoretical underpinnings of these
methods are still emerging. In this paper, we focus on the above problems when
the signals follow a restrictive signal model. In particular, we consider
signals living on a union of smooth band-limited surfaces in high dimensions.
We show that an exponential mapping transforms the data to a union of low-
dimensional subspaces. Using this relation, we introduce a sampling
theoretical framework for the recovery of smooth surfaces from few samples and
the learning of functions living on smooth surfaces. The low-rank property of
the features is used to determine the number of measurements needed to recover
the surface. Moreover, the low-rank property of the features also provides an
efficient approach, which resembles a neural network, for the local
representation of multidimensional functions on the surface. The direct
representation of such a function in high dimensions often suffers from the
curse of dimensionality; the large number of parameters would translate to the
need for extensive training data. The low-rank property of the features can
significantly reduce the number of parameters, which makes the computational
structure attractive for learning and inference from limited labeled training
data.

Keywords. level set; surface recovery; function representation; image
denoising; neural networks

##  1 Introduction

The learning and exploitation of the extensive structure in real-world
multidimensional datasets is a key problem in machine learning. For example,
surface embedding methods model the data as points on smooth surfaces in high-
dimensional spaces for visualization [36]. Deep-learned generative models are
widely used to learn and exploit the structure in applications [15]. A popular
hypothesis in machine learning is the manifold assumption [10, 40], which
assumes that the data as points are lying close to a smooth low dimensional
manifold. Another key problem in machine learning is to learn functions or
mappings of high dimensional signals (e.g. denoising [48] and classification
[7]). The representations of these functions using conventional linear models
suffer from the curse of dimensionality. Specifically, the extensive number of
parameters translate to the need for large amounts of training data, which
makes these representations impractical in applications. While deep learning
methods offer fast and efficient alternatives, the impact of the structure of
the signals as well as the complexity of the specific function on the
performance of the learning process are relatively less understood.

We focus on the representation of points and functions on smooth band-limited
surfaces in this work. The main focus of this work is to introduce a sampling
theory for (a) recovering surfaces in high dimensional spaces, and (b) local
representation of functions that live on surfaces, from few measurements.
Popular approaches to represent surfaces include (a) explicit representations
using triangle mesh/graphs [6, 46], and (b) implicit level set representations
[35, 22]. In this work, we follow the latter approach because of its
flexibility in dealing with arbitrary topologies. We model the surface or
union of surfaces as the zero level set of a multidimensional band-limited
function $\psi$. A more band-limited function will translate to a smoother
surface and hence the support or bandwidth of $\psi$ can be viewed as a
measure of the complexity of the surface. We define an exponential lifting of
each point on the surface to obtain a feature vector; we show that the
exponential features of any point on the surface will be annihilated by the
inner product with the coefficients of $\psi$, which implies that the
exponential features live in a subspace. This property provides a bridge
between the complex union of surfaces representation and the well-studied
union of subspaces model. Specifically, the non-uniform mapping translates the
surface model to a subspace model. We use the low-rank property of the feature
matrix to estimate the surface from few of its samples. Sampling conditions
which will guarantee the recovery of the surface with probability 1, when the
surface is irreducible or when it is the union of several irreducible
surfaces, are introduced.

We also use the above property to efficiently represent and recover band-
limited multidimensional functions living on the surfaces. We show that such
representations are linear in the lifted feature space. Specifically, the low-
rank nature of the exponential features provides an elegant approach to
locally represent the function using considerably lower number of parameters.
In particular, we show that the feature vectors of a few anchor points on the
surface span the space, which allows us to efficiently represent the function
as the interpolation of the function values at the anchor points using a
Dirichlet kernel. The significant reduction in the number of free parameters
offered by this local representation makes the learning of the function from
finite samples tractable. We note that the computational structure of the
representation is essentially a one-layer kernel network. Note that the
approximation is highly local; the true function and the local representation
match only on the surface, while they may deviate significantly on points
which are not on the surface. This work is related to kernel methods, which
are widely used for the approximation of functions [4, 9, 24]. It is well-
known that an arbitrary function can be approximated using kernel methods, and
the computational structure resembles a single hidden layer neural network.
Our work has two key distinctions with the above approaches (a) unlike most
kernel methods that choose infinite bandwidth kernels (e.g. Gaussians), we
restrict our attention to a band-limited kernel. (b) we focus on a restrictive
data model, where the data samples are localized or close to the zero set of a
band-limited function. These two restrictions enable us to represent band-
limited functions on such smooth surfaces exactly, as opposed to classical
methods that can only approximate functions. The analysis shows the quite
significant benefit in representing complex functions, when the data have a
simple structure. In addition, it also shows how the composition of lower
bandwidth functions as in a multi-layer network can enable the efficient
representation of high bandwidth function.

This work builds upon our prior work [31, 30, 50], where we considered the
recovery of planar curves from several of the samples. This work extends the
above results in three important ways (i). The planar results are generalized
to the high dimensional setting. (ii). The worst-case sampling conditions are
replaced by high-probability results, which are far less conservative, and are
in good agreement with experimental results. (iii). The sampling results are
extended to the local representation of functions.

###  1.1 Terminology and Notation

We introduce some commonly used terminologies and notations throughout the
paper. We term the zero level set of a trigonometric polynomial as a surface.
Usually, the lower-case Greek letters $\psi,\eta,$ etc. are used to represent
the trigonometric polynomials. The calligraphic letter $\mathcal{S}$ or
$\mathcal{S}[\psi]$ is used to represent the zero level sets of the
trigonometric polynomials and hence the surfaces. The bold lower-case letters
$\mathbf{x}$ denotes the real variable in $[0,1)^{n}$ and sometimes the points
on the surface. The indexed bold lower-case letters $\mathbf{x}_{i}$ represent
the samples on the surface. The upper-case Greek letters
$\Lambda,\Gamma\subset\mathbb{Z}^{n}$ are used to denote the support of the
trigonometric polynomials. In other words, the upper-case Greek letters are
the coefficients index set. The coefficients set is shown as
$\\{\mathbf{c_{k}}:\mathbf{k}\in\Lambda\\}$. The cardinality of support
$\Lambda$ is given by $|\Lambda|$, which will serve as a measure of the
complexity of the surface. The notation $\Gamma\ominus\Lambda$ indicates the
set of all the possible shifts of the set $\Lambda$ within the set $\Gamma$.
The specific Greek letter $\Phi$ (sometimes subscripts are used to identify
the corresponding support) is used to represent the lifting map (feature map)
of the point on the surface. The notation $\Phi(\mathbf{X})$ denotes the
feature matrix of the sampling set $\mathbf{X}$.

##  2 Implicit surface representation

In this work, we use the level set representation to describe a
(hyper-)surface. We model a (hyper-)surface $\mathcal{S}$ in $[0,1)^{n};n\geq
2$ as the zero level set of a function $\psi$:

$\mathcal{S}[\psi]=\\{\mathbf{x}\in\mathbb{R}^{n}|\psi(\mathbf{x})=0\\}.$ (1)

For example, when $n=2$, $\mathcal{S}$ is a (hyper-)surface of dimension 1,
which is a curve. We note that the level set representation is widely used in
image segmentation [22]. The normal practice in image segmentation is the non-
parametric level set representation of a time-dependent evolution function
$\psi$.

###  2.1 Parametric level set representation of surfaces

Several authors have recently proposed to represent the level set function
$\psi$ as a linear combination of basis functions
$\varphi_{\mathbf{k}}(\mathbf{x})$ [44, 5]. These schemes argue that the
reduced number of parameters translate to fast and efficient algorithms.
Motivated by these schemes, we represent $\psi(\mathbf{x})$ as

$\psi(\mathbf{x})=\sum_{\mathbf{k}\in\Lambda}\mathbf{c_{k}}~{}\varphi_{\mathbf{k}}(\mathbf{x}).$
(2)

Since the level set function is the linear combination of some basis
functions, we term the corresponding zero level set as parametric zero level
set. We now focus on specific parametric representations, depending on the
choices of the basis functions.

Note that lower dimensional (hyper-)surfaces can be represented as the
intersection of multiple (hyper-)surfaces. For example, we can represent an
$n-2$ dimensional surface as the intersection of two $n-1$ dimensional
surfaces $\mathcal{S}[\psi_{1}]$ and $\mathcal{S}[\psi_{2}]$:

$\mathcal{S}[\psi_{1},\psi_{2}]=\\{\mathbf{x}\in\mathbb{R}^{n}|\psi_{1}(\mathbf{x})=\psi_{2}(\mathbf{x})=0\\}.$
(3)

###  2.2 Shift invariant surface representation

A popular choice for the basis functions is the shift invariant
representation, where compactly supported basis functions such as B-splines
are used. Specifically, the basis functions are shifted copies of a template
$\varphi$, denoted by:

$\varphi_{\mathbf{k}}(\mathbf{x})=\varphi\left(\frac{\mathbf{x}}{T}-\mathbf{k}\right)$
(4)

Here $T$ is the grid spacing, which controls the degrees of freedom of the
representation. The number of B-splines in the above representation is
$1/T^{n}-1$. One may also choose a multi-resolution or sparse wavelet surface
representation, when the basis functions are shifted and dilated copies of a
template. This approach allows the surface to have different smoothness
properties at different spatial regions.

###  2.3 band-limited surface representation

In this paper, we focus on the band-limited surface representation [27], where
the level set function $\psi$ is expressed as a Fourier series:

$\psi(\mathbf{x})=\sum_{\mathbf{k}\in\Lambda}\mathbf{c_{k}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x}),\quad\mathbf{x}\in[0,1)^{n}.$
(5)

This choice is motivated by the close connection of the above representation
with polynomials, as discussed below, as well as the existence of fast
algorithms using the fast Fourier transform. In the above representation, the
set $\Lambda$ denotes the support of the Fourier coefficients
$\mathbf{c}=\\{\mathbf{c_{k}}:\mathbf{k}\in\Lambda\\}$; its cardinality
$|\Lambda|$ is the number of free parameters in the surface representation.
The band-limited representation (5) is an important model in both mathematical
theories and engineering applications. It is commonly used in digital image
processing [42, 28, 50], biomedical image processing [43, 33, 27], and
geophysics [34].

As mentioned above, there is a close connection between the representation (5)
and polynomials [27]. We will make use of the one-to-one correspondence
between the representation (5) and the complex polynomials. We consider the
mappings $\nu_{i}:[0,1)\to\\{z\in\mathbb{C}:|z|=1\\}$:

$\nu_{i}(x_{i})=\exp(j2\pi x_{i})=:z_{i}$ (6)

to study the properties of the zero sets of (5). With this transformation, the
representation (5) simplifies to the complex polynomial denoted as
$\mathcal{P}[\psi]$, which is of the form

$\mathcal{P}[\varphi](\mathbf{z})=\sum_{\mathbf{k}\in\Lambda}c_{\mathbf{k}}\prod_{i=1}^{n}z_{i}^{k_{i}}.$
(7)

Since the mapping involves powers of $z_{i}$, where $z_{i}$ are specified by
the trigonometric mapping (6), we term the expansion in (5) as a trigonometric
polynomial.

####  2.3.1 Minimal/Irreducible band-limited surface

We note that the function $\psi(\mathbf{x})$ with support $\Lambda$ in (5) can
be expressed with a larger support $\Gamma\supset\Lambda$ in the Fourier
domain as

$\psi(\mathbf{x})=\sum_{\mathbf{k}\in\Gamma}\tilde{\mathbf{c}}_{\mathbf{k}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x}),\quad\mathbf{x}\in[0,1)^{n},$
(8)

where in this case the coefficients set $\tilde{\mathbf{c}}$ is the zero-
filled version of the vector $\mathbf{c}$, denoted by
$\tilde{\mathbf{c}}\in\mathbb{C}^{|\Gamma|}$:

$\tilde{\mathbf{c}}_{\mathbf{k}}=\begin{cases}\mathbf{c_{k}}&\text{if
}\mathbf{k}\in\Lambda\\\ 0&\text{else}\end{cases}.$ (9)

When the support $\Gamma$ is larger than the minimal one, we observe that the
representation of the surface is not unique. In particular, multiple level-set
functions of the form (8) can correspond to the same surface
$\mathcal{S}[\psi]$ denoted by (1) and (5). For example, the zero set of
$\psi$ will remain unchanged when multiplied by any function, whose magnitude
is non-zero everywhere. An example of special importance is

$\varphi_{1}=\varphi\cdot\exp(j2\pi\mathbf{k}_{0}^{T}\mathbf{x});\quad\mathbf{k}_{0}\in\Gamma\ominus\Lambda,$
(10)

whose zero sets are identical to ${\mathcal{S}}[\psi]$. Note that multiplying
$\psi$ by an exponential results in shifting its coefficients; any shift of
the coefficients in (9) such that the non-zero coefficients stay within
$\Gamma$ will yield the same surface. The notation $\Gamma\ominus\Lambda$ in
(10) denotes the set of shifts of $\mathbf{k}_{0}$ such that the shifted
support $\Lambda+\mathbf{k}_{0}$ is still within $\Gamma$. See Fig. 3 for an
illustration. We note that the shifting process translates to multiplying the
polynomial in (7) by a monomial, resulting in an increased degree. The
following result shows that for each $\mathcal{S}$, there exists a unique
$\psi$ satisfying (1), whose coefficient set
$\\{\mathbf{c_{k}}:\mathbf{k}\in\Lambda\\}$ has the smallest support, which we
term it as the _minimal trigonometric polynomial_.

######  Proposition 1.

For every (hyper-)surface $\mathcal{S}$ given by the zero level set of (8),
there is a unique (up to scaling) minimal trigonometric polynomial $\psi$,
which satisfies $\psi(\mathbf{x})=0;\forall\mathbf{x}\in\mathcal{S}$. Any
other trigonometric polynomial $\psi_{1}$ that also satisfies
$\psi_{1}(\mathbf{x})=0;\forall\mathbf{x}\in\mathcal{S}$ will have
$supp(\psi_{1})\supseteq supp(\psi)$. Here, $supp(\psi)$ denotes the support
of the Fourier coefficients of $\psi$.

The proof of this result is given in Appendix 8.1.

We refer to the $\psi$ of the form (5) with the minimal support $\Lambda$ that
satisfy

$\psi(\mathbf{x})=0;~{}~{}\forall\mathbf{x}\in\mathcal{S}$ (11)

as the _minimal trigonometric polynomial_ of the surface $\mathcal{S}$. In
other words, when $\psi$ is the minimal trigonometric polynomial of a surface
$\mathcal{S}$, it does not have a factor that never vanishes on $[0,1)^{n}$ or
only vanishes at some isolated points. Specifically, if a polynomial has a
factor that never vanishes on $[0,1)^{n}$, one can remove this factor and
obtain a polynomial with a smaller support. Note from (7) that the minimal
trigonometric polynomial will correspond to $\mathcal{P}[\psi]$ being a
polynomial with the minimal degree.

We note that the support $\Lambda$ of the minimal polynomial of the surface
$\mathcal{S}$ grows with the complexity of $\mathcal{S}$; a more oscillatory
surface with a lot of details correspond to a high support minimal polynomial,
while a simple and highly smooth surface correspond to a low support minimal
polynomial. We hence consider $|\Lambda|$ as a _complexity measure_ of the
surface. Furthermore, we note that the surface model can approximate an
arbitrary closed surface with any degree of accuracy, as long as the support
is large enough [27]. See Fig. 1 for an illustration in 3D.

(a) $17\times 17\times 17$ coefficients (b) $25\times 25\times 25$
coefficients

(c) $33\times 33\times 33$ coefficients

Figure 1: Illustration the fertility of our level set representation model in
3D. The three examples show that our model is capable to capture the geometry
of the shape even though the shape has complicated topologies, which
demonstrated that the representation is not restrictive.

We note that the mapping $\nu=(\nu_{1},\cdots,\nu_{n})$ defined by (6) is a
bijection from $[0,1)^{n}$ onto the complex unit torus
$\mathbb{CT}^{n}=\\{(z_{1},\cdots,z_{n}):|z_{i}|=1,i=1,\cdots,n\\}$. Hence,

$\psi(\mathbf{x})=0~{}~{}\Leftrightarrow~{}~{}\mathcal{P}[\psi][\mathbf{z}]=0~{}\text{on}~{}\mathbb{CT}^{n},~{}~{}\mbox{where}~{}~{}z_{i}=\nu_{i}(x_{i}),~{}~{}i=1,\cdots,n,$
(12)

which implies that there is a one-to-one correspondence between the zero sets
of $\psi$ and the zeros of $\mathcal{P}[\psi]$ on the unit torus. Accordingly,
we can study the algebraic properties of trigonometric polynomials and their
zero sets by studying their corresponding complex polynomials under the
mapping $\nu$. Following the steps of [27], we can then define the
irreducibility of a trigonometric polynomial.

######  Definition 2.

A trigonometric polynomial $\psi(\mathbf{x})$ is said to be irreducible, if
the corresponding polynomials $\mathcal{P}[\psi]$ is irreducible as a
polynomial in $\mathbb{C}[z_{1},\cdots,z_{n}]$. A polynomial $p$ is
irreducible over a field of complex numbers, if it cannot be expressed as the
product of two or more non-constant polynomials with complex coefficients.

######  Definition 3.

A surface is termed as irreducible, if it is the zero set of an irreducible
trigonometric polynomial.

When $\psi$ can be written as the product of several irreducible components
$\psi=\prod_{i=1}^{m}\psi_{i}$, then $\mathcal{S}[\psi]$ is essentially the
union of irreducible surfaces:

$\mathcal{S}[\psi]=\bigcup_{i=1}^{m}\mathcal{S}[\psi_{i}].$ (13)

##  3 Lifting and feature spaces

We now use the parametric level set representation of surfaces to show that
the points on these surfaces satisfy certain annihilation properties. Consider
the non-linear lifting mapping
$\Phi_{\Gamma}:\mathbb{R}^{n}\rightarrow\mathbb{C}^{|\Gamma|}$:

$\Phi_{\Gamma}(\mathbf{x})=\begin{bmatrix}\varphi_{\mathbf{k}_{1}}(\mathbf{x})\\\
\vdots\\\ \quad\varphi_{\mathbf{k}_{|\Gamma|}}(\mathbf{x}))\end{bmatrix}$ (14)

of each point in the original domain $[0,1]^{n}$. We can view
$\Phi_{\Gamma}(\mathbf{x})$ as the feature vector of the point $\mathbf{x}$,
analogous to the ones used in kernel methods [38]. Here, $|\Gamma|$ denotes
the cardinality of the set $\Gamma$. We denote the set

$\mathcal{V}_{\Gamma}(\mathcal{S})=\\{\Phi_{\Gamma}(\mathbf{x})|\mathbf{x}\in\mathcal{S}\\}$
(15)

as the feature space of the surface $\mathcal{S}$. Since any point on a
surface $\mathcal{S}$ satisfies (1), the feature vectors of points from
$\mathcal{S}$ satisfy

${\mathbf{c}}^{T}\Phi_{\Gamma}(\mathbf{x})=0,~{}~{}\forall\mathbf{x}\in\mathcal{S},$
(16)

where $\mathbf{c}$ is the coefficients vector in the representation of $\psi$
in (5). The above relation implies that $\mathbf{c}$ is orthogonal to all the
feature vectors of points living on $\mathcal{S}$.

Figure 2: Illustration of the annihilation relations (16) in 2D. We assume
that the curve is the zero level set of a band-limited function
$\psi(\mathbf{x})$, shown in the top left. The Fourier coefficients of $\psi$,
denoted by $\mathbf{c}$, are support limited in $\Lambda$, denoted by the red
square in the figure in the bottom right. Each point on the curve satisfies
$\psi(\mathbf{x}_{i})=0$. Using the representation of the curve specified by
(1), we thus have $\mathbf{c}^{T}\Phi_{\Lambda}(\mathbf{x}_{i})=0$. Note that
$\Phi_{\Lambda}(\mathbf{x}_{i})$ is the exponential feature map of the point
$\mathbf{x}_{i}$, whose dimension is specified by the cardinality of the set
$\Lambda$. This means that the feature map will lift each point in the level
set to a $|\Lambda|$ dimensional subspace whose normal vector is specified by
$\mathbf{c}$, as illustrated by the plane and the red vector $\mathbf{c}$ in
the top right. Note that if more than one closed curve are presented, each
curve will be lifted to a lower dimensional subspace in the feature space, as
shown by the two lines in the plane, and the lower dimensional spaces will
span the $|\Lambda|$ dimensional subspace. (Figure courtesy of Q. Zou, reprint
from [50] with permission from IEEE).

###  3.1 Shift invariant representation

We consider
$\varphi_{\mathbf{k}}=\beta^{p}\left(\frac{\mathbf{x}}{T}-\mathbf{k}\right)$
to be the $p^{\rm th}$-degree tensor-product B-spline function. Note that
$\beta^{p}(\mathbf{x})$ is support limited in $[-(p+1)/2,(p+1)/2]^{n}$. If the
support of $\varphi_{k}(\mathbf{x})$ does not overlap with $\mathcal{S}$, we
have $\varphi_{k}(\mathbf{x})=0;\forall\mathbf{x}\in\mathcal{S}$. Hence, we
have

$\mathbf{i}_{\mathbf{k}}^{T}~{}\Phi_{\Gamma}(\mathbf{x})=0,~{}~{}\forall\mathbf{x}\in\mathcal{S}$
(17)

where $\mathbf{i}_{\mathbf{k}}$ is the indicator vector. Note that all of
these indicator vectors are linearly independent. The number of basis vectors
whose support does not overlap with $\mathcal{S}$ is dependent on the area of
$\mathcal{S}$ as well as the support of $\psi$. Thus, the dimension of
$\mathcal{V}_{\Gamma}(\mathcal{S})$ is a measure of the area of the surface
$\mathcal{S}$, and satisfies

${\rm{dim}}\left(\mathcal{V}_{\Gamma}(\mathcal{S})\right)\leq|\Gamma|-(P+1)=A,$
(18)

where $P$ is the number of basis functions whose support does not overlap with
$\mathcal{S}$.

###  3.2 band-limited surface representation

We now consider the case of an arbitrary point $\mathbf{x}$ on the zero level
set of $\psi(\mathbf{x})$. Using (8), the lifting is specified by:

$\Phi_{\Gamma}(\mathbf{x})=\begin{bmatrix}\exp(j2\pi\mathbf{k}_{1}^{T}\mathbf{x})\\\
\exp(j2\pi\mathbf{k}_{2}^{T}\mathbf{x})\\\ \vdots\\\
\quad\exp(j2\pi\mathbf{k}_{|\Gamma|}^{T}\mathbf{x})\end{bmatrix}.$ (19)

We note from (19) that the lifting $\Phi_{\Gamma}$ can be evaluated with
different supports $\Gamma\supset\Lambda$. When the set $\Gamma$ corresponds
to the support of the minimal polynomial ($\Gamma=\Lambda$), we term the
corresponding lifting as the _minimal lifting_. We now show that the feature
space $\mathcal{V}_{\Lambda}(\mathcal{S})$ is a lower dimensional space, or
equivalently contained in a subspace of $\mathbb{C}^{|\Lambda|}$.

####  3.2.1 Irreducible surface with minimal lifting ($\Gamma=\Lambda$)

We first focus on the case where $\psi$ is an irreducible trigonometric
polynomial and the support of the lifting is specified by $\Lambda$, which is
the support of the minimal polynomial. The annihilation relation (16) implies
that $\mathbf{c}$ is orthogonal to the feature vectors
$\Phi_{\Lambda}(\mathbf{x})$. This implies that

${\dim}(\mathcal{V}_{\Lambda})\leq|\Lambda|-1.$ (20)

####  3.2.2 Irreducible surface with non-minimal lifting
($\Gamma\supset\Lambda$)

We now consider the setting where the non-linear lifting is specified by
$\Phi_{\Gamma}(\mathbf{x})$, where $\Gamma\supset\Lambda$. Because of the
annihilation relation, we have

$\tilde{\mathbf{c}}^{T}~{}\Phi_{\Gamma}(\mathbf{x})=0,$

where $\tilde{\mathbf{c}}$ is the zero filled coefficients in (8). Since the
zero set of the function
$\psi_{\mathbf{k}_{0}}(\mathbf{x})=\psi(\mathbf{x})\cdot\exp(j2\pi\mathbf{k}_{0}^{T}\mathbf{x})$
is exactly the same as that of $\psi$, we have

$\sum_{\mathbf{k}}\mathbf{c_{\mathbf{k}-\mathbf{k}_{0}}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x})=0;~{}~{}\forall\mathbf{x}\in\mathcal{S}[\psi].$
(21)

This implies that any shift of $\tilde{\mathbf{c}}$ within
$\Gamma\ominus\Lambda$, denoted by
$\tilde{\mathbf{d}}_{\mathbf{k}}=\mathbf{c}_{\mathbf{k}-\mathbf{k}_{0}}$ will
satisfy $\tilde{\mathbf{d}}^{T}~{}\Phi_{\Gamma}(\mathbf{x})=0.$ It is
straightforward to see that $\tilde{\mathbf{d}}$ and $\tilde{\mathbf{c}}$ are
linearly independent for all values of $\mathbf{k}_{0}$. We denote the number
of possible shifts such that the shifted set $\Lambda+\mathbf{k}_{0}$ is still
within $\Gamma$ (i.e, $\Lambda+\mathbf{k}_{0}\subseteq\Gamma$ ) by
$|\Gamma\ominus\Lambda|$:

$\Gamma\ominus\Lambda=\\{\mathbf{l}\in\Gamma\mid\mathbf{l}-\mathbf{k}\in\Gamma,\forall\mathbf{k}\in\Lambda\\}.$
(22)

This set is illustrated in Fig. 3 along with $\Gamma$ and $\Lambda$.

Figure 3: The non-minimal filter support $\Gamma$ (green) is illustrated along
with the minimal filter support $\Lambda$ (red). The set
$\Gamma\ominus\Lambda$ (blue) contains all indices at which $\Lambda$ can be
centered, while remaining inside $\Gamma$.

Since the vectors $\mathbf{c}_{\mathbf{k}-\mathbf{k}_{0}}$ are linearly
independent and are orthogonal to any feature vector
$\Phi_{\Gamma}(\mathbf{x})$ on $\mathcal{S}[\psi]$, the dimension of the
subspace is bounded by

${\dim}(\mathcal{V}_{\Gamma})\leq|\Gamma|-|\Gamma\ominus\Lambda|.$ (23)

####  3.2.3 Union of irreducible surfaces with $\Gamma\supset\Lambda_{i}$

When $\psi=\prod_{i=1}^{m}\psi_{i}$, each irreducible surface
$\mathcal{S}[\psi_{i}]$ will be mapped to a subspace of dimension
$|\Gamma|-|\Gamma\ominus\Lambda_{i}|$. This implies that the non-linear
lifting transforms the union of irreducible surfaces to the well-studied union
of subspace model [23, 14, 25].

##  4 band-limited surface recovery from samples

We now consider the recovery of the surface $\mathcal{S}$ from its samples
$\mathbf{x}_{i};i=1,\cdots,N$. As discussed in the introduction, this is an
important problem in the recovery of shapes from point clouds as well as in
segmentation problems. According to the analysis in the previous section, if
the sampling point $\mathbf{x}_{i}$ is located on the zero level set of
$\psi(\mathbf{x})$, we will then have the annihilation relation specified by
(16). Notice that equation (16) is a linear equation with $\mathbf{c}$ as its
unknowns. Since all the samples $\mathbf{x}_{i};i=1,..,N$ satisfy the
annihilation relation (16), we have

$\mathbf{c}^{T}\underbrace{\left[\Phi_{\Gamma}(\mathbf{x}_{1})\quad\cdots\quad\Phi_{\Gamma}(\mathbf{x}_{N})\right]}_{\Phi_{\Gamma}(\mathbf{X})}=0.$
(24)

We call $\Phi_{\Gamma}(\mathbf{X})$ the feature matrix of the sampling set
$\mathbf{X}=\\{\mathbf{x}_{1},\cdots,\mathbf{x}_{N}\\}$. We propose to
estimate the coefficients $\mathbf{c}$, and hence the surface
$\mathcal{S}[\psi]$ using the above linear relation (24). Note that
$\mathcal{S}[\psi]$ is invariant to the scale of $\mathbf{c}$; without loss of
generality, we reformulate the estimation of the surface as the solution to
the system of equations

$\mathbf{c}^{T}~{}\bm{\Phi}_{\Gamma}(\mathbf{X})=0;~{}~{}\|\mathbf{c}\|_{F}=1.$
(25)

###  4.1 Sampling theorems

The above estimation scheme yields a unique solution, if the matrix
$\Phi_{\Lambda}(\mathbf{X})$ has a unique null-space basis vector. We will now
focus on the number of samples $N$ and its distribution on
$\mathcal{S}[\psi]$, which will guarantee the unique recovery of
$\mathcal{S}[\psi]$. We will consider different lifting scenarios introduced
in Section 3 separately. As we will see, in some cases considered below, the
null-space has a large dimension. However, the minimal null-space vector
(coefficients with the minimal support) will still uniquely identify the
surface, provided the sampling conditions are satisfied.

####  4.1.1 Irreducible surfaces with minimal lifting

Suppose $\psi(\mathbf{x})$ is an irreducible trigonometric polynomial with
support $\Lambda$. Consider the lifting which is specified by the minimal
support $\Lambda$. We see from (20) that
$\rm{rank}\left(\Phi_{\Lambda}(\mathbf{X})\right)\leq|\Lambda|-1$. The
following result shows when the inequality is replaced by an equality.

######  Proposition 4.

Let $\\{\mathbf{x}_{1},\cdots,\mathbf{x}_{N}\\}$ be $N$ independent and random
samples on the surface $\mathcal{S}[\psi]$, where $\psi(\mathbf{x})$ is an
irreducible (minimal) trigonometric polynomial with support $\Lambda$. The
feature matrix $\Phi_{\Lambda}(\mathbf{X})$ will have rank $|\Lambda|-1$ with
probability $1$, if

$N\geq|\Lambda|-1.$

The above proposition guarantees that the solution to the system of equations
specified by (25) is unique (up to scaling) when the number of samples exceed
$N=|\Lambda|-1$ with unit probability. The proof of this proposition can be
found in Appendix 8.2. With Proposition 4, we obtain the following sampling
theorem.

######  Theorem 5 (Irreducible surfaces of any dimension).

Let $\psi(\mathbf{x}),\mathbf{x}\in[0,1]^{n},n\geq 2$ be an irreducible
trigonometric polynomial whose support is given by $\Lambda$. The zero level
set of $\psi(\mathbf{x})$ is denoted as $\mathcal{S}[\psi]$. If we are
randomly given $N\geq|\Lambda|-1$ samples on $\mathcal{S}[\psi]$, then with
probability $1$ that we can uniquely recover the surface $\mathcal{S}[\psi]$.

The above theorem generalizes the results in [50] to any dimension $n\geq 2$.
When $n=2$, then $\mathcal{S}$ is a planar curve. In this setting, if the
support of $\psi$ is given by $\Lambda=k_{1}\times k_{2}$. Then by this
sampling theorem, we get perfect recovery with probability one, when the
number of random samples on the curve exceeds $k_{1}\cdot k_{2}-1$. Note that
the degrees of freedom in the representation (5) is $k_{1}\cdot k_{2}-1$, when
we constrain $\|\mathbf{c}\|_{F}=1$. This implies that if the number of
samples exceed the degrees of freedom, we get perfect recovery. Note that
these results are significantly less conservative than the ones in [50], which
required a minimum of $(k_{1}+k_{2})^{2}$ samples. We note that the results in
[50] were the worst case guarantees, and will guarantee the recovery of the
curve from any $(k_{1}+k_{2})^{2}$ samples. By contrasts, our current results
are high probability results; there may exist a set of $N\geq k_{1}\cdot
k_{2}-1$ samples from which we cannot get unique recovery.

We note that the current work is motivated by the phase transition experiments
(Fig. 5) in [50], which shows that one can recover the curve in most cases
when the number of samples exceeds $k_{1}\cdot k_{2}-1$ rather than the
conservative bound of $(k_{1}+k_{2})^{2}$. We also note that it is not
straightforward to extend the proof in [50] to the cases beyond $n=2$.
Specifically, we relied on Bezout’s inequality in [50], which does not
generalize easily to high dimensional cases.

(a) Curve (b) Recovery #1 (c) Recovery #2

(d) original $\psi(x,y)$ (e) $\psi(x,y)$ with 7 samples (f) $\psi(x,y)$ with 8
samples

Figure 4: Illustration of Theorem 5 in 2D. The irreducible curve given by (a)
is the original curve, which is obtained by the zero level set of a
trigonometric polynomial whose support is $3\times 3$. According to Theorem 5,
we will need at least 8 samples to recover the curve. In (b), we randomly
choose 7 samples (the red dots) on the original curve (the gray curve). The
blue dashed curve shows the recovered curve from this 7 samples. Since the
sampling condition is not satisfied, the recovery failed. In (c), we randomly
choose 8 points (the red dots). From (c), we see that the blue dashed curve
(recovered curve) overlaps the gray curve (the original curve), meaning that
we recover the curve perfectly. In (d) - (f), we showed the original
trigonometric polynomial, the polynomial obtained from 7 samples and the
polynomial obtained from 8 samples.

(a) Surface (b) Recovery #1

(c) Recovery #2

Figure 5: Illustration of Theorem 5 in 3D. The irreducible surface given by
(a) is the original surface, which is given by the zero level set of a
trigonometric polynomial whose support is $3\times 3\times 3$. According to
Theorem 5, we will need at least 26 samples to recover the surface. In (b), we
randomly choose 25 samples (the blue dots) on the original surface (the gray
part). The red surface is what we recovered from the 25 samples. Since the
sampling condition is not satisfied, the recovery failed. In (c), we randomly
choose 26 points (the blue dots). From (c), we see that the red surface
(recovered surface) overlaps the gray surface (the original surface), meaning
that we recover the surface perfectly.

####  4.1.2 Union of irreducible surfaces with minimal lifting

We now consider the union of irreducible surfaces $\mathcal{S}[\psi]$, where
$\psi$ has several irreducible factors
$\psi(\mathbf{x})=\psi_{1}(\mathbf{x})\cdots\psi_{M}(\mathbf{x})$. Then we
have $\mathcal{S}[\psi]=\bigcup_{i=1}^{M}\mathcal{S}[\psi_{i}]$. Suppose the
support of $\psi(\mathbf{x})$ is given by $\Lambda$ and the support of each
factor $\psi_{i}(\mathbf{x})$ is given by $\Lambda_{i}$. We have the following
result for this setting.

######  Proposition 6.

Let $\psi(\mathbf{x})$ be a trigonometric polynomial with $M$ irreducible
factors, i.e.,

$\psi(\mathbf{x})=\psi_{1}(\mathbf{x})\cdots\psi_{M}(\mathbf{x}).$ (26)

Suppose the support of each factor $\psi_{i}(\mathbf{x})$ is given by
$\Lambda_{i}$ and the support of $\psi$ is $\Lambda$. Assume that
$\\{\mathbf{x}_{1},\cdots,\mathbf{x}_{N}\\}$ are $N$ random samples on
$\mathcal{S}[\psi]$ that are chosen independently. Then with probability 1
that the feature matrix $\Phi_{\Lambda}(\mathbf{X})$ will be of rank
$|\Lambda|-1$ if

  1. 1.

each irreducible factor is sampled with $N_{i}\geq|\Lambda_{i}|-1$ points, and

  2. 2.

the total number of samples satisfy $N\geq|\Lambda|-1$.

The proof of this result can be seen in Appendix 8.2.3. Based on this
proposition, we have the following sampling conditions.

######  Theorem 7 (Union of irreducible surfaces of any dimension).

Let $\psi(\mathbf{x})$ be a trigonometric polynomial with $M$ irreducible
factors as in (27). If the samples $\mathbf{x}_{1},.,\mathbf{x}_{N}$ satisfy
the conditions in Proposition 6, then the surface can be uniquely recovered by
the solution of (25) with probability 1.

Unlike the sampling conditions in Theorem 5 that does not impose any
constraints on the sampling, the above result requires each component to be
sampled with a minimum rate specified by the degrees of freedom of that
component. We illustrate the above result in Fig. 6 in 2D ($n=2$), where
$\mathcal{S}$ is the union of two irreducible curves with support of $3\times
3$, respectively. The above results show that if each of the curves are
sampled with at least eight points and if the total number of samples is no
less than 24, we can uniquely identify the union of curves. The results show
that if any of the above conditions are violated, the recovery fails; by
contrast, when the number of randomly chosen points satisfy the conditions, we
obtain perfect recovery.

(a) Curve (b) $7+17$

(c) $8+16$

(d) $17+7$ (e) $16+8$

(f) $8+8$

Figure 6: Illustration of Theorem 7. The original curve (a) is given by the
zero set of a reducible trigonometric polynomial with support $5\times 5$,
which is the product of two trigonometric polynomials with support $3\times
3$. According to the sampling theorem, we totally need at least 24 samples and
each of the component need to be sampled for at least 8 samples. We first
choose 7 samples (red dots) on the first component and 17 samples (red
circles) on the second one. The gray curve in (b) is the original curve and
the blue dashed curve is what we recovered from the $7+17=24$ samples. Since
the sampling condition is not satisfied, the recovery failed. In (c), we
choose 8 samples (red dots) on the first component and 16 samples (red
circles) on the second one. From (c), we see that the gray curve (the original
curve) overlaps the blue dashed curve (recovered curve), meaning that we
recovered the curve successfully. In (d), we choose 17 samples on the first
component and 7 samples on the other one. From (d), we see that the recovery
is not successful. In (e), we have 16 samples on the first component and 8
samples on the second one. The original curve overlaps the recovered one. So
we recovered it perfectly. Lastly, we choose 8 samples on each of the
component and we failed to recover the curve as shown in (f). Note that the
recovered curves pass through the samples in all cases.

####  4.1.3 Non-minimal lifting

In Section 4.1.1, we introduced theoretical guarantees for the perfect
recovery of the surface in any dimensions. The sampling theorems introduced in
Section 4.1.1 and 4.1.2 assume that we know exactly the support of the surface
or the union of surfaces. However, in practice, the true support of the
surface is usually unknown. We now consider the recovery of the surface, when
the support is over-estimated, or equivalently the lifting is performed
assuming $\Gamma\supset\Lambda$. As discussed in Section 3.2.2, the dimension
of $\mathcal{V}_{\Gamma}$ is upper bounded by
$|\Gamma|-|\Gamma\ominus\Lambda|$, which implies that

$\operatorname{rank}(\Phi_{\Gamma}(\mathbf{X}))\leq|\Gamma|-|\Gamma\ominus\Lambda|,$

where $\Gamma\ominus\Lambda$ represents the number of valid shifts of
$\Lambda$ within $\Gamma$ as discussed in Section 3.2.2.

The following two propositions show when the inequality in the rank relation
above can be an equality and hence we can recover the surface.

######  Proposition 8 (Irreducible surface with non-minimal lifting).

Let $\\{\mathbf{x}_{1},\cdots,\mathbf{x}_{N}\\}$ be $N$ random samples on the
surface $\mathcal{S}[\psi]$, chosen independently. The trigonometric
polynomial $\psi(\mathbf{x})$ is irreducible whose true support is $\Lambda$.
Suppose the lifting mapping is performed using support
$\Gamma\supset\Lambda$.Then
$\operatorname{rank}(\Phi_{\Gamma}(\mathbf{X}))=|\Gamma|-|\Gamma\ominus\Lambda|$
with probability 1, if

$N\geq|\Gamma|-|\Gamma\ominus\Lambda|.$

The proof of this proposition can be found in Appendix 8.2.4.

######  Proposition 9 (Union of irreducible surfaces with non-minimal
lifting).

Let $\psi(\mathbf{x})$ be a trigonometric polynomial with $M$ irreducible
factors, i.e.,

$\psi(\mathbf{x})=\psi_{1}(\mathbf{x})\cdots\psi_{M}(\mathbf{x}).$ (27)

Suppose the support of each factor $\psi_{i}(\mathbf{x})$ is given by
$\Lambda_{i}$ and the support of $\psi$ is $\Lambda$. Let
$\Gamma_{i}\supset\Lambda_{i}$ be the non-minimal support of each factor
$\psi_{i}(\mathbf{x})$ and $\Gamma\supset\Lambda$ is the support of the non-
minimal lifting. Assume $\\{\mathbf{x}_{1},\cdots,\mathbf{x}_{N}\\}$ are $N$
random samples on $\mathcal{S}[\psi]$ that are chosen independently. Then with
probability 1 that the feature matrix $\Phi_{\Lambda}(\mathbf{X})$ will be of
rank $|\Gamma|-|\Gamma\ominus\Lambda|$ if

  1. 1.

each irreducible factor is sampled with
$N_{i}\geq|\Gamma_{i}|-|\Gamma_{i}\ominus\Lambda_{i}|$ points, and

  2. 2.

the total number of samples satisfy $N\geq|\Gamma|-|\Gamma\ominus\Lambda|$.

We prove this result in Appendix 8.2.5. Note that in practice, when non-
minimal lifting mapping is performed, we then randomly sample approximately
$|\Gamma|-|\Gamma\ominus\Lambda|$ positions on $\mathcal{S}$. This random
strategy ensures that the samples are distributed to the factors, roughly
satisfying the conditions in Proposition 9.

####  4.1.4 Surface recovery algorithm for the non-minimal setting

The two propositions in Section 4.1.3 show that $\Phi_{\Gamma}(\mathbf{X})$
has $|\Gamma\ominus\Lambda|$ null space vectors
$\mathbf{n}_{i}\leftrightarrow\mu_{i}$, when the non-minimal lifting with
support $\Gamma$ is performed. The following result from [50] shows the
property of these null space vectors.

######  Proposition 10 (Proposition 9 in [50]).

The coefficients of the trigonometric polynomials of the form

$\theta_{\mathbf{k}}(\mathbf{x})=\exp(j2\pi\mathbf{l}^{T}\mathbf{x})\psi(\mathbf{x}),\quad\forall\mathbf{k}\in\Gamma\ominus\Lambda.$

is a null space vector of $\Phi_{\Gamma}(\mathbf{X})$.

Note that the coefficients of $\theta_{{\mathbf{k}}}(\mathbf{x})$ correspond
to the shifted versions of the coefficients of $\psi$ and hence are linearly
independent. We also note that any such function is a valid annihilating
functions for points on $\mathcal{S}$. When the dimension of the null-space is
$\Gamma\ominus\Lambda$, these corresponding coefficients form a basis for the
null-space. Therefore, we have that any function in the null-space can be
expressed as

$\displaystyle\eta(\mathbf{x})$ $\displaystyle=$
$\displaystyle\sum_{\mathbf{k}\in\Gamma\ominus\Lambda}\alpha_{\mathbf{k}}~{}\psi(\mathbf{x})\exp(j2\pi\mathbf{k}^{T}\mathbf{x})$
(28) $\displaystyle=$
$\displaystyle\psi(\mathbf{x})\underbrace{\sum_{\mathbf{k}\in\Gamma\ominus\Lambda}\alpha_{\mathbf{k}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x})}_{\gamma(\mathbf{x})}=\psi(\mathbf{x})\gamma(\mathbf{x}),$
(29)

where $\alpha_{\mathbf{k}}$ and $\gamma$ are arbitrary coefficients and
function, respectively. Note that all of the functions obtained by the null-
space vectors have $\psi$ as a common factor.

Accordingly, we have that $\psi(\mathbf{x})$ is the greatest common divisor of
the polynomials $\mu_{i}(\mathbf{x})\leftrightarrow\mathbf{n}_{i}$, where
$\mathbf{n}_{i}$ are the null-space vectors of $\Phi_{\Gamma}(\mathbf{X})$,
which can be estimated using singular value decomposition (SVD). Since we
consider polynomials of several variables, it is not computationally efficient
to find the greatest common divisor. We note that we are not interested in
recovering the minimal polynomial, but are only interested in finding the
common zeros of $\mu_{i}(\mathbf{x})$. We hence propose to recover the
original surface as the zeros of the the sum of squares polynomial

$\gamma(\mathbf{x})=\sum_{i=1}^{|\Gamma\ominus\Lambda|}|\mu_{i}(\mathbf{x})|^{2}.$

Note that rank guarantees in Propositions 8.2.4 and 8.2.5 ensure that the
entire null-space will be fully identified by the feature matrix. Coupled with
Proposition 10, we can conclude that the recovery using the above algorithm
(SVD, followed by the sum of squares of the inverse Fourier transforms of the
coefficients) will give perfect recovery of the surface under noiseless
conditions.

(a) The original curve (b) The sampling points

(c) 1st null-space function (d) 2nd null-space function

(e) sum-of-squares polynomial

Figure 7: Illustration of the sampling fashion for non-minimal Fourier
support. We consider the curve as shown in (a), which is given by the zero
level set of a trigonometric polynomial of support $5\times 5$. We choose the
non-minimal support $\Gamma$ as $11\times 11$. According to the sampling
condition for non-minimal Fourier support, we sampled on 72 random locations.
We randomly chose two null-space vectors for the feature matrix of the
sampling set, which gave us functions (c) and (d). We can see that all of
these functions have zeros on the original zero set, in addition to processing
several other zeros. The sum of squares function is shown in (e), showing the
common zeros, which specifies the original curve.

###  4.2 Surface recovery from noisy samples

The analysis in Section 4.1.3 shows that when the support of the surface is
small, the feature matrix is low rank. In practice, the sampling points are
usually corrupted with some noise. We denote the noisy sampling set by
$\mathbf{Y}=\mathbf{X}+\mathbf{N}$, where $\mathbf{N}$ is the noise. We
propose to exploit the low-rank nature of the feature matrix to recover it
from noisy measurements. Specifically, when the sampling set $\mathbf{X}$ is
corrupted by noise, the points will deviate from the original surface, and
hence the features will cease to be low rank.

We pose the denoising as:

$\mathbf{X}^{*}=\arg\min_{\mathbf{X}}||\mathbf{X}-\mathbf{Y}||^{2}+\lambda||\Phi(\mathbf{X})||_{*}$
(30)

where we use the nuclear norm of the feature matrix of the sampling set as a
regularizer. Unlike traditional convex nuclear norm formulations, the above
scheme is non-convex.

We adapt the kernel low-rank algorithm in [31] to the high dimensional setting
to solve (30). This algorithm relies on an iteratively reweighted least
squares (IRLS) approach [11, 26] which alternates between the following two
steps:

$\mathbf{X}^{(n)}=\arg\min_{\mathbf{X}}||\mathbf{X}-\mathbf{Y}||^{2}+\lambda{\mathrm{trace}}[\mathcal{K}(\mathbf{X})\mathbf{P}^{(n-1)}],$
(31)

and

$\mathbf{P}^{(n)}=\left[\mathcal{K}(\mathbf{X}^{(n)})+\gamma^{(n)}\mathbf{I}\right]^{-1/2}$
(32)

where $\gamma^{(n)}=\frac{\gamma^{(n-1)}}{\eta}$ and $\eta>1$ is a constant.
Here,
$\mathcal{K}(\mathbf{X})=\Phi_{\Gamma}(\mathbf{X})^{T}\Phi_{\Gamma}(\mathbf{X})$.
We use the _kernel-trick_ to evaluate $\mathcal{K}\left(\mathbf{X}\right)$,
thus eliminating the need to explicitly evaluating the features. Each entry of
the matrices $\mathcal{K}\left(\mathbf{X}\right)$ correspond to inner-products
in feature space:

$\displaystyle\left(\mathcal{K}\left(\mathbf{X}\right)\right)_{(i,j)}$
$\displaystyle=$
$\displaystyle\underbrace{\Phi(\mathbf{x}_{i})^{H}\Phi(\mathbf{x}_{j})}_{\kappa(\mathbf{x}_{i},\mathbf{x}_{j})}$
(33)

which can be evaluated efficiently using the nonlinear function $\kappa$
(termed as kernel function) of their inner-products in $\mathbb{R}^{n}$. The
dependence of the kernel function on the lifting is detailed in Section 5.3.

The details of the above minimization algorithms can be found in [50] in the
context of curve denoising; the present work extends the above result to the
multidimensional setting. Once the denoised null-space matrix is obtained, we
use the sum of square approach described in Section 4.1.3 to recover the
original surface.

We illustrate this approach in the context of recovering 3D shapes from noisy
point clouds in Fig. 8. The data sets are obtain from AIM@SHAPE [1]. We note
that the direct approach, where the null-space of the original noisy feature
matrix, often results in perturbed shapes. By contrast, the nuclear norm prior
is able to regularize the recovery, thus yielding shapes that are closer to
the original shapes. The parameters of the algorithm are tabulated in Table 1.

(a) Torus (b) 200 samples

(c) Recovery #1 (d) Recovery #2

(e) Cat

(f) 1500 samples (g) Recovery #1

(h) Recovery #2

(i) Duck

(j) 1800 samples (k) Recovery #1

(l) Recovery #2

Figure 8: Illustration of the points cloud denoising algorithm and surface recovery algorithm with unknown support. The first column shows the three original surfaces. The second column indicates the number of samples that we used to recover the surfaces. The third column corresponds to the surface recovery without using the points cloud denoising algorithm (30) beforehand. The last column gives the recovered surfaces using the denoising algorithm (30) beforehand. The parameters for the recovery of these three surfaces, including the number of iterations used in the denoising algorithm, the balancing parameter $\lambda$ in (30), and etc., in the recovery algorithm are reported in Table 1. Name of surface | # samples | # denoise iteration |  $\lambda$ in (30) | chosen support
---|---|---|---|---
Torus | 200 | 3 | 0.8 | $7\times 7\times 7$
Cat | 1500 | 15 | 1.4 | $31\times 31\times 31$
Duck | 1800 | 20 | 1.2 | $27\times 27\times 27$
Table 1: The parameters chosen for the surfaces recovery in Fig. 8.

##  5 Recovery of functions on surfaces

One of the main tasks in machine learning is to learn a function using given
input and output data pairs [19]. It is of paramount importance in control
[8], communications [49], dynamical system identification [20] and has wide
applications in language processing [18], image classification [12] and speech
recognition [37].

The problem can be formulated as estimating a nonlinear function
$y=f(\mathbf{x}),\,\,\,\mathbf{x}=(x_{1},x_{2},\cdots,x_{n})$, given input and
output data pairs $(\mathbf{x}_{i},y_{i})$. A challenge in the representation
of such high dimensional function is the large number of parameters and the
increasing of the dimensions, which is also termed as the curse of dimensions.
Kernel methods [29], random forests [39] and neural networks [47] provide a
powerful class of machine learning models that can be used in learning highly
nonlinear functions. These models have been widely used in many machine
learning tasks [17].

We now show that the results shown in the previous sections provide an
attractive option to compactly represent functions, when the data lies on a
smooth surface or manifold in high dimensional spaces. We note that the
manifold assumption is widely assumed in a range of machine learning problems
[13, 10]. We now show that if the data lie on a smooth surface in high
dimensional space, one can represent the multidimensional functions very
efficiently using few parameters. We now consider the representation of a
band-limited multidimensional function:

$f(\mathbf{x})=\sum_{\mathbf{k}\in\Gamma}\beta_{\mathbf{k}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x})=\bm{\beta}^{T}\Phi_{\Gamma}(\mathbf{x}),$
(34)

where $\mathbf{x}\in\mathbb{R}^{n}$. The number of free parameters in the
above representation is $|\Gamma|$, where $\Gamma\subset\mathbb{Z}^{n}$ is the
support of the function. Note that $|\Gamma|$ grows rapidly with the dimension
$n$. The large number of parameters needed for such a representation makes it
difficult to learn such functions from few labeled data points. Fortunately,
natural datasets often lie on simpler constructs in high dimensional space. We
assume that we are only interested in the points on the union of irreducible
surfaces as in (13), where the support of $\psi$ is given by
$\Lambda\subset\Gamma$.

###  5.1 Compact representation of features using anchor points

We use the upper bound of the dimension of the feature matrix in (23) to come
up with an efficient representation of functions of the form 34. The dimension
bound (23) implies that the features of points on $S[\psi]$ lies in a subspace
of dimension $r=|\Gamma|-|\Gamma\ominus\Lambda|$, which is far smaller than
$|\Gamma|$ especially when the dimension $n$ is large; one can then find $r$
basis vectors in the feature space that span
$\mathcal{V}_{\Gamma}(\mathcal{S})$. However, there is no guarantee that these
basis vectors are mappings of some points on $\mathcal{S}$.

We now show that it is possible to find a set of $N\geq r$ anchor points
$\mathbf{a}_{1},\cdots,\mathbf{a}_{N}\in\mathcal{S}[\psi]$, such that the
feature space $\mathcal{V}_{\Gamma}(\mathcal{S})$ is in
${\rm{span}}\\{\Phi_{\Gamma}(\mathbf{a}_{1}),\cdots,\Phi_{\Gamma}(\mathbf{a}_{N})\\}$.
This result is a Corollary of Proposition 9.

######  Corollary 11.

Let $\psi(\mathbf{x})$ be a trigonometric polynomial with $M$ irreducible
factors as in (27). Suppose $\Gamma_{i}\supset\Lambda_{i}$ is the non-minimal
support of each factor $\psi_{i}(\mathbf{x})$ and $\Gamma\supset\Lambda$ is
the total support. Let $\\{\mathbf{a}_{1},\cdots,\mathbf{a}_{N}\\}$ be $N$
randomly chosen anchor points on $\mathcal{S}[\psi]$ satisfying

  1. 1.

each irreducible factor $\mathcal{S}[\psi_{i}]$ is sampled with
$N_{i}\geq|\Gamma_{i}|-|\Gamma_{i}\ominus\Lambda_{i}|$ points, and

  2. 2.

the total number of samples satisfy $N\geq|\Gamma|-|\Gamma\ominus\Lambda|$.

Then,

$\mathcal{V}_{\Gamma}(\mathcal{S})\subseteq{\rm
span}\left\\{\Phi_{\Gamma}(\mathbf{a}_{i});i=1,\cdots,N\right\\}$ (35)

with probability 1.

As discussed in Section 4.1.3, if we randomly choose
$N\geq|\Gamma|-|\Gamma\ominus\Lambda|=r$ points on $\mathcal{S}[\psi]$, the
feature matrix will satisfy the conditions in Corollary 11 and hence (35) with
high probability. This relation implies that the feature vector of any point
$\mathbf{x}\in\mathcal{S}[\psi]$ can be expressed as the linear combination of
the features of the anchor points
$\Phi_{\Gamma}(\mathbf{a}_{i});i=1,\cdots,N$:

$\displaystyle\Phi_{\Gamma}(\mathbf{x})$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\mathbf{\alpha}_{i}(\mathbf{x})~{}\Phi_{\Gamma}(\mathbf{a}_{i})$
(36) $\displaystyle=$
$\displaystyle\underbrace{\begin{bmatrix}\Phi_{\Gamma}(\mathbf{a}_{1})&\cdots&\Phi_{\Gamma}(\mathbf{a}_{N})\end{bmatrix}}_{\Phi(\mathbf{A})}\underbrace{\begin{bmatrix}\alpha_{1}(\mathbf{x})\\\
\vdots\\\ \alpha_{N}(\mathbf{x})\end{bmatrix}}_{\bm{\alpha}(\mathbf{x})}$ (37)

Here, $\alpha_{i}(\mathbf{x})$ are the coefficients of the representation.
Note that the complexity of the above representation is dependent on $N$,
which is much smaller than $|\Gamma|$, when the surface is highly band-
limited. We note that the above compact representation is exact only for
$\mathbf{x}\in\mathcal{S}[\psi]$ and not for arbitrary
$\mathbf{x}\in\mathbb{R}^{n}$; the representation in (37) will be invalid for
$\mathbf{x}\notin\mathcal{S}[\psi]$.

However, this direct approach requires the computation of the high dimensional
feature matrix, and hence may not be computationally feasible for high
dimensional problems. We hence consider the normal equations and solve for
$\bm{\alpha}(\mathbf{x})$ as

$\bm{\alpha}(\mathbf{x})=\left(\underbrace{\Phi(\mathbf{A})^{H}\Phi(\mathbf{A})}_{\mathcal{K}\left(\mathbf{A}\right)}\right)^{{\dagger}}\underbrace{\left(\Phi(\mathbf{A})^{H}\Phi_{\Gamma}(\mathbf{x})\right)}_{\mathbf{k}_{\mathbf{A}}(\mathbf{x})},$
(38)

where $(\cdot)^{\dagger}$ denotes the pseudo-inverse.

We note that Corollary 11 guarantees $\Phi(\mathbf{A})$ to have full column
rank as $N=r$. However, the condition number of this matrix may be poor,
depending on the choice of the anchor points. It may be worthwhile to choose
the anchors such that the condition number of $\Phi(\mathbf{A})$ is low, which
will reduce the noise amplification in (38). One may pose this as an
optimization problem; we will discuss it in detail in a later section.

###  5.2 Representation and learning of functions

Using (34), (37), and (38), the function
$\mathbf{f}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{m}$ can be written as

$\displaystyle\mathbf{f}(\mathbf{x})$ $\displaystyle=$
$\displaystyle\bm{\beta}^{T}~{}\Phi({\mathbf{A}})~{}\mathcal{K}\left(\mathbf{A}\right)^{{\dagger}}~{}\mathbf{k}_{\mathbf{A}}(\mathbf{x})$
(39) $\displaystyle=$
$\displaystyle\underbrace{\left[\overbrace{\bm{\beta}^{T}~{}\Phi_{\Gamma}(\mathbf{a}_{1})}^{\mathbf{f}(\mathbf{a}_{1})},\ldots,\overbrace{\bm{\beta}^{T}~{}\Phi_{\Gamma}(\mathbf{a}_{N})}^{\mathbf{f}(\mathbf{a}_{N})}\right]}_{\mathbf{F}}~{}\underbrace{\mathcal{K}\left(\mathbf{A}\right)^{{\dagger}}~{}\mathbf{k}_{\mathbf{A}}(\mathbf{x})}_{\bm{\alpha}(\mathbf{x})}$
(40)

Here, $\mathbf{f}(\mathbf{x})$ is an $M\times 1$ vectors, while $\mathbf{F}$
is an $M\times N$ matrix. $\mathcal{K}(\mathbf{A})$ is an $N\times N$ matrix
and $\mathbf{k_{A}}(\mathbf{x})$ is an $N\times 1$ vector. Thus, if the
function values at the anchor points, specified by
$\mathbf{f}(\mathbf{a}_{i});i=1,\cdots,N$ are known, one can compute the
function for any point $\mathbf{x}\in\mathcal{S}[\psi]$.

We note that the direct representation of a function
$f:\mathbb{R}^{n}\rightarrow\mathbb{R}$ in (34) requires $|\Gamma|$
parameters, which can be viewed as the area of the green box in Fig. 3. By
contrast, the above representation only requires
$|\Gamma|\ominus|\Gamma:\Lambda|$ anchor points, which can be viewed as the
area of the grey region in Fig. 3. The more efficient representation allows
the learning of complex functions from few data points, especially in high
dimensional applications.

We demonstrate the above local function representation result in a 2D setting
in Fig. 9. Specifically, the original band-limited function is with support
$13\times 13$. The direct representation of the function has $13\times 13=169$
degrees of freedom. Now, if we only care about points on a curve which is with
support $3\times 3$, then the same function living on the curve can be
represented exactly using 48 anchor points, thus significantly reducing the
degrees of freedom. However, note that the above representation is only exact
on the curve. We note that the function goes to zero as one moves away from
the curve.

The choice of anchor points depends on the geometry of the surface, including
the number of irreducible components. While a random choice of anchor points
satisfies (35) with high probability; the choice of the anchor points controls
the condition number of $\Phi(\mathbf{A})$ and needs to be chosen carefully
for robust representation of functions. By contrast, one often does not have
perfect control on the data pairs
$(\mathbf{x}_{i},\mathbf{y}_{i}=\mathbf{f}(\mathbf{x}_{i}))$ themselves. For
arbitrary training samples, we can estimate the unknowns $\mathbf{F}$ in (40)
from the linear relations

$\displaystyle\underbrace{\left[\mathbf{y}_{1},..\mathbf{y}_{P}\right]}_{\mathbf{Y}}$
$\displaystyle=$
$\displaystyle\mathbf{F}~{}\underbrace{[\bm{\alpha}(\mathbf{x}_{1}),\ldots,\bm{\alpha}(\mathbf{x}_{P})]}_{\mathbf{Z}}$
(41)

as
$\mathbf{F}=\mathbf{Y}\mathbf{Z}^{H}\left(\mathbf{Z}\mathbf{Z}^{H}\right)^{{\dagger}}$.
The above recovery is exact when we have $N=r$ achor points because
$\mathbf{Z}$ has full column rank in this case. The reason why $\mathbf{Z}$
has full column rank is due to (37) and (38). Equation (37) suggests that
${\rm{rank}}(\mathbf{Z})\geq N$, while equation (38) shows
${\rm{rank}}(\mathbf{Z})\leq N$. Therefore, we have
${\rm{rank}}(\mathbf{Z})=N$, indicating that $\mathbf{Z}$ has full rank in
this case. When $N>r$, the $\mathbf{F}$ is obtained using the pseudo-inverse,
which is based on the least square approximation.

(a) Curve (b) band-limited function

(c) Function on curve

(d) Anchor points (e) Approximation

(f) Approx on curve

Figure 9: Illustration of the local representation of functions in 2D. We
consider the local approximation of the band-limited function in (b) with a
support of $13\times 13$, living on the band-limited curve shown in (a). The
support of the curve is $3\times 3$. The curve is overlaid on the function in
(b) in yellow. The restriction of the function to the vicinity of the curve is
shown in (c). Our results suggest that the local function approximation
requires $13^{2}-11^{2}=48$ anchor points. We randomly select the points on
the curve, as shown in (d). The interpolation of the function values at these
points yields the global function shown in (e). The restriction of the
function to the curve in (f) shows that the approximation is good.

###  5.3 Efficient computation using _kernel trick_

We use the _kernel-trick_ to evaluate $\mathcal{K}\left(\mathbf{A}\right)$ and
$\mathbf{k}_{\mathbf{A}}(\mathbf{x})$, thus eliminating the need to explicitly
evaluating the features of the anchor points and $\mathbf{x}$. Each entry of
the matrix $\mathcal{K}\left(\mathbf{A}\right)$ is computed as in (33), while
the vector $\mathbf{k}_{\mathbf{A}}(\mathbf{x})$ is specified by:

$\displaystyle(\mathbf{k}_{\mathbf{A}}(\mathbf{x}))_{i}$ $\displaystyle=$
$\displaystyle\underbrace{\Phi_{\Gamma}(\mathbf{a}_{i})^{H}\Phi_{\Gamma}(\mathbf{x})}_{\kappa(\mathbf{a}_{i},\mathbf{x})},$
(42)

which can be evaluated efficiently as nonlinear function $\kappa$ (termed as
kernel function) of their inner-products in $\mathbb{R}^{n}$. We now consider
the kernel function $\kappa$ for specific choices of lifting.

Using the lifting in (19), we obtain the kernel as

$\kappa(\mathbf{x},\mathbf{y})=\sum_{\mathbf{k}\in\Gamma}\exp(j2\pi\mathbf{k}^{T}(\mathbf{y}-\mathbf{x})).$

Note that the kernel is shift invariant in this setting. Since
$\kappa:\mathbb{R}^{n}\rightarrow\mathbb{R}$ is an $n$ dimensional function,
evaluating and storing it is often challenging in multidimensional
applications. We now focus on approximating the kernel efficiently for fast
computation.

We consider the impact of the shape of the support set $\Gamma$ on the shape
of the kernel. Specifically, we consider sets of the form

$\Gamma=\\{\mathbf{k}\in\mathbb{Z}^{n},||\mathbf{k}||_{q}\leq d\\}$ (43)

where $d$ denotes the size of the support. The integer $q$ specifies the shape
of $\Gamma$ [45], which translates to the shape of the kernel

$k_{d,n}^{q}(\mathbf{x}):=\sum_{\mathbf{k}\in\mathbb{Z}^{n},||\mathbf{k}||_{q}\leq
d}\exp(j2\pi\mathbf{k}^{T}\mathbf{x}),$ (44)

We term the $q=1$ case as the diamond Dirichlet kernel. If $q=2$, we call it
the circular Dirichlet kernel. We call the Dirichlet kernel the triangular
Dirichlet kernel if $q=\infty$. See Figure 10 for the support and Figure 11 to
see the associated kernel.

(a) $q=1$

(b) $q=2$

(c) $q=\infty$

Figure 10: Support of the set $\Lambda$ with different $q$ values.

(a) Gaussian kernel (b) Dirichlet with $q=2$

(c) Dirichlet with $q=\infty$ (d) Plot of $\gamma$

Figure 11: Visualization of kernels in $\mathbb{R}^{2}$ and the non-linear
function $\gamma$ with some commonly used activation functions.

We note from the above figures that the circular Dirichlet kernel ($q=2$) is
roughly circularly symmetric, unlike the triangular or diamond kernels. This
implies that we can safely approximate it as

$\kappa(\mathbf{x},\mathbf{y})\approx g(\|\mathbf{x}-\mathbf{y}\|^{2})$ (45)

where $g:\mathbb{R}_{+}\rightarrow\mathbb{R}$. We note that this approximation
results in significantly reduced computation in the multidimensional case. The
function $g$ may be stored in a look-up table or computed analytically. We use
this approach to speed up the computation of multidimensional functions in
Section 6.

An additional simplification is to assume that $\mathbf{x}$ and $\mathbf{y}$
are unit-norm vectors. In this case, we can approximate

$g(||\mathbf{x}_{i}-\mathbf{y}_{i}||_{2}^{2})=g(||\mathbf{x}_{i}||_{2}^{2}+||\mathbf{y}_{i}||_{2}^{2}-2\langle\mathbf{x}_{i},\mathbf{y}_{i}\rangle)\approx
g(2-2\langle\mathbf{x},\mathbf{y}\rangle)=:\gamma(\langle\mathbf{x},\mathbf{y}\rangle),$
(46)

where $\gamma(z)=g(1-z/2)$. Here, we term $\gamma$ as the activation function.
While we do not make this simplifying assumption in our computations, it
enables us to show the similarity of the computational structure of (39) to
current neural network. The plot of this activation function, along with
commonly used activation functions, is shown in Figure 11 (d).

##  6 Relation to neural networks

We now focus on the computation of an arbitrary function in (39) and will show
its equivalence with neural networks.

$\displaystyle\mathbf{f}(\mathbf{x})$ $\displaystyle=$
$\displaystyle\underbrace{\left[\mathbf{f}_{1},\ldots,\mathbf{f}_{n}\right]}_{\mathbf{F}}~{}\mathcal{K}\left(\mathbf{A}\right)^{{\dagger}}~{}\underbrace{\begin{bmatrix}g(\|\mathbf{x}-\mathbf{a}_{1}\|^{2})\\\
\vdots\\\
g(\|\mathbf{x}-\mathbf{a}_{n}\|^{2})\end{bmatrix}}_{\mathbf{k}_{\mathbf{A}}(\mathbf{x})}$
(47) $\displaystyle\approx$
$\displaystyle\underbrace{\mathbf{F}~{}\mathcal{K}(\mathbf{A})^{{\dagger}}}_{\tilde{\mathbf{F}}}~{}\begin{bmatrix}\gamma\left(\left\langle\mathbf{x},\mathbf{a}_{1}\right\rangle\right)\\\
\vdots\\\
\gamma\left(\left\langle\mathbf{x},\mathbf{a}_{n}\right\rangle\right)\end{bmatrix}$
(48)

In the second step, we used the approximation in (46). The computation
involves the inner product of the input signal $\mathbf{x}$ with templates
$\mathbf{a}_{i};i=1,..,n$, followed by the non-linear activation function
$\gamma$ to obtain $\mathbf{k}_{\mathbf{A}}(\mathbf{x})$. These terms are then
weighted by the fully connected layer
$\mathcal{K}\left(\mathbf{A}\right)^{{\dagger}}$, followed by weighting by the
second fully connected layer $\tilde{\mathbf{F}}$. See Fig. 12 for the visual
illustration.

(a) One layer network

(b) Two layers network

Figure 12: Computational structure of function evaluation. (a) corresponds to
(39) to compute the band-limited multidimensional function $\mathbf{f}$ on
$\mathcal{S}[\psi]$. The inner-product between the input vector $\mathbf{x}$
and the anchor templates on the surface are evaluated, followed by non-linear
activation functions $\gamma$ to obtain the coefficients
$\alpha_{i}(\mathbf{x})$. These coefficients are operated with the fully
connected linear layers $\mathbf{K}_{\mathbf{A}}^{{\dagger}}$ and
$\mathbf{F}(\mathbf{A})$. The fully connected layers can be combined to obtain
a single fully connected layer $\tilde{\mathbf{F}}$. Note that this structure
closely mimics a neural network with a single hidden layer. (b) uses an
additional quadratic layer, which combines functions of a lower support to
obtain a function of a higher support.

###  6.1 Illustration in a denoising application

We now illustrate the proposed network in the context of a toy-problem: the
denoising of MNIST data [21]. Specifically, we consider the learning of a
function $f:\mathbb{R}^{p^{2}}\rightarrow\mathbb{R}$, which predicts the
denoised center pixel of a patch from the noise $p\times p$ patch. The
function $f$ in $p^{2}$ dimensional space is associated with a large number of
free parameters; learning of these unknowns are challenging due to the curse
of dimensionality. In this work, we assume that the points (image patches) are
lying on a band-limited surface. Then the result in the previous section
offers a work-around, which suggests that the function can be expressed as the
linear combination of the features of “anchor-patches”, weighted by
$\mathbf{p}$. In our experiments, we choose the patch size to be $7\times 7$.

We propose to learn the anchor patches $\mathbf{a}_{i}$ and the function
values $f(\mathbf{a}_{i})$ from exemplar data using stochastic gradient
descent. Note that the learned representation is valid for any patch, and
hence the proposed scheme is essentially a convolutional neural network. The
difference of our structure in (48) with the commonly used convolutional
neural networks (CNN) structure is the activation function $\gamma$. We
replaced the ReLU non-linearity in a network with a single hidden layer with
the proposed function $\gamma$ and trained the algorithm on the MNIST dataset
with added noise. We also trained a ReLU network with the same parameters for
comparison. Besides, we applied our two-layer network and the two-layer ReLU
network on the MNIST dataset as well. The comparison of the testing results is
shown in Figure 13.

Figure 13: Comparison of our learned denoiser using Dirichlet activation
function and the ReLU network. The testing results show that the denoising
performance using the proposed activation function is comparable to the
performance using ReLU. The six rows in the figure correspond to the original
images, the noisy images, the denoised images using the proposed one-layer
network, the denoised images using one layer ReLU network, the denoised images
using the proposed two-layer network and the denoised images using two-layer
ReLU network. The averaged PSNR of the denoised images using the proposed one-
layer network, one layer ReLU network, the proposed two-layer network and the
two-layer ReLU network are $19.6753$ dB, $20.0254$ dB, $20.8576$ dB and
$17.4818$ dB respectively. From the quantitative results, we can see that our
proposed one-layer network performs almost the same as the one-layer ReLU
network. For the proposed two-layer network, the performance is getting better
from both quantitative and visual points of view. For the two-layer ReLU
network, visually the performance is better than that of the one-layer ReLU
network. But the PSNR is getting worse. The main reason that causes the low
PSNR for the two-layer ReLU network is the change of the pixel values on each
hand-written digits.

The convergence of the training loss of the proposed models are shown in
Figure 14. We note that the final training errors of both single layer
networks are comparable. We observe that both of the proposed schemes (one-
layer and two-layer) converged rapidly to a low error, while the ReLU network
took more iterations to converge. Note that our proposed activation function
has support size $d$ as a parameter, indicated by (44). We reported the impact
of the support on the denoising performance in Fig. 14 (d). From the table, we
can see that the higher support usually produces better performance. But
support will not significantly affect the performance.

(a) Proposed vs ReLU (one layer) (b) Proposed vs ReLU (two layers)

(c) One layer vs two layers (proposed)

(d) Variation of the PSNR of the denoised results with the support of the
kernel

Figure 14: Comparison of the training loss for the proposed networks and ReLU
network, and the variation of the PSNR of the denoised results with the
support of the kernel. (a) shows the training loss of the proposed one-layer
network and the one-layer ReLU network. We can see that the training loss of
the proposed model decreases faster than the ReLU network. After certain
number of epochs, the training loss curves of both the proposed one-layer
network and the ReLU network are very similar. (b) compares the training loss
of the proposed two-layer network with the two-layer ReLU network. We can find
that the proposed two-layer network performs slightly better than the two-
layer ReLU network. In (c), the training losses of the proposed one-layer
network and the proposed two-layer network are compared. We can see that the
two-layer network is better than the one-layer network, which demonstrates the
necessity for going in-depth. In (d), we reported the variation of the PSNR
(dB) of the denoised results with the support of the kernel. In the
experiments, we choose the support $d$ to be 7, 11 and 15 respectively and we
then compared the performance. We did the testing on the first 100 testing
images in the MNIST dataset obtained from TensorFlow. From the reported
average PSNR, we can find that higher support usually gives better PSNR, while
the effective is not significant.

###  6.2 Relation to denoising auto-encoders

We note that the space of band-limited functions of the form (34) can
reasonably approximate lower order polynomials in $\mathbb{R}^{n}$ for
sufficiently high support $\Gamma$ [41]. In particular, let us assume that
there exists a set of coefficients $\bm{\beta}$ such that

$\mathbf{x}\approx\tilde{\mathbf{x}}=\sum_{\mathbf{k}\in\Gamma}\beta_{\mathbf{k}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x})$
(49)

In this case, the above results imply that one can exactly represent any point
on the surface $\mathcal{S}[\psi]$ as

$\displaystyle\mathbf{x}$ $\displaystyle\approx$
$\displaystyle\underbrace{\left[\mathbf{a}_{1},..,\mathbf{a}_{n}\right]}_{\mathbf{A}}~{}\underbrace{\mathcal{K}\left(\mathbf{A}\right)^{{\dagger}}~{}\mathbf{k}_{\mathbf{A}}(\mathbf{x})}_{\bm{\alpha}(\mathbf{x})}$
(50)

We note that the resulting network is hence essentially an auto-encoder.
Specifically, the inner-products between the feature vectors of $\mathbf{x}$
and the anchor point $\mathbf{a}_{i}$ denoted by
$\alpha(\mathbf{\mathbf{missing}}x)$ can be viewed as the latent features or
compact code. As described previously, the coefficients
$\bm{\alpha}=\mathcal{K}(\mathbf{A})^{{\dagger}}\mathbf{k_{A}}(\mathbf{x})$
captures the geometry of the surface, while the top layer $\mathbf{A}$ is the
decoder that recover the signal from its latent vectors.

We note that the surface recovery algorithms in Section 4 follow a null-space
approach, where we identify the null-space of the feature space or
equivalently the annihilation functions from the samples of the surface.
Specifically, the sum of squares of the null-space functions in Section 4.1.4
provides a measure of the error in projecting the feature vector to the null-
space of the feature matrix.

$\displaystyle\gamma(\mathbf{x})$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{|\Gamma\ominus\Lambda|}|\mu_{i}(\mathbf{x})|^{2}=\sum_{i=1}^{|\Gamma\ominus\Lambda|}|\mathbf{n}_{i}^{T}\Phi_{\Gamma}(\mathbf{x})|^{2}$
(51) $\displaystyle=$
$\displaystyle\|\mathbf{N}~{}\Phi_{\Gamma}(\mathbf{x})\|^{2}$ (52)

where $\mathbf{n}_{i}$ are the null-space vectors. The projection energy is
zero if the point $\mathbf{x}$ is on $\mathcal{S}$ and is high when it is far
from it.

By contrast, the auto-encoder approach can be viewed as a signal subspace
approach, where we project the samples to the basis vectors specified by the
feature vectors of the anchors $\Phi_{\Gamma}(\mathbf{a}_{i})$. We illustrate
this approach in Figure 15. Specifically, we trained a two layer neural
network described by (50). We assumed the non-linearity to be specified by
$\eqref{approx}$ and trained the network parameters ($\mathbf{A}$ as well as
the weights of the inner-products) using stochastic gradient descent. The
training data corresponds to randomly drawn points on the surface. To ensure
that the network learns a projection, we trained the network as a denoising
auto-encoder; the inputs correspond to samples on the surface corrupted with
Gaussian noise, while the labels are the true samples. Once the training is
complete, we plot the approximation error

$E(\mathbf{x})=\|\mathbf{x}-\mathbf{F}\mathcal{K}(\mathbf{A})^{{\dagger}}\mathbf{k}_{\mathbf{A}}(\mathbf{x})\|^{2}=\|\underbrace{\left(\mathbf{I}-\mathbf{F}\mathcal{K}(\mathbf{A})^{{\dagger}}\mathbf{k}_{\mathbf{A}}\right)}_{\mathcal{R}}(\mathbf{x})\|^{2}$
(53)

as a function of the input point in Fig. 15.

We trained the network using the exemplar curve shown in Fig. 7. We randomly
choose 1000 points on the curve as the training data and 250 features are
chosen in the middle layer. The support of the Dirichlet kernel is chosen to
be $15$. The trained network is then used to learn the curve. The learned
results are shown in Fig. 15. From which one can see that the proposed
learning framework performs well. We note that the projection error is close
to zero on the surface, while it is high away from it. Note that this closely
mimics the plot in Fig. 7. Once trained, the surface can be estimated in low-
dimensional settings as the zero set of the projection error as shown in Fig.
15.(b), which closely approximates the true curve in (c). We note that
$\mathcal{R}$ can be viewed as a residual denoising auto-encoder. Once
trained, this network can be used as a prior in inverse problems as in [2]. We
have used the null-space network in Section 6. We have also used the null-
space prior (51) in our prior work [32], where the null-space basis was
learned as described in Section 4.2.

(a) Learned curve (b) Contour line

(c) Original curve

Figure 15: Illustration of the surface learning network using the curve in
Fig. 7. (a) and (b) are the learned results. We compared the learned curve
(blue curve) with the original curve (red curve) in (c). From which we see
that the two curves are almost the same, indicating that the learned network
performs well.

###  6.3 Benefit of depth and pooling

As noted above, the representation using anchor points reduces the degrees of
freedom significantly compared to the direct representation. However, we note
that the number of parameters needed to represent a high bandwidth function in
high dimensions is still high. We now provide some intuition on how the low-
rank tensor approximation of functions and composition can explain the benefit
of common operations in deep networks.

We now consider the cases when the band-limited multidimensional function
$f:\mathbb{R}^{n}\rightarrow\mathbb{R}$ in (34) can be approximated as

$f(\mathbf{x})=\left(\sum w_{i}~{}f_{i}(\mathbf{x})\right)^{2}.$ (54)

Clearly, the Fourier support of $f$ is almost twice that of
$f_{i}:\mathbb{R}^{n}\rightarrow\mathbb{R}$, showing the benefit of adding
layers. While an arbitrary function with the same support as $f$ cannot be
represented as in (54), one may be able to approximate it closely. The new
layer will have a quadratic non-linearity $Q$, if the function has the form
(54). Note that one may use arbitrary non-linearity in place of the quadratic
one in (54).

Similarly, one may perform a low-rank tensor approximation of an arbitary $N$
dimensional function $f:\mathbb{R}^{n}\rightarrow\mathbb{R}$. Specifically,
the approximation involves the sum of products of 1-D functions.

$f(x_{1},..,x_{N})\approx\sum_{i=1}^{r}h_{1}^{(i)}(x_{1})\cdot~{}h_{2}^{(i)}(x_{2})\ldots~{}h_{N}^{(i)}(x_{N}),$
(55)

where $h_{i}:\mathbb{R}\rightarrow\mathbb{R}$. The above sum of products can
also be realized by taking weighted linear combination of 1-D functions,
followed by a non-linearity as in (54). This allows one to have a hierarchical
structure, where lower dimensional functions are pooled together to represent
a multidimensional function.

In image processing applications, the functions to be learned are shift-
invariant. This allows one to learn functions of small image patches (e.g.
$3\times 3$) of a specified dimension at each layer. The functions on nearby
pixels in the output thus correspond to information from different $3\times 3$
neighborhoods. The low-dimensional functions from non-overlapping $3\times 3$
neighborhoods could be combined with downsampling as in (55) to represent a
high dimensional function (e.g. $9\times 9$) neighborhoods. The process can be
repeated to improve the efficiency of representation.

##  7 Conclusion

In this work, we considered a data model, where the signals are localized to a
surface that is the zero level-set of a band-limited function $\psi$. The
support of the function can be seen as a complexity measure of the surface. We
show that the non-linear features of the samples, obtained by an exponential
lifting, satisfies an annihilation relation. Using the annihilation relation,
we developed theoretical sampling guarantees for the unique recovery of the
surface. Our main contribution here is to prove that with probability 1, the
surface can be uniquely recovered using a collection of samples, whose number
is equal to the degrees of freedom of the representation. When the true
support of the surface is unknown, which is usually the case, we introduced a
method using the SoS polynomial to specify the surface. We also introduced the
way to get back the samples when the original samples are corrupted by noise.

We then use this model to efficiently represent arbitrary band-limited
functions $f$ of points living on the surface. We show that the exponential
features of the points on the surface live in a low-dimensional subspace. This
subspace structure is used to represent the $f$ efficiently using very few
parameters. We note that the computational structure of the function
evaluation mimics a single-layer neural network. We applied the proposed
computational structure to the context of image denoising.

##  8 Appendix

###  8.1 Proof of Proposition 1

As we mentioned in Section 2.3.1, if we have a (hyper-)surface $\mathcal{S}$
which is given by the zero level set of a trigonometric polynomial, then there
will be a minimal polynomial which defines $\mathcal{S}$ (Proposition 1). To
prove this result, we need to following famous result.

######  Lemma 12 (Hilbert’s Nullstellensatz [3]).

Let $\mathbb{K}$ be an algebraically closed field (for example $\mathbb{C}$).
Suppose $I\subset\mathbb{K}[x_{1},\cdots,x_{n}]$ is an ideal of polynomials,
and $\mathcal{Z}(I)$ denotes the set of common zeros of all the polynomials in
$I$. Let $\mathcal{I}(\mathcal{Z}(I))$ represents the ideal of polynomials in
$\mathbb{K}[x_{1},\cdots,x_{n}]$ vanishing on $\mathcal{Z}(I)$. Then, we have

$\mathcal{I}(\mathcal{Z}(I))=\sqrt{I},$

where $\sqrt{I}$ denotes the radical of $I$, specified by the set

$\sqrt{I}=:\\{p|p^{n}\in I,~{}\mbox{for some}~{}n\in\mathbb{Z}^{+}\\}$ (56)

Remark 1. We say a set $I\subset K[x_{1},\cdots,x_{n}]$ is an ideal, if $I$ is
closed under the addition operation (e.g. addition“$+$”), satisfies the
associative property, has a unit element $0$, and a valid inverse for every
element in $I$. For the operation multiplication (e.g. “$\cdot$”), we have
$r\cdot p\in I$ and $p\cdot r\in I$ for any $r\in K[x_{1},\cdots,x_{n}]$

Remark 2. An important property of the radical of the ideal $I$ is that
$I\subset\sqrt{I}$. Note that setting $n=1$ in (56) will yield $I$.

Remark 3. The above lemma states that the set of all polynomials that vanish
on the common zeros $\mathcal{Z}(I)$ of the polynomials in $I$ is given by
$\sqrt{I}\supset I$. Specifically, if we are given another polynomial
$\eta(\mathbf{x})$ that also vanishes on the common zero set $\mathcal{Z}(I)$,
then there must be positive integer $n$ such that $\eta^{n}(\mathbf{x})\in I$.

We denote the ideal generated by a function $f$ by
$(f)=\\{\mu|\mu=f\gamma\\}$, where $\gamma$ is an arbitrary polynomial. The
identity in this ideal is the zero polynomial. In particular, $(f)$ is the
family of all functions that have $f$ as a factor. We note that the set of
common zeros of all the functions in $(f)$, denoted by $\mathcal{Z}[(f)]$ is
the same as the zero set of $f$, denoted by $Z[f]$.

######  Lemma 13.

Let $f,g$ be two polynomials in $\mathbb{C}[x_{1},\cdots,x_{n}]$ with the same
zero set. Then the two polynomials must have (up to scaling) the same factors.

###### Proof.

Suppose $Z[f]=Z[g]=Z$ is the zero set of $f$ and $g$. Since
$Z[f]=\mathcal{Z}[(f)]$, we have $\mathcal{Z}[(f)]=\mathcal{Z}[(g)]=Z$. By the
Hilbert’s Nullstellensatz, we have

$\mathcal{I}(\mathcal{Z}(f))=\sqrt{(f)},\qquad\mathcal{I}(\mathcal{Z}(g))=\sqrt{(g)}.$

Since $Z(f)=Z(g)$, we then have
$\mathcal{I}(\mathcal{Z}(f))=\mathcal{I}(\mathcal{Z}(g))$ and hence
$\sqrt{(f)}=\sqrt{(g)}$. As mentioned above, we have $I\subset\sqrt{I}$ for
any ideal $I$. Therefore, we have $(f)\subset\sqrt{(f)}$ and
$(g)\subset\sqrt{(g)}$. This implies that $f\in\sqrt{(f)}$ and
$g\in\sqrt{(g)}$. Because we have $\sqrt{(f)}=\sqrt{(g)}$, we can obtain that
$f\in\sqrt{(g)}$ and $g\in\sqrt{(f)}$. By which we have that there exist
$m,n\in\mathbb{Z}$ and $p,q\in\mathbb{C}[x_{1},\cdots,x_{n}]$ such that

$f^{n}=p\cdot g,\qquad g^{m}=q\cdot f.$

Therefore, we can obtain that the irreducible factors of $g$ are of $f$ as
well and vice versa, which proves the desired conclusion. ∎

With this conclusion, we can now prove Proposition 1.

###### Proof of Proposition 1.

The proof of the existence and uniqueness about $\psi$ is same as the proof of
Proposition A.3 in [27] and thus we omit them here.

In this proof, we show that $supp(\psi)\subseteq supp(\psi_{1})$. Note that
the algebraic surface $X=\\{p=\mathcal{P}[\psi]=0\\}$ is the union of
irreducible surfaces $X_{j}=\\{p_{i_{j}}=0\\}\subset\mathbb{C}^{n}$. Define

$\nu(x_{1},\cdots,x_{n})=(e^{j2\pi x_{1}},\cdots,e^{j2\pi x_{n}}).$

Let $\mathcal{S}_{j}=\nu^{-1}(X_{j}\cap\mathbb{CT}^{n})$. Then we have a
decomposition of $\mathcal{S}$ as the union of surfaces $\mathcal{S}_{j}$. If
$\psi_{1}$ is another trigonometric polynomial with $\mathcal{S}$ as the zero
level set as well. Then $\psi_{1}$ vanishes on each $\mathcal{S}_{j}$. Let
$q=\mathcal{P}[\psi_{1}]$. Then we have $q=0$ on the infinite set
$\nu(\mathcal{S}_{j})$, by which we can infer that $q$ and $p$ will have the
same zero set using Theorem 14. Then by Lemma 13, we have $p\mid q$, which
implies that $supp(\psi)\subseteq supp(\psi_{1})$. ∎

###  8.2 Proof of results in Section 4

The key property of surfaces that we exploit is that the dimension of the
intersection of two band-limited surfaces of dimension $k$ is strictly lower
than $k$, provided their level set functions do not have any common factors.
Hence, if we randomly sample one of the surfaces, the probability that the
samples fall on the intersection of the two surfaces is zero. This result
enables us to come up with the sampling guarantees. We will now show the
results about the intersections of the zero sets of two trigonometric
surfaces.

####  8.2.1 Intersection of surfaces

We will first state a known result about the intersection of the zero sets of
two polynomials (non-trigonometric) whose level set functions do not have a
common factor.

######  Theorem 14 ([16],pp.115, Theorem 14).

Let $\mathcal{S}[\psi]$ and $\mathcal{S}[\eta]$ be two surfaces of dimension
$n-1$ over a field $\mathbb{K}$, which are the zero sets of the polynomials
$\psi:\mathbb{K}^{n}\rightarrow\mathbb{K}$ and
$\eta:\mathbb{K}^{n}\rightarrow\mathbb{K}$, respectively. If $\psi$ and $\eta$
do not have a common factor, then

$\dim\big{(}\mathcal{S}[\psi]\cap\mathcal{S}[\eta]\big{)}<n-1.$

The above result is a generalization of the two dimensional case
($\mathbb{C}^{2}$) in [27], where Bézout’s inequality was used to prove the
result. Specifically, the result in [27] suggests that the intersection of two
curves consists of a set of isolated points, if their potential function does
not have any common factor. Theorem 14 generalizes the above result to $n>2$;
it suggests that the intersection of two surfaces with dimension $k$ is
another surface, whose dimension is strictly less than $k$. For instance, the
intersection of two 3-D surfaces which are given by the zero level set of some
polynomials, could yield 2D curves or isolated points. We now extend Theorem
14 to trigonometric polynomials using the mapping $\nu$ specified by (6).

######  Lemma 15.

Let $\mathcal{S}[\psi]$ and $\mathcal{S}[\eta]$ within
$[0,1]^{n}\subset\mathbb{R}^{n}$ be two surfaces of dimension $n-1$ over
$\mathbb{R}$, which are the zero level sets of the trigonometric polynomials
$\psi$ and $\eta$. Suppose $\psi$ and $\eta$ do not have a common factor, then

$\dim(\mathcal{S}[\psi]\cap\mathcal{S}[\eta])<n-1.$

###### Proof.

Let $\nu=(\nu_{1},\cdots,\nu_{n})$ be defined by (6). We now would like to
prove the result by way of contradiction. Suppose

$\dim(\mathcal{S}[\psi]\cap\mathcal{S}[\eta])=\dim(\mathcal{S}[\psi])=\dim(\mathcal{S}[\eta])=n-1.$

This implies that $\nu(\mathcal{S}[\psi]\cap\mathcal{S}[\eta])$ will have the
same dimension of $\nu(\mathcal{S}[\psi])$ and $\nu(\mathcal{S}[\eta])$.
However, this is impossible according to Theorem 14. Therefore, we have the
desired result. ∎

Based on this lemma, we can directly have the following Corollary.

######  Corollary 16.

Suppose $\psi(\mathbf{x}),\eta(\mathbf{x}),\mathbf{x}\in[0,1]^{n}$ are two
trigonometric polynomials as in Lemma 15. Consider the $n-1$ dimensional
Lebesgue measure on $\mathcal{S}[\psi]$. Then this Lebesgue measure of the
intersection of the zero level sets of the trigonometric polynomials is zero,
i.e.,

$m(\mathcal{S}[\psi]\cap\mathcal{S}[\eta])=0.$

The Lebesgue measure can be viewed as the area of the $n-1$ dimensional
surface. For example, when $n=3$, $\mathcal{S}[\psi]$ and $\mathcal{S}[\eta]$
are 2-D surfaces, while their intersection is a 1-D curve or a set of isolated
points with zero area.

####  8.2.2 Proof of Proposition 4

###### Proof.

We note that $N\geq|\Lambda|-1$ is a necessary condition for the matrix to
have a rank of $|\Lambda|-1$. We now assume that the surface is sampled with
$N\geq|\Lambda|-1$ random samples, chosen independently, denoted by
$\mathbf{x}_{i};i=1,\cdots,N\in\mathcal{S}[\psi]$. Since
$\mathbf{c}\leftrightarrow\psi$ is a valid non-trivial null-space vector for
the feature matrix $\Phi_{\Lambda}(\mathbf{X})$ formed from these samples, we
have ${\rm rank}\left(\Phi_{\Lambda}(\mathbf{X})\right)\leq|\Lambda|-1$. The
polynomial $\psi(\mathbf{x})=\mathbf{c}^{T}\Phi_{\Lambda}(\mathbf{x})$ is the
minimal irreducible polynomial that defines the surface.

We now prove the desired result by contradiction. Assume that these exists
another linearly independent null-space vector $\mathbf{d}\leftrightarrow\eta$
or equivalently the rank of $\Phi_{\Lambda}(\mathbf{X})$ is strictly less than
$|\Lambda|-1$. Since $\mathbf{c}$ and $\mathbf{d}$ are linearly independent
and $\psi(\mathbf{x})$ is the minimal polynomial, we know that
$\psi(\mathbf{x})$ and $\eta(\mathbf{x})$ will not share a common factor. Also
note that $\mathbf{x}_{i}\in\mathcal{S}[\psi]\cap\mathcal{S}[\eta]$. However,
since $\psi(\mathbf{x})$ and $\eta(\mathbf{x})$ do not share a common factor,
the probability of each sample to be at the intersection of the two
polynomials $(\mathbf{x}_{i}\in\mathcal{S}[\psi]\cap\mathcal{S}[\eta])$ is
zero by Corollary 16. Therefore, with probability 1 that such $\mathbf{d}$
does not exist, meaning that with probability 1 that the feature matrix will
be of rank $|\Lambda|-1$ when $N\geq|\Lambda|-1$. ∎

####  8.2.3 Proof of Proposition 6

###### Proof.

We note that $N\geq|\Lambda|-1$ is a necessary condition for the matrix to
have a rank of $|\Lambda|-1$. We now assume that the surface is sampled with
$N$ random samples $\mathbf{x}_{i};i=1,\cdots,N$ satisfying the conditions in
Proposition 6. The minimal polynomial
$\psi(\mathbf{x})=\mathbf{c}^{T}\Phi_{\Lambda}(\mathbf{x})$ that defines the
surface can be factorized as
$\psi(\mathbf{x})=\psi_{1}(\mathbf{x})\cdot\psi_{2}(\mathbf{x})\cdots\psi_{M}(\mathbf{x})$.

We will prove the result by contradiction. Assume that these exists another
linearly independent null-space vector $\mathbf{d}\leftrightarrow\eta$, or
equivalently the rank of $\Phi_{\Lambda}(\mathbf{X})$ is less than
$|\Lambda|-1$. Since $\mathbf{c}$ and $\mathbf{d}$ are linearly independent,
$\psi$ and $\eta$ should differ by at least one factor. Without loss of
generality, let us assume that
$\eta(\mathbf{x})=\mu(\mathbf{x})\prod_{i=1}^{M-1}\psi_{i}(\mathbf{x})$, where
$\mu$ is an arbitrary polynomial of support $\Lambda_{M}$. Besides, $\mu$ and
$\psi_{M}$ does not share a factor. Using the result of Proposition 4, we see
that the probability of $\mu$ and an irreducible $\psi_{M}$ vanish at
$|\Lambda_{i}|-1$ independently drawn random locations is zero. If multiple
factors are shared, the same argument can be extended to each one of the
factors independently. ∎

####  8.2.4 Proof of Proposition 8

###### Proof.

We note that $N\geq|\Gamma|-|\Gamma\ominus\Lambda|$ is a necessary condition
for the matrix to have the specified rank. We now assume that the surface is
sampled with $N\geq|\Gamma|-|\Gamma\ominus\Lambda|$ random samples, chosen
independently. We note that $\mathbf{c}\leftrightarrow\psi$ specified by (8),
as well as the $|\Gamma\ominus\Lambda|$ translates of $\mathbf{c}$ within
$\Gamma$, are valid linearly independent null-space vectors of
$\Phi_{\Lambda}(\mathbf{X})$. We thus have

${\rm
rank}\left(\Phi_{\Lambda}(\mathbf{X})\right)\leq|\Gamma|-|\Gamma\ominus\Lambda|$
(57)

We will show that the rank condition can be satisfied with probability 1 by
contradiction. Assume that these exists another linearly independent null-
space vector $\mathbf{d}\leftrightarrow\eta$ or equivalently the rank of
$\Phi_{\Lambda}(\mathbf{X})$ is less than $|\Gamma|-|\Gamma\ominus\Lambda|$.
Since $\mathbf{d}$ are linearly independent with $\mathbf{c}$ and its
translates within $\Gamma$, we cannot express $\mathbf{d}$ as the linear
combinations of the the other null-space vectors. Specifically, we have

$\displaystyle\eta(\mathbf{x})$ $\displaystyle\neq$
$\displaystyle\sum_{\mathbf{k}\in\Gamma\ominus\Lambda}\alpha_{\mathbf{k}}~{}\psi(\mathbf{x})\exp(j2\pi\mathbf{k}^{T}\mathbf{x})$
(58) $\displaystyle=$
$\displaystyle\psi(\mathbf{x})\underbrace{\sum_{\mathbf{k}\in\Gamma\ominus\Lambda}\alpha_{\mathbf{k}}\exp(j2\pi\mathbf{k}^{T}\mathbf{x})}_{\gamma(\mathbf{x})}=\psi(\mathbf{x})\gamma(\mathbf{x}).$
(59)

Here $\alpha_{\mathbf{k}}$ is an arbitrary coefficients and hence $\gamma$ is
an arbitrary polynomial. The linear independence property implies that
$\eta(\mathbf{x})$ cannot have $\psi(\mathbf{x})$ as a factor. Since
$\psi(\mathbf{x})$ is the minimal polynomial, this also means that $\eta$ and
$\psi$ does not have any common factor.

Consider now the random sampling set
$\mathbf{x}_{i};i=1..|\Gamma|-|\Gamma\ominus\Lambda|$. We have

$\mathbf{c}^{T}\Phi_{\Lambda}(\mathbf{x}_{i})=\mathbf{d}^{T}\Phi_{\Lambda}(\mathbf{x}_{i})=0,\,\,\,i=1,\cdots,|\Lambda|-1.$

This implies that $\mathbf{x}_{i}\in\mathcal{S}[\psi]\cap\mathcal{S}[\eta]$.
However, since $\psi(\mathbf{x})$ and $\eta(\mathbf{x})$ do not share a common
factor, the probability of each sample to be at the intersection of the two
polynomials $(\mathbf{x}_{i}\in\mathcal{S}[\psi]\cap\mathcal{S}[\eta])$ is
zero by Corollary 16. Therefore, we have ${\rm
rank}\left(\Phi_{\Lambda}(\mathbf{X})\right)=|\Gamma|-|\Gamma\ominus\Lambda|$
with probability one.

∎

####  8.2.5 Proof of Proposition 9

###### Proof.

We note that $N\geq|\Gamma|-|\Gamma\ominus\Lambda|$ is a necessary condition
for the matrix to have the specified rank. We now assume that the surface is
sampled with $N$ random samples satisfying the sampling conditions in
Proposition 9. The minimal polynomial
$\psi(\mathbf{x})=\mathbf{c}^{T}\Phi_{\Lambda}(\mathbf{x})$ that defines the
surface can be factorized as
$\psi(\mathbf{x})=\psi_{1}(\mathbf{x})\cdot\psi_{2}(\mathbf{x})\cdots\psi_{M}(\mathbf{x})$.

Assume that there exists another linearly independent null-space vector
$\mathbf{d}\leftrightarrow\eta$ or equivalently the rank of
$\Phi_{\Lambda}(\mathbf{X})$ is less than $|\Gamma|-|\Gamma\ominus\Lambda|$.
Similar to the above arguments, if $\eta$ and $\psi$ does not have any common
factors, the rank condition is satisfied with probability 1. Similar to
Section 8.2.3, linear independence implies that $\eta(\mathbf{x})$ cannot be a
factor of $\psi$; there is at least one factor $\psi_{i}$ that is distinct.
Based on Proposition 8, these factors cannot vanish on more than
$|\Gamma_{i}|-|\Gamma_{i}\ominus\Lambda_{i}|$ common samples. ∎

## Acknowledgments

The authors would like to thank Dr. Greg Ongie for the valuable discussion
about the proofs of the results in this paper. The first author would also
like to thank Prof. Theodore Shifrin for his help during the discussion on the
intersection of two surfaces and Mr. Biao Ma for the discussion on Lemma 15
and Corollary 16.

## References

  * [1] Aim@shape, digital shape workbench.  http://www.infra-visionair.eu/.
  * [2] H. K. Aggarwal, M. P. Mani, and M. Jacob, Modl: Model-based deep learning architecture for inverse problems, IEEE transactions on medical imaging, 38 (2018), pp. 394–405.
  * [3] M. Atiyah and I. MacDonald, Introduction To Commutative Algebra, Addison-Wesley series in mathematics, Avalon Publishing, 1994.
  * [4] M. Belkin, P. Niyogi, and V. Sindhwani, Manifold regularization: A geometric framework for learning from labeled and unlabeled examples, Journal of machine learning research, 7 (2006), pp. 2399–2434.
  * [5] O. Bernard, D. Friboulet, P. Thévenaz, and M. Unser, Variational b-spline level-set: a linear filtering approach for fast deformable model evolution, IEEE Transactions on Image Processing, 18 (2009), pp. 1179–1191.
  * [6] S. Campagna, L. Kobbelt, and H.-P. Seidel, Directed edges—a scalable representation for triangle meshes, Journal of Graphics tools, 3 (1998), pp. 1–11.
  * [7] T.-H. Chan, K. Jia, S. Gao, J. Lu, Z. Zeng, and Y. Ma, Pcanet: A simple deep learning baseline for image classification?, IEEE Trans. Image Process., 24 (2015), pp. 5017–5032.
  * [8] S. Chen, S. Billings, and P. Grant, Non-linear system identification using neural networks, International journal of control, 51 (1990), pp. 1191–1214.
  * [9] Y. Cho and L. K. Saul, Kernel methods for deep learning, in Advances in neural information processing systems, 2009, pp. 342–350.
  * [10] C. Fefferman, S. Mitter, and H. Narayanan, Testing the manifold hypothesis, Journal of the American Mathematical Society, 29 (2016), pp. 983–1049.
  * [11] M. Fornasier, H. Rauhut, and R. Ward, Low-rank matrix recovery via iteratively reweighted least squares minimization, SIAM Journal on Optimization, 21 (2011), pp. 1614–1640.
  * [12] A. Frome, Y. Singer, F. Sha, and J. Malik, Learning globally-consistent local distance functions for shape-based image retrieval and classification, in 2007 IEEE 11th International Conference on Computer Vision, IEEE, 2007, pp. 1–8.
  * [13] V. Gallese, The roots of empathy: the shared manifold hypothesis and the neural basis of intersubjectivity, Psychopathology, 36 (2003), pp. 171–180.
  * [14] K. Gedalyahu and Y. C. Eldar, Time-delay estimation from low-rate samples: A union of subspaces approach, IEEE Transactions on Signal Processing, 58 (2010), pp. 3017–3031.
  * [15] I. Goodfellow, J. Pouget-Abadie, M. Mirza, B. Xu, D. Warde-Farley, S. Ozair, A. Courville, and Y. Bengio, Generative adversarial nets, in Advances in neural information processing systems, 2014, pp. 2672–2680.
  * [16] R. C. Gunning and H. Rossi, Analytic functions of several complex variables, vol. 368, American Mathematical Soc., 2009.
  * [17] T. Hastie, R. Tibshirani, J. Friedman, and J. Franklin, The elements of statistical learning: data mining, inference and prediction, The Mathematical Intelligencer, 27 (2005), pp. 83–85.
  * [18] F. Huang, A. Ahuja, D. Downey, Y. Yang, Y. Guo, and A. Yates, Learning representations for weakly supervised natural language processing tasks, Computational Linguistics, 40 (2014), pp. 85–120.
  * [19] C. Koulamas and G. J. Kyparisis, Single-machine and two-machine flowshop scheduling with general learning functions, European Journal of Operational Research, 178 (2007), pp. 402–407.
  * [20] J. G. Kuschewski, S. Hui, and S. H. Zak, Application of feedforward neural networks to dynamical system identification and control, IEEE Transactions on Control Systems Technology, 1 (1993), pp. 37–49.
  * [21] Y. LeCun, C. Cortes, and C. Burges, Mnist handwritten digit database, ATT Labs [Online]. Available: http://yann. lecun. com/exdb/mnist, 2 (2010).
  * [22] C. Li, C. Xu, C. Gui, and M. D. Fox, Distance regularized level set evolution and its application to image segmentation, IEEE transactions on image processing, 19 (2010), pp. 3243–3254.
  * [23] Y. M. Lu and M. N. Do, A theory for sampling signals from a union of subspaces, IEEE transactions on signal processing, 56 (2008), pp. 2334–2345.
  * [24] J. Mairal, P. Koniusz, Z. Harchaoui, and C. Schmid, Convolutional kernel networks, in Advances in neural information processing systems, 2014, pp. 2627–2635.
  * [25] M. Mishali, Y. C. Eldar, and A. J. Elron, Xampling: Signal acquisition and processing in union of subspaces, IEEE Transactions on Signal Processing, 59 (2011), pp. 4719–4734.
  * [26] K. Mohan and M. Fazel, Iterative reweighted algorithms for matrix rank minimization, The Journal of Machine Learning Research, 13 (2012), pp. 3441–3473.
  * [27] G. Ongie and M. Jacob, Off-the-grid recovery of piecewise constant images from few fourier samples, SIAM Journal on Imaging Sciences, 9 (2016), pp. 1004–1041.
  * [28] H. Pan, T. Blu, and P. L. Dragotti, Sampling curves with finite rate of innovation, IEEE Transactions on Signal Processing, 62 (2013), pp. 458–471.
  * [29] G. Pillonetto, F. Dinuzzo, T. Chen, G. De Nicolao, and L. Ljung, Kernel methods in system identification, machine learning and function estimation: A survey, Automatica, 50 (2014), pp. 657–682.
  * [30] S. Poddar and M. Jacob, Recovery of noisy points on bandlimited surfaces: Kernel methods re-explained, in 2018 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), IEEE, 2018, pp. 4024 – 4028\.
  * [31] , Recovery of point clouds on surfaces: Application to image reconstruction, in Biomedical Imaging (ISBI 2018), 2018 IEEE 15th International Symposium on, IEEE, 2018, pp. 1272–1275.
  * [32] S. Poddar, Y. Q. Mohsin, D. Ansah, B. Thattaliyath, R. Ashwath, and M. Jacob, Manifold recovery using kernel low-rank regularization: application to dynamic imaging, IEEE Transactions on Computational Imaging, 5 (2019), pp. 478–491.
  * [33] D. Potts and G. Steidl, Fourier reconstruction of functions from their nonstandard sampled radon transform, Journal of Fourier Analysis and Applications, 8 (2002), pp. 513–534.
  * [34] M. Rauth and T. Strohmer, Smooth approximation of potential fields from noisy scattered data, Geophysics, 63 (1998), pp. 85–94.
  * [35] M. Rousson and N. Paragios, Shape priors for level set representations, in European Conference on Computer Vision, Springer, 2002, pp. 78–92.
  * [36] S. T. Roweis and L. K. Saul, Nonlinear dimensionality reduction by locally linear embedding, science, 290 (2000), pp. 2323–2326.
  * [37] H. Sak, A. Senior, K. Rao, O. Irsoy, A. Graves, F. Beaufays, and J. Schalkwyk, Learning acoustic frame labeling for speech recognition with recurrent neural networks, in 2015 IEEE international conference on acoustics, speech and signal processing (ICASSP), IEEE, 2015, pp. 4280–4284.
  * [38] B. Schölkopf and A. J. Smola, Learning with kernels: support vector machines, regularization, optimization, and beyond, MIT press, 2002.
  * [39] F. Schroff, A. Criminisi, and A. Zisserman, Object class segmentation using random forests., in BMVC, 2008, pp. 1–10.
  * [40] H. Shao, A. Kumar, and P. Thomas Fletcher, The riemannian geometry of deep generative models, in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops, 2018, pp. 315–323.
  * [41] T. Sørevik and M. A. Nome, Trigonometric interpolation on lattice grids, BIT Numerical Mathematics, 56 (2016), pp. 341–356.
  * [42] T. Strohmer, Computationally attractive reconstruction of bandlimited images from irregular samples, IEEE Transactions on image processing, 6 (1997), pp. 540–548.
  * [43] T. Strohmer, T. Binder, and M. Sussner, How to recover smooth object boundaries in noisy medical images, in Proceedings of 3rd IEEE International Conference on Image Processing, vol. 1, IEEE, 1996, pp. 331–334.
  * [44] S. Wang and M. Y. Wang, Radial basis functions and level set method for structural topology optimization, International journal for numerical methods in engineering, 65 (2006), pp. 2060–2090.
  * [45] F. Weisz, Summability of multi-dimensional trigonometric fourier series, Surv. Approx. Theory, 7 (2012), pp. 1–179.
  * [46] J. Yu, C. Wojtan, G. Turk, and C. Yap, Explicit mesh surfaces for particle based fluids, in Computer Graphics Forum, vol. 31, Wiley Online Library, 2012, pp. 815–824.
  * [47] J. Zhang, G. G. Walter, Y. Miao, and W. N. W. Lee, Wavelet neural networks for function learning, IEEE transactions on Signal Processing, 43 (1995), pp. 1485–1497.
  * [48] K. Zhang, W. Zuo, Y. Chen, D. Meng, and L. Zhang, Beyond a gaussian denoiser: Residual learning of deep cnn for image denoising, IEEE Trans. Image Process., 26 (2016), pp. 3142–3155.
  * [49] Z. Zhu and H. Leung, Adaptive identification of nonlinear systems with application to chaotic communications, IEEE Transactions on Circuits and Systems I: Fundamental Theory and Applications, 47 (2000), pp. 1072–1080.
  * [50] Q. Zou, S. Poddar, and M. Jacob, Sampling of planar curves: Theory and fast algorithms, IEEE Transactions on Signal Processing, 67 (2019), pp. 6455–6467.