#  Approximation of Functions on Manifolds in High Dimension from Noisy
Scattered Data

Shira Faigenbaum-Golovin1,∗ David Levin1
1 School of Mathematical Sciences, Tel Aviv University, Israel
∗ Corresponding author, E-mail address<EMAIL_ADDRESS>
###### Abstract

In this paper, we consider the fundamental problem of approximation of
functions on a low-dimensional manifold embedded in a high-dimensional space,
with noise affecting both in the data and values of the functions. Due to the
curse of dimensionality, as well as to the presence of noise, the classical
approximation methods applicable in low dimensions are less effective in the
high-dimensional case. We propose a new approximation method that leverages
the advantages of the Manifold Locally Optimal Projection (MLOP) method [14]
and the strengths of the method of Radial Basis Functions (RBF) [13]. The
method is parametrization free, requires no knowledge regarding the manifold’s
intrinsic dimension, can handle noise and outliers in both the function values
and in the location of the data, and is applied directly in the high
dimensions. We show that the complexity of the method is linear in the
dimension of the manifold and squared-logarithmic in the dimension of the
codomain of the function. Subsequently, we demonstrate the effectiveness of
our approach by considering different manifold topologies and show the
robustness of the method to various noise levels.

keywords: Manifold learning, Approximation of functions, High dimensions,
Dimensional reduction, Noisy data

MSC classification: 65D99
(Numerical analysis - Numerical approximation and computational geometry)

##  1 Introduction

In this paper, we consider the following formulation of the problem of
approximation of functions in a high-dimensional space. Let
$f:\mathbb{R}^{n}\to\mathbb{R}^{s}$ be a function, and let
$\\{f(x_{i})\\}_{i=1}^{K}$ be its values on a given sample set of points
$\\{x_{i}\\}_{i=1}^{K}\subset\mathbb{R}^{n}$ with noise present both in the
domain of the function and in its codomain. The goal of the approximation is
to estimate the values of the function at a new set of points. While in low
dimensions numerous methods were suggested to solve this problem (e.g.,
splines, or Moving Least-Squares [18]), in high dimensions this is a
challenging task due to the presence of noise and the curse of dimensionality.
For instance, with respect to the latter challenge, if one merely assumes that
the function is smooth, then approximation rates deteriorate severely with the
growth of the dimension, the reason being that the amount of sampled data
should grow exponentially with respect to the dimension if one wishes to
maintain the same order of approximation.

We categorize high-dimensional approximation methods according to whether the
domain of the function to be approximated is a manifold or not. If no
assumptions are made on the data domain, several methods were suggested. For
example, solutions which treat non-smooth multivariate functions, [1], are
based on sparse occupancy trees [5], Radial Basis Functions [13] (which we
will discuss in detail below), or address the problem in the case where the
values of the function lie on a manifold [16].

In many situations, the high-dimensional data reside on a low-dimensional
manifold, and this information can be exploited to improve the approximation
via one of the following two approaches: approximating in low dimension after
dimension reduction, or alternatively approximating in high dimension. At
times, reducing the dimension (e.g., in PCA [23], Multidimensional Scaling
[11], Linear Discriminant Analysis [15], Locality Preserving Projections [17],
Locally Linear Embedding [24], ISOMAP [27], Diffusion Maps [9], and Neural
Networks in their general form, [21]) can lead to a better approximation (in
terms of handling the challenge of the dimensionality, as well as the noise in
the data). However, it may be non-efficient if the data volume is very large,
and in addition, may result in information loss (due to some assumptions that
need to be made on the data, e.g., regarding the data geometry, or the
intrinsic dimension).

On the other hand, the assumption that the data reside on a manifold can be
utilized to improve the approximation in high dimensions. Approximation of
functions on manifolds is studied using local polynomials [4], wavelets [10],
local linear regression [4], or neural networks [2, 25, 8]. For smooth
functions on $[0,1]^{N}$ which depend on a much smaller number $l$ of
variables, a solution was suggested in [12]. In addition, a recent paper,
[26], proposed a solution based on Moving Least-Squares (MLS) that was
designed to deal with noisy data with good rates of approximation.

In this paper, we propose a method of approximation of functions that
leverages the advantages of the Manifold Locally Optimal Projection (MLOP)
algorithm [14] to complement the strengths of the method of Radial Basis
Functions (RBF) [13, 6]. We introduce this duet for approximation in high
dimensions under noisy conditions (both in the domain and in the codomain of
the function). In what follows we will provide a short introduction to the RBF
method as well as to the Locally Weighted Average Approximation method (which
will be used to contrast the RBF numerical results). In the next section, we
will explain the proposed methodology, and demonstrate the improvement in
approximation via several numerical examples.

##  2 Preliminaries

###  2.1 Preliminaries — the MLOP framework

The Locally Optimal Projection (LOP) method was introduced in [22] to
approximate two-dimensional surfaces in $\mathbb{R}^{3}$ from point-set data.
The procedure does not require the estimation of local normals and planes, or
parametric representations. Its main advantage is that it performs well in the
case of noisy samples. In [14] the LOP mechanism was generalized to devise the
so-called Manifold Locally Optimal Projection (MLOP) method. Here, we give a
concise overview of the MLOP method and its key properties.

First, we introduce the $h$-$\rho$ condition, defined for scattered-data
approximation of functions (which is an adaptation of the condition in [18]
for low-dimensional data), to handle finite discrete data on manifolds.

######  Definition 1.

$h$-$\rho$ sets of fill distance $h$ and density $\leq\rho$ with respect to a
manifold $\mathcal{M}$. Let $\mathcal{M}$ be a manifold in $\mathbb{R}^{n}$
and consider sets of data points sampled from $\mathcal{M}$. We say that such
a set $P=\\{{P_{j}}\\}_{i=1}^{J}$ is an $h$-$\rho$ set if:
1\. $h$ is the fill distance, i.e., $h=\text{median}_{p_{j}\in
P}\mathop{}\\!\mathrm{min}_{p_{j}\in P\backslash\\{p_{i}\\}}\|p_{i}-p_{j}\|$.
2\. $\\#\\{P\cap\bar{B}(y,kh)\\}\leq\rho k^{n},\quad k\geq 1,\quad
y\in\mathbb{R}^{n}$.
Here $\\#Y$ denotes the number of elements in a set $Y$ and $\bar{B}(x,r)$
denotes the closed ball of radius $r$ centered at $x$.

Note that the last condition regarding the point separation $\delta$ defined
in [18], which states that there exists $\delta>0$ such that
$\|p_{i}-p_{j}\|\geq\delta,\quad 1\leq i\leq j\leq J$, is redundant in the
case of finite data. We also note that the vanilla definition of the fill
distance uses the supremum $\sup$ in its expression; here we use the median in
order to deal with the presence of outliers.

The setting of the high-dimensional reconstruction problem is the following:
Let $\mathcal{M}$ be a manifold in $\mathbb{R}^{n}$ of unknown intrinsic
dimension $d\ll n$. There is given a noisy point-cloud
$P=\\{p_{j}\\}_{j=1}^{J}\subset\mathbb{R}^{n}$ situated near the manifold
$\mathcal{M}$ such that $P$ is a $h$-$\rho$ set. We wish to find a new point-
set $Q=\\{q_{i}\\}_{i=1}^{I}\subset\mathbb{R}^{n}$ which will serve as a
noise-free approximation of $\mathcal{M}$. We seek a solution in the form of a
new, quasi-uniformly distributed point-set $Q$ that will replace the given
data $P$ and provide a noise-free approximation of $\mathcal{M}$. This is
achieved by leveraging the well-studied weighted $L_{1}$-median [28] used in
the LOP algorithm and requiring a quasi-uniform distribution of points
$q_{i}\in Q$. These ideas are encoded by the cost function

$G(Q)=E_{1}(P,Q)+\Lambda E_{2}(Q)=\sum\limits_{q_{i}\in
Q}\sum\limits_{p_{j}\in
P}\|q_{i}-p_{j}\|_{H_{\epsilon}}w_{i,j}+\sum\limits_{q_{i}\in
Q}\lambda_{i}\sum\limits_{q_{i^{\prime}}\in
Q\backslash\\{q_{i}\\}}\eta(\|q_{i}-q_{i}^{\prime}\|)\widehat{w}_{i,i^{\prime}}\,,$
(1)

where the weights $w_{i,j}$ are rapidly decreasing smooth functions. The MLOP
implementation uses
$w_{i,j}=\exp\big{\\{}\\!-\|q_{i}-p_{j}\|^{2}/{h_{1}^{2}}\big{\\}}$ and
$\widehat{w}_{i,i^{\prime}}=\exp\big{\\{}\\!-\|q_{i}-q_{i}^{\prime}\|^{2}/{h_{2}^{2}}\big{\\}}$.
The $L_{1}$-norm used in [22] is replaced by the “norm”
$\|\cdot\|_{H_{\epsilon}}$ introduced in [19] as
$\|v\|_{H_{\epsilon}}=\sqrt{v^{2}+\epsilon}$, where $\epsilon>0$ is a fixed
parameter (in our case we take $\epsilon=0.1$). As shown in [19], using
$\|\cdot\|_{H_{\epsilon}}$ instead of $\|\cdot\|_{1}$ has the advantage that
one works with a smooth cost function and outliers can be removed. In
addition, $h_{1}$ and $h_{2}$ are the support size parameters of $w_{i,j}$ and
$\widehat{w}_{i,i^{\prime}}$ which guarantee a sufficient amount of $P$ or $Q$
points for the reconstruction (for more details, see the subsection “Optimal
Neighborhood Selection” in [14]). Also, $\eta(r)$ is a decreasing function
such that $\eta(0)=\infty$; in our case we take $\eta(r)=\frac{1}{3r^{3}}$.
Finally, $\\{\lambda_{i}\\}_{i=1}^{I}$ are constant balancing parameters.

In order to solve the problem with the cost function (1), we look for a point-
set $Q$ that minimizes $G(Q)$. The solution $Q$ is found via the gradient
descent iterations

$q_{i^{\prime}}^{(k+1)}=q_{i^{\prime}}^{(k)}-\gamma_{k}\nabla
G(q_{i^{\prime}}^{(k)}),\qquad i^{\prime}=1,\dots,I\,,$ (2)

where the initial guess $\\{q_{i}^{(0)}\\}_{i-1}^{I}=Q^{(0)}$ consists of
points sampled from $P$.
The gradient of $G$ is

$\nabla
G(q_{i^{\prime}}^{(k)})=\sum\limits_{j=1}^{J}{\big{(}q_{i^{\prime}}^{(k)}-p_{j}\big{)}\alpha_{j}^{i^{\prime}}}-\lambda_{i^{\prime}}\sum\limits_{\begin{subarray}{c}i=1\\\
i\neq
i^{\prime}\end{subarray}}^{I}{\big{(}q_{i^{\prime}}^{(k)}-q_{i}^{(k)}\big{)}\beta_{i}^{i^{\prime}}}\,,$
(3)

with the coefficients $\alpha_{j}^{i^{\prime}}$ and $\beta_{j}^{i^{\prime}}$
given by the formulas

$\alpha_{j}^{i^{\prime}}=\frac{w_{i,j}}{\|q_{i}-p_{j}\|_{H_{\epsilon}}}\left(1-\frac{2}{h_{1}^{2}}\|q_{i}-p_{j}\|_{H_{\epsilon}}^{2}\right)$
(4)

and

$\beta_{i}^{i^{\prime}}=\frac{\widehat{w}_{i,i^{\prime}}}{\|q_{i}-q_{i^{\prime}}\|}\left(\left|{\frac{\partial\eta\left(\|q_{i}-q_{i^{\prime}}\|\right)}{\partial
r}}\right|+\frac{2\eta\left(\|q_{i}-q_{i^{\prime}}\|\right)}{h_{2}^{2}}\|q_{i}-q_{i^{\prime}}\|\right),$
(5)

for $i=1,...,I$, $i\neq i^{\prime}$. In order to balance the two terms in
$\nabla G(q_{i^{\prime}}^{(k)})$, the factors $\lambda_{i^{\prime}}$ are
initialized in the first iteration as

$\lambda_{i^{\prime}}=-\,\frac{\bigg{\|}\sum\limits_{j=1}^{J}{\big{(}q_{i^{\prime}}^{(k)}-p_{j}\big{)}\alpha_{j}^{i^{\prime}}}\bigg{\|}}{\bigg{\|}\sum\limits_{i=1}^{I}{\big{(}q_{i^{\prime}}^{(k)}-q_{i}^{(k)}\big{)}\beta_{i}^{i^{\prime}}}\bigg{\|}}\,.$
(6)

Balancing the contribution of the two terms is important in order to maintain
equal influence of the attraction and repulsion forces in $G(Q)$. The step
size in the direction of the gradient $\gamma_{k}$ is calculated as indicated
in [3]:

$\displaystyle\gamma_{k}=\frac{\langle\bigtriangleup
q_{i^{\prime}}^{(k)},\bigtriangleup
G_{i^{\prime}}^{(k)}\rangle}{\langle\bigtriangleup
G_{i^{\prime}}^{(k)},\bigtriangleup G_{i^{\prime}}^{(k)}\rangle}\,,$ (7)

where $\bigtriangleup
q_{i^{\prime}}^{(k)}=q_{i^{\prime}}^{(k)}-q_{i^{\prime}}^{(k-1)}$ and
$\bigtriangleup G_{i^{\prime}}^{(k)}=\nabla G_{i^{\prime}}^{(k)}-\nabla
G_{i^{\prime}}^{(k-1)}$.

######  Remark 2.1.

The reasoning in terms of Euclidean distances, which is the cornerstone of the
MLOP method, works well in low dimensions, e.g., for the reconstruction of
surfaces in 3D, but breaks down in high dimensions once noise is present. To
deal with this issue, a dimension reduction is performed via random linear
sketching [29]. It should be emphasized that the dimension reduction procedure
is utilized solely for the calculation of norms, and the manifold
reconstruction is performed in the high-dimensional space. Thus, given a point
$x\in\mathbb{R}^{n}$, we project it to a lower dimension $m\ll n$ using a
random matrix, $S$, with certain properties. Subsequently, the norm of
$\|S^{t}x\|$ will approximate $\|x\|$. The construction of $S$ is carried out
in the following steps:

  1. 1.

Sample $G\in\mathbb{R}^{J\times m}$ with $G\sim N(0,1)$.

  2. 2.

Compute $B\in\mathbb{R}^{n\times m}$ as $B:=P^{\rm t}G$.

  3. 3.

Calculate the QR decomposition of $B$ as $B=SR$, and use $S$ as the dimension
reduction matrix.

Preliminaries — Optimal Neighborhood Selection

The support sizes $h_{1}$ and $h_{2}$ of the functions
$w_{i,j}=\exp\big{\\{}\\!-\|q_{i}-p_{j}\|^{2}/{h_{1}^{2}}\big{\\}}$ and
$\widehat{w}_{i,i^{\prime}}=\exp\big{\\{}\\!-\|q_{i}-q_{i}^{\prime}\|^{2}/{h_{2}^{2}}\big{\\}}$,
respectively, are closely related to the fill distance of the $P$-points and
the $Q$\- points. Due to the importance of optimal selection of these
parameters, we quote here several definitions and results from [14]. Unlike
the standard definition of the notion of fill distance in scattered-data
function approximation [18], we introduce

######  Definition 2.

The fill distance of the set $P$ is

$\displaystyle h_{0}=\text{median}_{p_{i}\in
P}\mathop{}\\!\mathrm{min}_{p_{j}\in
P\backslash\\{p_{i}\\}}\|p_{i}-p_{j}\|\,.$ (8)

Note that the vanilla definition of fill distance uses the supremum in the
definition (instead of the median). However, as mentioned above, in our case
we replace the supremum with the median so as to deal with the presence of
outliers.

######  Definition 3.

Given two point-clouds, $P=\\{p_{j}\\}_{j=1}^{J}\subset\mathbb{R}^{n}$ and
$Q=\\{q_{i}\\}_{i=1}^{I}\subset\mathbb{R}^{n}$, situated near a manifold
$\mathcal{M}$ in $\mathbb{R}^{n}$, such that their sizes obey the constraint
$I\leq J$, denote $\nu=\left\lfloor\frac{J}{I}\right\rfloor$. Then we say that
the radius that guarantees approximately $\nu$ points from $P$ in the support
of each point $q_{i}$ is $\widehat{h}_{0}=c_{1}h_{0}$, with $c_{1}$ given by

$\displaystyle\widehat{h}_{0}=c_{1}h_{0},~{}{\rm
with}~{}~{}c_{1}=\text{argmin}\\{c:\\#(\bar{B}_{ch_{0}}(q_{i})\cap
P)\geq\nu,\,\forall q_{i}\in Q\\}\,.$ (9)

######  Remark 2.2.

Let $\sigma$ be the variance of the Gaussian
$w(r)=\exp\\{-{r^{2}}/{\sigma^{2}}\\}$. For the normal distribution, four
standard deviations away from the mean account for $99.99\%$ of the set. In
our case, by the definition of $w_{i,k}$, since $h$ is the square root of the
variance, $4\sigma=4\frac{h}{\sqrt{2}}=2\sqrt{2}h_{1}$ covers $99.99\%$ of the
support size of $w_{i,k}$.

The following theorem, proved in [14], indicates how the parameters $h_{1}$
and $h_{2}$ should be selected.

######  Theorem 2.3.

Let $\mathcal{M}$ be a $d$-dimensional manifold in $\mathbb{R}^{n}$. Suppose
given two point-clouds, $P=\\{p_{j}\\}_{j=1}^{J}\subset\mathbb{R}^{n}$ and
$Q=\\{q_{i}\\}_{i=1}^{I}\subset\mathbb{R}^{n}$, situated near $\mathcal{M}$,
such that their sizes obey the constraint $I\leq J$, and let
$\nu=\left\lfloor\frac{J}{I}\right\rfloor$. Let $w_{i,j}$ be the locally
supported weight function given by
$w_{i,j}=\exp\\{-{\|q_{i}-p_{j}\|^{2}}/{h^{2}}\\}$. Then a neighborhood size
of $h=2\sqrt{2}\widehat{h}_{0}$ guarantees $2^{1.5d}\nu$ points in the support
of $w_{i,j}$, where $\widehat{h}_{0}=c_{1}h_{0}$, with $c_{1}$ given by (9).

Theoretical Analysis of the MLOP Method

For the sake of completeness, we mention here several important results
regarding the convergence of the MLOP method, its order of approximation, rate
of convergence and complexity (for more details, see [14]).

######  Theorem 2.4 (Convergence to a stationary point).

Let $\mathcal{M}$ be a in $\mathbb{R}^{n}$ of unknown intrinsic dimension $d$.
Suppose that the scattered data points $P=\\{{p_{j}}\\}_{j=1}^{J}$ were
sampled near the manifold $\mathcal{M}$, $h_{1}$ and $h_{2}$ are set as in
Theorem 2.3, and the $h$-$\rho$ set condition is satisfied with respect to
$\mathcal{M}$. Let the points $Q^{(0)}=\\{q_{i}^{(0)}\\}_{i=1}^{I}$ be sampled
from $P$. Then the gradient descent iterations (6) converge almost surely to a
local minimizer $Q^{*}$.

######  Theorem 2.5 (Order of approximation).

Let $P=\\{p_{j}\\}_{j=1}^{J}$ be a set of points that are sampled (without
noise) from a $d$–dimensional $C^{2}$ manifold $\mathcal{M}$, and satisfy the
$h$-$\rho$ condition. Then for a fixed $\rho$ and a finite support of size $h$
of the weight functions $w_{i,j}$, the set $Q$ defined by the MLOP algorithm
has an order of approximation $O(h^{2})$ to $\mathcal{M}$.

######  Definition 4.

A differentiable function $f(\cdot)$ is called $L$-smooth if for any
$x_{1},x_{2}$

$\displaystyle\|\nabla f(x_{1})-\nabla f(x_{2})\|\leq L\|x_{1}-x_{2}\|\,.$

######  Theorem 2.6 (Rate of convergence).

Suppose the point-set $P=\\{{p_{j}}\\}_{j=1}^{J}$ is sampled near a
$d$-dimensional manifold in $\mathbb{R}^{n}$ and the assumptions of Theorem
2.4 are satisfied. Suppose the cost function $G$ defined in (1) is an
$L$-smooth function. For $\epsilon>0$, let $Q^{*}$ be a local fixed-point
solution of the gradient descent iterations, with step size
$\gamma=\epsilon^{-1}$. Set the termination condition as $\|\nabla
f(x)\|\leq\epsilon$. Then $Q^{*}$ is an $\epsilon$-first-order stationary
point that will be reached after $k=L(G(Q^{(0)})–G(Q^{*}))\epsilon^{-2}$
iterations, where $L=l^{2}$ and $l<\infty$ is a bounded parameter.

######  Theorem 2.7 (Complexity).

Given a point-set $P=\\{{p_{j}}\\}_{j=1}^{J}$ sampled near a $d$-dimensional
manifold $\mathcal{M}\in\mathbb{R}^{n}$, let $Q=\\{q_{i}\\}_{i=1}^{I}$ be a
set of points that will provide the desired manifold reconstruction. Then the
complexity of the MLOP algorithm is $O(nmJ+kI(nm\widehat{I}+\widehat{J}))$,
where the number of iterations $k$ is bounded as in Theorem 2.6, $m\ll n$ is
the smaller dimension to which we reduce the dimension of the data, and
$\widehat{I}$ and $\widehat{J}$ are the numbers of points in the support of
the weight functions $\widehat{w}_{i,i^{\prime}}$, $w_{i,j}$ that belong to
the $Q$-set and $P$-set, respectively. Thus, the approximation is linear in
the ambient dimension $n$, and does not depend on the intrinsic dimension $d$.

###  2.2 Preliminaries — Radial Basis Functions

Radial Basis Functions constitute a very useful and convenient multivariate
interpolation tool [13, 6]. Given the values of a function
$f:\mathbb{R}^{n}\to\mathbb{R}^{s}$ at center points $x_{i}\in\mathbb{R}^{n}$,
$i=1,...,K$, we approximate the value of $f$ at a new point $x$ by the formula

$\widetilde{f}(x)=\sum\limits_{i=1}^{K}\lambda_{i}\phi(\|x-x_{i}\|)\,,$ (10)

where $\phi:\mathbb{R}_{+}\to\mathbb{R}$ is a radial basis function and
$\lambda_{i}$ are scalar parameters chosen to maintain interpolation at the
center points, i.e., $\widetilde{f}(x_{i})=f(x_{i})$. For examples of possible
choices of radial basis functions, see [7, 30]. In the numerical examples
presented below we choose to use the following Gaussian RBF with local
support:

$\displaystyle\phi_{1}(r)=\exp\\{-({r}/{h})^{2}\\},$
$\displaystyle\phi_{2}(r)=\exp\\{-({r}/{h})^{2}\\}\big{(}1+r/h\big{)},$
$\displaystyle\phi_{3}(r)=\exp\\{-({r}/{h})^{2}\\}\big{(}15+15(r/h)+6(r/h)^{2}+(r/h)^{3}\big{)}.$

In what follows, we will state a theorem proved in [31] on the order of
approximation in a Sobolev space of a method that uses a RBF of general form,
namely

$\widetilde{f}(x)=\sum\limits_{j=1}^{K}\lambda_{j}\phi_{w}(\|x-x_{j}\|)+\sum\limits_{i=1}^{I}\alpha_{i}p_{i}(x)\,,$
(11)

where $\phi_{w}:=\phi(\cdot/w)$, with $w$ depending on the fill-distance $h$
of the points $X=\\{x_{j}\\}_{j=1}^{K}$, $p_{1},...,p_{I}$ is a polynomial
basis for $\Pi_{m}$ (the space of polynomials of degree $\leq m$), and the
coefficients $\lambda_{j}$ and $\alpha_{j}$ are chosen to satisfy the linear
system $\widehat{f}(x_{j})=f(x_{j})$ for $j=1,...,K$ and
$\sum\limits_{j=1}^{K}\lambda_{j}p_{i}(x_{j})=0$, $i=1,...,I$.

To state the theorem on the order of approximation of the RBF method we nee
some additional definitions.

######  Definition 5.

For $k\in\mathbb{N}$ $p$, the Sobolev space $W_{p}^{k}(\Omega)$ is defined as

$W_{p}^{k}(\Omega):=\bigg{\\{}f:\|f\|_{k,L_{p}(\Omega)}:=\bigg{(}\sum_{|\alpha|_{1}\leq
k}\|D^{\alpha}f\|_{L_{p}(\Omega)}^{p}\bigg{)}^{1/p}<\infty\bigg{\\}},$

for $p<\infty$, and as

$W_{\infty}^{k}(\Omega):=\bigg{\\{}f:\|f\|_{k,L_{\infty}(\Omega)}:=\sum_{|\alpha|_{1}\leq
k}\|D^{\alpha}f\|_{L_{\infty}(\Omega)}<\infty\bigg{\\}},$

for $p=\infty$.

######  Definition 6.

Let $X=\\{x_{j}\\}_{j=1}^{K}$ be a set of points with fill-distance $h$ and
separation distance $\delta=\underset{1\leq i\neq j\leq
N}{\mathop{}\\!\mathrm{min}}\|x_{i}-x_{j}\|/2$. Then we say that $X$ is quasi-
uniformly distributed if there exists a constant $\eta>0$ independent of $X$
such that

$2\delta\leq h\leq\eta\delta\,.$ (11)

######  Definition 7.

Let $\phi_{w}$ be a radial basis function, and let $\widehat{\phi}_{w}$ be its
Fourier transforms. Define the supremum of the norm of $\widehat{\phi}_{w}$ as

$M_{\phi,w}(r):=\underset{\theta\in
B(0,r)}{\sup}\|\widehat{\phi}_{w}(\theta)\|^{-1/2}\,.$ (12)

######  Definition 8.

Let $f$ be as defined in (11), and $\phi_{w}$ be a radial basis function. Then
the norm of the corresponding error functional is defined as

$P_{\phi,X}(x)=\underset{\|f\|_{\phi}\neq
0}{\sup}\frac{\|f(x)-\widetilde{f}(x)\|}{\|f\|_{\phi}},\quad{\rm
where}~{}\|f\|_{\phi}=\int_{\mathbb{R}^{d}}\frac{\|\widehat{f}(\theta)\|^{2}}{\widehat{\phi}_{w}(\theta)}d\theta\,.$
(13)

Finally, we are ready to state the promised theorem, proven in [31]

######  Theorem 2.8.

Let $X=\\{x_{j}\\}_{j=1}^{K}$ be a set of quasi-uniformly distributed
scattered points (see Definition 6), and let $\widetilde{f}(x)$, defined as in
(10), be an interpolant to $f$ on $X$ using the radial basis function
$\phi_{w}=\phi(\cdot/w)$. Let $M_{\phi,w}(r)$ with $r>0$, be defined as in
(12). Assume that there exists a constant $\delta_{0}>0$ such that

$P_{\phi,X/w}(x/w)M_{\phi,w}(\delta_{0}/h)\leq o(h^{k})\,.$ (14)

Then, for every function $f\in W_{\infty}^{k}(\Omega)$, with $k\in\mathbb{N}$,
the error of the RBF method is estimated as

$\|f-\widetilde{f}\|_{L_{\infty}(\Omega)}=o(h^{k})\,.$ (15)

######  Remark 2.9.

As a result, an important key advantage of the RBF method is that it performs
better on quasi-uniform samples. This property will be utilized in the next
section.

We will also give a short introduction to the Locally Weighted Average
Approximation, which will be used as a reference in the section devoted to our
numerical examples. Given the values of the function
$f:\mathbb{R}^{n}\to\mathbb{R}^{s}$ at the points $x_{i}\in\mathbb{R}^{n}$,
the locally weighted average approximation of $f$ at a point $x$ is defined as

$\displaystyle f(x)=\frac{\sum_{i}w_{i}f(x_{i})}{\sum_{i}w_{i}}\,,$

where $w_{i}=\exp\\{-{\|x_{i}-x\|^{2}}/{h^{2}}\\}$, and $h$ is the fill-
distance of the points $\\{x_{i}\\}_{i=1}^{I}$. Concerning the accuracy of the
method, we note that the locally weighted average approximation reconstructs
constant functions.

##  3 Extending the MLOP Method to Approximation of Functions on a Manifold

Let $\mathcal{M}$ be a smooth $d$-dimensional manifold in $\mathbb{R}^{n}$,
where $d\ll n$. Let $P=\\{p_{j}\\}_{j=1}^{J}$ be a set of points which were
sampled from $\mathcal{M}$ and are affected by noise, and let
$f:\mathcal{M}\to\mathbb{R}^{s}$ be a smooth function. Given noisy
measurements of $f$ at the points in $P$, the approximation problem has two
steps:

  1. 1.

Find a noise-free representation of the manifold $\mathcal{M}$ and the noise-
free values of $f$.

  2. 2.

Estimate the value of the function $f$ at a new given point $x$.

Our solution for the first step is based on generalizing the MLOP method
designed for manifold denoising to the case of function denoising. The key
idea of the solution is to embed the approximation problem in a higher-
dimensional space, and denoise the data there. Given the input data, which
consists of the set of points $\\{p_{j}\\}_{j=1}^{J}\subset\mathbb{R}^{n}$ and
the set of values $\\{f(p_{j})\\}_{j=1}^{J}\subset\mathbb{R}^{s}$, we define a
new point-set $\widehat{P}=\\{\widehat{p}_{j}\\}_{j=1}^{J}$ to be the graph of
a function $f$, i.e. as the set of ordered pairs, where
$\widehat{p_{j}}=(p_{j},f(p_{j}))$ is the pairing of $p_{j}$ with the value
$f(p_{j})$. The points in $\widehat{P}$ are now considered as data points in
$\mathbb{R}^{n+s}$, taken from the $\widehat{\mathcal{M}}={\rm
Graph}\,f=\\{(x,f(x)):\,x\in\mathcal{M}\\}$. The newly defined set
$\widehat{\mathcal{M}}$, being the graph of a smooth function $f$ defined on a
smooth $d$-dimensional manifold, is itself a smooth $d$-dimensional manifold.
It should be noted that prior to the embedding in the $(n+s)$-dimensional
space, the values of $f$ should be normalized to the maximum value of the
$p_{j}$ coordinates, in order to avoid them dominating the $p_{j}$ entries
during the norm calculations of the MLOP algorithm.

In this setting, we are now denoising a $d$-dimensional manifold embedded in
$n+s$ dimensions. We apply the MLOP method on the new data set $\widehat{P}$,
and look for a clean dataset
$\widehat{Q}=\\{\widehat{q}_{i}\\}_{i=1}^{I}\subset\mathbb{R}^{n+s}$ which
will serve as a noise-free approximation of $\widehat{\mathcal{M}}$. The main
advantage of this approach is that with a single MLOP execution on the
$\\{\widehat{p_{j}}=(p_{j},f(p_{j}))\\}$ data in $\mathbb{R}^{n+s}$ we produce
a noise-free set $\widehat{Q}\subset\mathbb{R}^{n+s}$, which in fact consist
of a noise-free set $Q$, that reconstructs the manifold $\mathcal{M}$, and an
estimate of the clean value of the function evaluated at these points,
$\widetilde{f}(Q)$.

In the second step, we address the problem of evaluating the function at a new
point $z\in\mathcal{M}$. The outcome of the first step is a set of points
which is not only noise-free, but also quasi-uniformly distributed on the
manifold $\widehat{\mathcal{M}}$. This key idea paves the way towards
estimating the value of the function at a new given point on $\mathcal{M}$, or
near $\mathcal{M}$, with a good order of approximation. Our solution is based
on utilizing the RBF approximation, as defined in (10), while setting the
centers at the cleaned points $Q$ and appending the corresponding cleaned
function values $\widetilde{f}(Q)$. These steps are summarized in Algorithm 1:

Algorithm 1 Function Approximation on a Manifold in High Dimensions

1:Input: $P=\\{p_{j}\\}_{j=1}^{J}\subset\mathbb{R}^{n}$, $\\{f(p_{j})\\}$,
$\\{z_{k}\\}_{k=1}^{K}$ $\triangleright$ where $\\{z_{k}\\}_{k=1}^{K}$ is a
set of new points for function approximation

2:Output: $Q=\\{q_{i}\\}_{i=1}^{I}\subset\mathbb{R}^{n}$, $\widetilde{f}(Q)$,
$\\{\widetilde{f}(z_{k})\\}_{k=1}^{K}$

3:Denoise the input data by running MLOP with
$\widehat{P}=(P,f(P))\subset\mathbb{R}^{n+s}$
$\to\widehat{Q}=(Q,\widetilde{f}(Q))$

4:for each $z_{k}\in\\{z_{k}\\}_{k=1}^{K}$ do

5: Approximate $\widetilde{f}(z_{k})$ via RBF, with centers set at $Q$ and the
corresponding $\widetilde{f}(Q)$ values

6:end for

##  4 Theoretical Analysis of the Method

In this section we discuss some of the theoretical aspects of the proposed
approach to approximation of functions. Our analysis relies on the theory of
the MLOP method as well as the RBF method. Specifically, we can use Theorems
2.4, 2.5, 2.6, 2.7 regarding the convergence of the MLOP algorithm to a
stationary point, its rate of convergence, and its complexity, estimated for
the problem at hand. We also build on the results about the order of
approximation of the MLOP. Thus, we can state the following theorem on the
order of approximation for the our approximation problem.

###  4.1 Order of Approximation

######  Theorem 4.1 (Order of approximation).

Let $P=\\{p_{j}\\}_{j=1}^{J}$ be a set of points sampled from a
$d$-dimensional $C^{2}$ manifold $\mathcal{M}$ without noise that satisfy the
$h$-$\rho$ condition. Let $f:\mathcal{M}\to\mathbb{R}^{s}$ be a smooth
multivariate function given at points of $P$. Suppose that $f\in
W_{\infty}^{k}(\Omega)$ and fulfills all the conditions of Theorem 2.8. Then:

  1. (i)

For fixed $\rho$ and $\delta$, there is a set of points
${Q}=\\{{q}_{i}\\}_{i=1}^{I}\subset\mathbb{R}^{n}$ which approximates
$\mathcal{M}$ with order $O(\widehat{h}^{2})$, where
$\widehat{h}=\max\\{\widehat{h}_{1},\widehat{h}_{2}\\}$ and $\widehat{h}_{1}$
and $\widehat{h}_{2}$ are as in Definition 3 with respect to the high-
dimensional data $\widehat{P}=\\{(p_{j},f(p_{j}))\\}\subset\mathbb{R}^{n+s}$.
Moreover, the order of approximation of $f$ on the set $Q$ is also
$O(\widehat{h}^{2})$.

  2. (ii)

The order of approximation of $f$ at a new point is less than
$C_{1}\widehat{h}^{2}+C_{2}h_{2}^{k}$ when using the RBF approximation with
centers at $Q$, and $\widetilde{f}(Q)$, where $h_{2}$ is the optimal support
of $\widehat{w}_{i,i^{\prime}}$ as in Definition 3, and $C_{1}$ and $C_{2}$
are constants.

###### Proof.

Given the input data, which consists of the points $\\{p_{j}\\}_{j=1}^{J}$ and
the function values $\\{f(p_{j})\\}_{j=1}^{J}$, we define a new point-set,
which is the graph of the function $f$,
$\widehat{P}=\\{\widehat{p}_{j}\\}_{j=1}^{J}$ sampled from a new manifold
$\widehat{\mathcal{M}}\in\mathbb{R}^{n+s}$, where $\widehat{p_{j}}$ are
defined as the pairing of the $p_{j}$-data with the corresponding values
$f(p_{j})$. In this setting, we use the MLOP algorithm to denoise a
$d$-dimensional manifold embedded in $(n+s)$-dimension. Thus, by applying the
MLOP method to the new data we obtain the set $\widehat{Q}$ of points that
reconstruct $\widehat{\mathcal{M}}$. From Theorem 2.5 it follows that
$\widehat{Q}$ approximates $\widehat{\mathcal{M}}$ with the order
$O(\widehat{h}^{2})$, where $\widehat{h}$ is the representative distance
introduced in Definition 3 for the extended data set
$\widehat{P}=\\{(p_{j},f(p_{j}))\\}\in\mathbb{R}^{n+s}$. Recall that
$\widehat{Q}$ is in fact a combination of a noise-free set $Q$, which
reconstructs the manifold $\mathcal{M}$, and an estimate of the clean values
of the function at these points in $Q$, $\widetilde{f}(Q)$. Therefore, the
order of approximation of $Q$ to $\mathcal{M}$, and of the estimated values
$\widehat{f}$ to $f$ is $O(\widehat{h}^{2})$.

Subsequently, we use an RBF approximation with theoretical order of
approximation $O(h^{k})$ (as stated in Theorem 2.8). In our case the relevant
$h$ is $h_{2}$, the support size of $\widehat{w}_{i,i^{\prime}}$ as in
Definition 3. Therefore, the overall order of approximation for a new point is
a combination of the two orders, namely $\leq
C_{1}\widehat{h}^{2}+C_{2}h_{2}^{k}$, with $C_{1}$, and $C_{2}$ constants. ∎

###  4.2 Complexity of the Approximation of Functions

######  Theorem 4.2.

Let $P=\\{{p_{j}}\\}_{j=1}^{J}$ be a set of points sampled near a
$d$-dimensional manifold $\mathcal{M}\subset\mathbb{R}^{n}$ and let
$Q=\\{q_{i}\\}_{i=1}^{I}$ be a set of points which will provide the desired
noise-free manifold reconstruction. Let $f:\mathcal{M}\to\mathbb{R}^{s}$ be a
multivariate function given at the points of $P$. Then the complexity of the
approximation of $f$ via the MLOP algorithm for the denoising step is
$O((n+s)mJ+kI((n+s)m\widehat{I}+\widehat{J})+nmI+I\,\textup{log}_{2}^{2}s)$,
and for evaluating $f$ at a new point is $O(nm+I)$, where the number of
iterations $k$ is bounded as in Theorem 2.6, $m$, with $m\ll n$, is the
smaller dimension to which we reduce the dimension of the data, and
$\widehat{I}$ and $\widehat{J}$ are the numbers of $Q$-points and $P$-points,
respectively, in the support of the weight function $\hat{w}_{i,i^{\prime}}$
and $w_{i,j}$.

###### Proof.

The estimate of the complexity of the algorithm can be separated into two
steps: pre-processing and evaluating the function at a new point. The pre-
processing step consists of applying the MLOP algorithm in $\mathbb{R}^{n+s}$,
as well as finding the RBF coefficients $\lambda_{i}$ by solving the Least-
Squares problem. By Theorem 2.7, the complexity of applying the MLOP in the
higher dimension $n+s$ is $O((n+s)mJ+kI((n+s)m\widehat{I}+\widehat{J}))$. The
output of this stage is a points set $Q$ of size $I$, and the corresponding
set of values $\widetilde{f}(Q)$. From this point on, all the approximation
operations are performed on $Q$, and if $I\ll J$ then we can increase
sufficiently the efficiency. The next part of the pre-processing step is to
evaluate the radial basis function $\phi$ for each $q_{i}\in Q$, which costs
$O(nmI)$, and then to find the $\lambda_{i}$ by solving the Least-Squares
problem, which takes $O(I\,\textup{log}_{2}^{2}s)$ (as shown in [20]). As a
result, the complexity of the pre-processing step is
$O((n+s)mJ+kI((n+s)m\widehat{I}+\widehat{J})+nmI+I\,\textup{log}_{2}^{2}s)$.

Finally, using the $\lambda_{i}$ already found, we evaluate the function at a
new point in time $O(nm+I)$. It should be stressed that although the pre-
processing steps are cost-effective, they are executed once before the
function is approximated at a new points set. Thus, if the number of new
points for which the approximation needs to be found is large, then the pre-
processing steps have less effect on the runtime. ∎

##  5 Numerical Examples

In what follows we present several numerical experiments to demonstrate the
advantages of our methodology. We can point out two strengths of the proposed
approximation approach. On the one hand, denoising the data domain as well as
the function codomain plays an important role in the approximation of
functions. On the other hand, sampling the manifold quasi-uniformly improves
significantly the approximation of functions by means of classical
approximation methods on new data.

Given data sampled from a manifold with noise, and the noisy values of a
function $f$ at these points, we follow the function approximation procedure
described above. Specifically, we define a new problem in $\mathbb{R}^{n+s}$,
and apply $k$ MLOP iterations to clean the newly defined manifold, which
results in a new point-set $Q^{(k)}$, and the corresponding cleaned value set
$\widetilde{f}(Q^{(k)})$. Next, we randomly select 100 points,
$\\{z_{i}\\}_{i=1}^{100}$, from a clean reference dataset, and estimate the
values of the function $f$ at these points using both the RBF approximation,
where the centers of the RBF function are taken at the points of $Q^{(k)}$
(with the radial basis function set to either $\phi_{1}$, $\phi_{2}$, or
$\phi_{3}$, as defined below equation (10)), and the locally weighted average
approximation (defined by formula 2.2). As a result, the approximation at the
new points relies on the clean quasi-uniformly distributed $Q^{(k)}$ points,
as well as on the clean values $\widetilde{f}(Q^{(k)})$. In all the stages
above we evaluate the accuracy of the approximation as the relative maximum
error of the $L_{1}$ norm of the difference between the value of
$\widetilde{f}$ at the new point and the value of $f$ at the closest point in
the reference dataset, as well as the root-mean-square error and the standard
deviation.

We start with two examples of functions, one smooth and the other non-smooth,
both on a one-dimensional manifold embedded in a high-dimensional space.
Although in principle the approximation requires a smooth function, it can
still be applied to a non-smooth function, provided that we end up with a
smoothed result. Specifically, we consider the case of the manifold $O(2)$ of
orthogonal matrices, embedded in a 60-dimensional linear space by using the
parameterization

$p=[\cos(\theta),-\sin(\theta),\sin(\theta),\cos(\theta),0,\dotsc,0]\,,$ (16)

where $\theta\in[-\pi,\pi]$. The input dataset $\widehat{P}$ was constructed
by sampling 500 equally distributed points in the parameter space. Next, we
randomly sampled an orthogonal matrix $A\in\mathbb{R}^{60\times 60}$, and
created a new point-set via the non-trivial vector embedding

$P=A\widehat{P}\,.$ (17)

Subsequently, we added a uniform noise $U(-0.1,0.1)$, and initialized the set
$Q$ by selecting $55$ points from $P$. Figure 1 (A) illustrates the first two
coordinates of the points in our set (after multiplication by the matrix
$A^{-1}$). The noisy sample points are shown in green, while the initial
reconstruction points are shown in red.

We start by approximating the smooth function
$f(x)=\frac{1}{4}(1+\sin(10\,\theta))$, where $\theta$ corresponds to the
value used in the expression (16) of $p$. Next, a uniform noise $U(-0.1,0.1)$
was added in the codomain (see Figure 1 (C)), and then we applied the MLOP
algorithm, which reconstructed the manifold (Figure 1 (B)), as well as
denoised function values (Figure 1 (D)). Table 1 summarizes the errors that
correspond to different scenarios. We first notice that due to the denoising
effect on $Q$ points the maximum relative error decreased from $0.31$ to
$0.13$ for noisy data $Q^{(0)}$ as opposed to the clean data $Q^{(150)}$.
Next, we can also see the benefits of using the MLOP algorithm prior to
approximating the function with the RBF method on the new data. In the present
example the maximum relative error of the best RBF execution decreased
dramatically from $0.66$ when running the RBF with centers at $Q^{(0)}$, to an
$0.12$ when running on the quasi-uniform data, and with very low error
variance. A quick comparison between the RBF method and the approximation via
locally weighted average shows that the latter loses the battle to RBF (even
though it produce better results on the clean data versus the noisy one).

(A)

(B)

(C)

(D)

Figure 1: Manifold of orthogonal matrices embedded in a 60-dimensional space.
Shown are the first two coordinates of the point-set (after multiplication by
$A^{-1}$). (A) Scattered data with uniformly distributed noise $U(-0.1;0.1)$
(green), and the initial point-set $Q^{(0)}$ (red). (B) The resulting point-
set of the MLOP algorithm after 150 iterations, $Q^{(150)}$ (red), overlaying
the noisy samples (green). (C) The initial function values evaluated at the
original point-set $Q^{(0)}$ with noise $U(-0.1,0.1)$. The black line shows
the noise-free reference data. (D) Smooth function approximation via the MLOP
algorithm at the data points $Q^{(150)}$.

We then applied the approximation procedure to the non-smooth function given
by $f(x)=\frac{1}{6}(1+\arccos(\cos(10\,\theta)))$. We evaluated the function
at the $P$-points and added the uniform noise $U(-0.1,0.1)$ (see Figure 2
(A)). Then, we applied the MLOP algorithm, which resulted in a reconstructed
manifold, as well as denoised function values; see Figure 2 (B). As this
figure shows, the non-smooth function $f$ is approximated reliably. This is
also reflected in the errors listed in Table 1, which shows that the
approximation error decreased from $0.2$ to $0.1$ after the denoising
procedure. The advantages of the MLOP approach are also demonstrated by the
error decrease, for the best choice of radial basis function, from $0.77$ to
$0.15$ for approximation on a new point-set. This shows the robustness of the
approximation process with respect to the clean data. Here again, we see that
the RBF method produces a better approximation then the weighted average.

(A)

(B)

Figure 2: Manifold of orthogonal matrices embedded in a 60-dimensional space. Shown are the first two coordinates of the point-set (after multiplication by $A^{-1}$). Scattered data with uniformly distributed noise $U(-0.1;0.1)$ (green), and the $Q$ point-set (red). Left: The initial function values evaluated at the original $Q^{(0)}$-points with noise $U(-0.1,0.1)$. The black line shows the noise-free reference data. Right: Approximation of our non-smooth function via MLOP at the data points $Q^{(150)}$. Table 1: Summary of the maximum and root mean squared errors with standard deviations errors of approximation of functions on the $O(2)$ manifold embedded into 60-dimensional space | $f(x)=\frac{1}{4}(1+\sin 10x)$ | $f(x)=\frac{1}{6}(1+\arccos(\cos 10x))$
---|---|---
|  Max relative error | RMSE $\pm$ var |  Max relative error |  RMSE $\pm$ var
Error over $Q^{(k)}$ |  |  |  |
$f(Q^{(0)})$ |  0.31 |  $0.06\pm 0.0011$ |  0.2 |  $0.04\pm 0.0007$
$f(Q^{(150)})$ |  0.13 |  $0.03\pm 0.0003$ |  0.1 |  $0.03\pm 0.0002$
Error over 100 new points
RBF, $\phi_{1}$, centers at $Q^{(0)}$, noisy $f$ |  0.66 |  $0.14\pm 0.0079$ |  0.77 |  $0.16\pm 0.010$
RBF, $\phi_{1}$, centers at $Q^{(150)}$, cleaned $\widetilde{f}$ |  0.12 |  $0.03\pm 0.0002$ |  0.15 |  $0.03\pm 0.0004$
RBF, $\phi_{2}$, centers at $Q^{(150)}$, cleaned $\widetilde{f}$ |  0.12 |  $0.03\pm 0.0002$ |  0.16 |  $0.03\pm 0.0003$
RBF, $\phi_{3}$, centers at $Q^{(150)}$, cleaned $\widetilde{f}$ |  0.24 |  $0.05\pm 0.0009$ |  0.21 |  $0.05\pm 0.0007$
Weighted average |  0.28 |  $0.10\pm 0.0018$ |  0.32 |  $0.09\pm 0.0024$

In what follows we demonstrate our function approximation methodology on
several examples of a low-dimensional manifold embedded in high-dimensional
space. Specifically, we embedded a two-dimensional cylindrical structure and
then a six-dimensional cylindrical structure in $\mathbb{R}^{60}$. We start
with the two-dimensional cylindrical structure. We sampled the structure using
the parameterization

$\displaystyle
p=tv_{1}+\frac{R}{\sqrt{2}}\big{(}\cos(u)v_{2}+\sin(u)v_{3}\big{)}\,,$

where $v_{1}=[1,1,1,1,1,\dotsc,1]$,
$v_{2}=[0,1,-1,0,0,\dotsc,0],v_{3}=[1,0,0,-1,0,\dotsc,0]$
$(v_{1},v_{2},v_{3}\in\mathbb{R}^{60})$, $t\in[0,2]$ and
$u\in[0.1\pi,1.5\pi]$. Using this representation, $800$ equally distributed
(in parameter space) points were sampled with uniformly distributed noise
(i.e., $U(-0.1,0.1)$). We evaluated the function
$f(t,u)=1.3(1+\textup{sin}(0.5u+1.5t))$ at these points, and constructed the
initial $Q$-set by randomly sampling $150$ points (see Figure 3 (A) for the
$P$ and $Q$ data, and Figure 3 (C) for the values of the function at the
$Q$-points). The representative distances of the $P$-set and $Q$-set were
$h_{1}=0.19$ and $h_{2}=0.27$, resepctively. Next, we applied the MLOP
algorithm on the new data of $(P,f(P))$, and extracted the sets $Q$ and
$\widetilde{f}(Q)$ (see Figure 3 (B) for the cleaned $Q$-points, and Figure 3
(D) for the cleaned values at the $Q$-points). In addition, we approximated
the values of $f$ at 100 new points, randomly selected from the reference
data. The evaluation error results are summarized in Table 2. The maximum
relative $L_{1}$ error and the RMSE accompanied with the variance are
summarized in Table 2. It should be noted that since we compare the maximum
error to clean data, the relative error can exceed 1. One can see that the new
data RBF with $\phi_{2}$ as well as $\phi_{3}$ achieve the lowest errors.

(A)

(B)

(C)

(D)

Figure 3: Two-dimensional cylindrical structure embedded in a 60-dimensional
space. The first three coordinates of the point-set are shown. (A) Scattered
data with uniformly distributed noise $U(-0.1;0.1)$ (green), and the initial
point-set $Q^{(0)}$ (red). (B) The resulting point-set of the MLOP algorithm
after 200 iterations, $Q^{(300)}$ (red), overlaying the noisy samples (green).
(C) The initial values of the function at the original $Q^{(0)}$-points with
noise $U(-0.1,0.1)$. (D) MLOP approximation at the data points $Q^{(300)}$.

### Six-dimensional cylindrical structure

Next, we tested our method on higher-dimensional manifolds by utilizing an
$n$-sphere to generate an $(n+1)$-dimensional cylinder (in the example of the
two-dimensional cylinder, we used a circle to generate the structure). Here,
we utilized a five-dimensional sphere to build a six-dimensional manifold,
using the parameterization

$\displaystyle x_{1}=R\cos(u_{1})\,,\quad
x_{2}=R\sin(u_{1})\cos(u_{2}),\quad\ldots,\quad
x_{6}=R\sin(u_{1})\sin(u_{2})\cdots\sin(u_{5})\sin(u_{6})\,.$

We then embedded the sampled data in a 60-dimensional space by the
parametrization

$\displaystyle
p=tv_{0}+R^{2}[x_{1},x_{2},x_{3},x_{4},x_{5},x_{6},0,\dotsc,0]\,,$ (18)

where $R=1.5$, $t\in[0,2]$, $u_{i}\in[0.1\pi,0.6\pi]$, and
$v_{0}\in\mathbb{R}^{60}$ is a vector with 1’s in positions $1,...,d+1$ and 0
in the remaining positions. We randomly sampled the $P$-points from the six-
dimensional cylindrical structure and embedded the sampled data in a
60-dimensional space. This process resulted in 1200 points in the $P$-set. We
evaluated the function $f(u_{1},\ldots,u_{6})=\sum\limits_{i=1}^{6}{u_{i}}$ at
these points, and constructed the initial $Q$-set by randomly sampling $460$
points. Next, uniformly distributed noise $U(-0.2;0.2)$ was added to the
points. To avoid trying to visualize a six-dimensional manifold, we plot here
the cross-section of the cylindrical structure in three dimensions. In Figure
4 (A) we present the $P$\- and $Q$-data, and in Figure 4 (C) the values of the
function at the $Q$-points. The representative distances of the $P$-set and
$Q$-set were $h_{1}=0.24$ and $h_{2}=0.37$, respectively. We then applied the
MLOP algorithm on the new data of $(P,f(P))$, and extracted the sets $Q$ and
$\widetilde{f}(Q)$ (see Figure 4 (B) for the cleaned $Q$-points, and Figure 4
(D) for the cleaned function values at the $Q$ points). In addition, we
approximated the values of the function at 100 new points, randomly selected
from the reference data. The maximum relative $L_{1}$ error and the RMSE
accompanied with the variance are summarized in Table 2. One can see that the
new data RBF with $\phi_{2}$ and the weighted average achieve lower errors.

(A)

(B)

(C)

(D)

Figure 4: Six-dimensional cylindrical structure embedded in a 60-dimensional space. Plot of the cross-section of the cylindrical structure in three dimensions. (A) Scattered data with uniformly distributed noise $U(-0.2;0.2)$ (green), and the initial point-set $Q^{(0)}$ (red). (B) The resulting point-set of the MLOP algorithm after 300 iterations, $Q^{(300)}$ (red) overlaying the noisy samples (green). (C) The initial function values evaluated at the original $Q^{(0)}$ points with noise $U(-0.2,0.2)$. (D) MLOP approximation at the data points $Q^{(300)}$. Table 2: Summary of the maximum and root mean squared error with standard deviations errors of approximation of a function on a two-dimensional and six-dimensional cylindrical manifold embedded in a 60-dimensional space | 2D cylinder in $\mathbb{R}^{60}$ | 6D cylinder in $\mathbb{R}^{60}$
---|---|---
|  Max relative error | RMSE $\pm$ var |  Max relative error |  RMSE $\pm$ var
Error over $Q^{(k)}$ |  |  |  |
$f(Q^{(0)})$ |  0.11 | $0.09\pm 0.0029$ |  0.066 |  $0.08\pm 0.0018$
$f(Q^{(300)})$ |  0.06 | $0.05\pm 0.0012$ |  0.054 |  $0.06\pm 0.0012$
Error over 100 new points
RBF with $\phi_{1}$, centers at $Q^{(0)}$, noisy $f$ |  1.37 | $0.3\pm 0.04$ |  0.42 |  $0.4\pm 0.063$
RBF with $\phi_{1}$, centers at $Q^{(300)}$, cleaned $\widetilde{f}$ |  0.38 | $0.13\pm 0.007$ |  0.26 |  $0.26\pm 0.026$
RBF with $\phi_{2}$, centers at $Q^{(300)}$, cleaned $\widetilde{f}$ |  0.11 | $0.05\pm 0.0009$ |  0.08 |  $0.07\pm 0.0018$
RBF with $\phi_{3}$, centers at $Q^{(300)}$, cleaned $\widetilde{f}$ |  0.10 | $0.04\pm 0.0006$ |  0.13 |  $0.04\pm 0.0007$
Weighted average |  0.2 | $0.1\pm 0.0024$ |  0.09 |  $0.06\pm 0.0027$

### Robustness to Noise

In the following example, we examine the effect of the noise level in the
target domain on the quality of the approximation. To do this numerically, we
sampled a function over a Swiss Roll using the parameterization

$\displaystyle p=\frac{1}{10}[x,y,z,0,\dotsc,0]\,,$

where $x=t\sin(t)$, $y$ is a random number in the range $[-6,6]$, and
$z=t\cos(t)$, with $t=8k/n+2$ and $k\in\mathbb{N}$. The approximated function
was $f(p)=t$. We created a Swiss Roll with 800 data points, and randomly
sampled 200 points as the initial $Q$-set. We added noise with various
magnitudes (0.1, 0.2, 0.5, and 0.7) to the $P$-points as well as to the values
of $f$ at the $P$-points. For example, Figure 5 left shows a case of
approximation with uniformly distributed noise $U(-0.2,0.2)$, while the right
plot presents the denoised version. We see that the data were cleaned both in
the domain and in the codomain of the function. In Figure 6, we plot the error
values under various noise scenarios in the codomain, both for the noisy
$Q^{(0)}$ data and the error of the RBF approximation on the $Q^{(300)}$ data.
One can see that although the approximation error increases on the noisy data
(from 0.03 to 0.22), the approximation error on clean and quasi-uniformly
distributed data error grows only moderately (from 0.02 to 0.18). We also see
that at high levels of noise (e.g., 0.7) the accuracy is good. This shows the
strengths of our approach, and justifies the need for data denoising as well
as uniform sampling before approximation algorithms are applied.

(A)

(B)

Figure 5: Swiss Roll embedded in $\mathbb{R}^{60}$. The figure depicts the
first tree coordinates. Left: The initial values of the function at the
original $Q^{(0)}$ points with noise $U(-0.2,0.2)$, with values indicated by
the color. Right: MLOP function approximation at the data points $Q^{(300)}$
cleaned via the MLOP. Figure 6: Effect of noise level on the accuracy of
function approximation for a Swiss Roll embedded in a 60-dimensional space.
The RMSE error evaluated on the original noisy data is shown in orange, while
the RMSE error on the cleaned data is presented in blue.

##  6 Discussion and Conclusions

In this paper, we consider the problem of approximating a function on a
manifold in high dimensions, with noise present both in the domain and in the
codomain of the function. Given a set of points and the values of an unknown
function evaluated at these points, the goal is to find the approximation of
the function at a new dataset. While in low dimensions this problem did
receive a lot of attention in approximation theory, in high dimensions the
solution is challenged by the curse of dimensionality. In our solution, we
propose to combine the best of both worlds, the Manifold Locally Optimal
Projection (MLOP) [14] and the method of Radial Basis Functions (RBF) methods
[13]. The approximation problem has two steps. First, find a noise-free
representation of the manifold in terms of new points, and the noise-free
values of the function at these points. Then, estimate the value of the
function at a new given point $x$ using the RBF method, with the centers set
at the cleaned points. The MLOP method is used here for noise removal and for
the generation of a quasi-uniform manifold sampling, as a pre-processing step
for the RBF mehod; this improves the approximation results dramatically. In
the paper, we showed that the order of approximation at the new data points is
less than $C_{1}h^{2}+C_{2}h_{2}^{k}$, where $h$ is the representative
distance of the graph of the function $f$ of the initial data, and $h_{2}$ is
the representative distance of the MLOP reconstruction.

A possible future direction would be to investigate a classification problem.
Thus, given data that lies on a manifold, and the corresponding values of a
non-smooth function which received only $k$ values, each representing a
different class. In addition, both the data and the labeling contain noise and
outliers. The research question to be addressed is whether new data can be
classified with high accuracy.

An additional possible future direction would be to deal with optimization of
functions on a manifold. In the past decade, optimization gained a lot of
attention, especially with the rise of Neural Networks computing (NN).
Optimization algorithms are the pillar stones of the NNs, as they are in
charge of constructing the networks, by learning from the training examples.
In our research, we propose introducing a new optimization process, that will
take into account the topology of the data. We would like to utilize the
manifold structure of the data to improve the optimization process, by
incorporating the manifold’s information into the NN optimization of a
function. We propose extending the MLOP framework to deal with this task, by
modifying the definition of the MLOP cost function to include the optimized
function. By extending the definition of the MLOP algorithm, the gradient
descent iterations will find not only the optimal manifold reconstruction, but
also minimize the function.

## Acknowledgments

We would like to thank Dr. Barak Sober for valuable discussions and comments.
This study was supported by a generous donation from Mr. Jacques Chahine, made
through the French Friends of Tel Aviv University, and was partially supported
by ISF grant 2062/18.

## References

  * [1] Amir, A., Levin, D.: High order approximation to non-smooth multivariate functions.  Computer Aided Geometric Design 63, 31–65 (2018)
  * [2] Andras, P.: High-dimensional function approximation with neural networks for large volumes of data.  IEEE Transactions on Neural Networks and Learning Systems 29(2), 500–508 (2017)
  * [3] Barzilai, J., Borwein, J.M.: Two-point step size gradient methods.  IMA Journal of Numerical Analysis 8(1), 141–148 (1988)
  * [4] Bickel, P.J., Li, B., et al.: Local polynomial regression on unknown manifolds.  In: Complex Datasets and Inverse Problems, pp. 177–186. Institute of Mathematical Statistics (2007)
  * [5] Binev, P., Dahmen, W., Lamby, P.: Fast high-dimensional approximation with sparse occupancy trees.  Journal of Computational and Applied Mathematics 235(8), 2063–2076 (2011)
  * [6] Buhmann, M.D.: Radial Basis Functions: Theory and Implementations, vol. 12.  Cambridge University Press, Cambridge Monographs on Applied and Computational Science (2003)
  * [7] Buhmann, M.D., De Marchi, S., Perracchione, E.: Analysis of a new class of rational RBF expansions.  IMA Journal of Numerical Analysis 40(3), 1972–1993 (2020)
  * [8] Chen, M., Jiang, H., Liao, W., Zhao, T.: Efficient approximation of deep ReLU networks for functions on low dimensional manifolds.  In: Advances in Neural Information Processing Systems, pp. 8174–8184 (2019)
  * [9] Coifman, R.R., Lafon, S., Lee, A.B., Maggioni, M., Nadler, B., Warner, F., Zucker, S.W.: Geometric diffusions as a tool for harmonic analysis and structure definition of data: Diffusion maps.  Proceedings of the National Academy of Sciences 102(21), 7426–7431 (2005)
  * [10] Coifman, R.R., Maggioni, M.: Diffusion wavelets.  Applied and Computational Harmonic Analysis 21(1), 53–94 (2006)
  * [11] Cox, T.F., Cox, M.A.: Multidimensional Scaling.  Chapman and Hall, London (2000)
  * [12] DeVore, R., Petrova, G., Wojtaszczyk, P.: Approximation of functions of few variables in high dimensions.  Constructive Approximation 33(1), 125–143 (2011)
  * [13] Dyn, N., Levin, D.: Iterative solution of systems originating from integral equations and surface interpolation.  SIAM Journal on Numerical Analysis 20(2), 377–390 (1983)
  * [14] Faigenbaum-Golovin, S., Levin, D.: Manifold reconstruction and denoising from scattered data in high dimension via a generalization of $l_{1}$-median (2020)
  * [15] Fisher, R.A.: The use of multiple measurements in taxonomic problems.  Annals of Eugenics 7(2), 179–188 (1936)
  * [16] Grohs, P., Sprecher, M.: Projection-based quasiinterpolation in manifolds.  SAM Report 23 (2013)
  * [17] He, X., Niyogi, P.: Locality preserving projections.  In: Advances in Neural Information Processing Systems, pp. 153–160 (2004)
  * [18] Levin, D.: The approximation power of moving least-squares.  Mathematics of Computation 67(224), 1517–1531 (1998)
  * [19] Levin, D.: Between moving least-squares and moving least-$\ell_{1}$.  BIT Numerical Mathematics 55(3), 781–796 (2015)
  * [20] Li, L.: A new complexity bound for the least-squares problem.  Computers & Mathematics with Applications 31(12), 15–16 (1996)
  * [21] Lin, T., Zha, H.: Riemannian manifold learning.  IEEE Transactions on Pattern Analysis and Machine Intelligence 30(5), 796–809 (2008)
  * [22] Lipman, Y., Cohen-Or, D., Levin, D., Tal-Ezer, H.: Parameterization-free projection for geometry reconstruction.  In: ACM Transactions on Graphics (TOG), vol. 26, p. 22. ACM (2007)
  * [23] Pearson, K.: LIII. On lines and planes of closest fit to systems of points in space.  The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science 2(11), 559–572 (1901)
  * [24] Roweis, S.T., Saul, L.K.: Nonlinear dimensionality reduction by locally linear embedding.  Science 290(5500), 2323–2326 (2000)
  * [25] Shaham, U., Cloninger, A., Coifman, R.R.: Provable approximation properties for deep neural networks.  Applied and Computational Harmonic Analysis 44(3), 537–557 (2018)
  * [26] Sober, B., Aizenbud, Y., Levin, D.: Approximation of functions over manifolds: A moving least-squares approach.  arXiv preprint arXiv:1711.00765 (2017)
  * [27] Tenenbaum, J.B., De Silva, V., Langford, J.C.: A global geometric framework for nonlinear dimensionality reduction.  Science 290(5500), 2319–2323 (2000)
  * [28] Vardi, Y., Zhang, C.H.: The multivariate L1-median and associated data depth.  Proceedings of the National Academy of Sciences 97(4), 1423–1426 (2000)
  * [29] Woodruff, D.P., et al.: Sketching as a tool for numerical linear algebra.  Foundations and Trends in Theoretical Computer Science 10(1–2), 1–157 (2014)
  * [30] Wu, Zong-min., Schaback, R.: Local error estimates for radial basis function interpolation of scattered data.  IMA journal of Numerical Analysis 13(1), 13–27 (1993)
  * [31] Yoon, J.: Spectral approximation orders of radial basis function interpolation on the Sobolev space.  SIAM Journal on Mathematical Analysis 33(4), 946–958 (2001)


# On the Expressibility of
the Reconstructional Color Refinement

V. Arvind, Johannes Köbler, Oleg Verbitsky† The Institute of Mathematical
Sciences (HBNI) and Chennai Mathematical Institute, Chennai, India.Institut
für Informatik, Humboldt-Universität zu Berlin, Germany.Supported by DFG grant
KO 1053/8–2. On leave from the IAPMM, Lviv, Ukraine.

###### Abstract

One of the most basic facts related to the famous Ulam reconstruction
conjecture is that the connectedness of a graph can be determined by the deck
of its vertex-deleted subgraphs, which are considered up to isomorphism. We
strengthen this result by proving that connectedness can still be determined
when the subgraphs in the deck are given up to equivalence under the color
refinement isomorphism test. Consequently, this implies that connectedness is
recognizable by Reconstruction Graph Neural Networks, a recently introduced
GNN architecture inspired by the reconstruction conjecture (Cotta, Morris,
Ribeiro 2021).

##  1 Introduction

For a vertex $v$ of a graph $G$, the _vertex-deleted subgraph_ $G\setminus v$
is obtained by removing $v$ along with all incident edges from $G$. The
isomorphism type of $G\setminus v$ is sometimes referred to as a _card_ , and
then the multiset of all cards is called the _deck_ of $G$. The famous Ulam
reconstruction conjecture [19] says that every graph $G$ with more than 2
vertices is, up to isomorphism, reconstructible from its deck. Though
published in 1960, the problem collection [19] is a descendant of the much
earlier famous _Scottish Book_ [12]. According to Harary [9], the conjecture
has already been discussed since 1929 and, about a hundred years later, it
remains a notoriously hard open problem in graph theory.

The subject has a rich literature, and one of the most frequently taken
approaches is examination of how much information about a graph $G$ can be
retrieved from its deck. One of the earliest results in this direction is
referred to in [9] as a “little theorem”: $G$ is connected if and only if at
least two subgraphs in its deck are connected. This implies the following fact
(see also [15, Theorem 2.2])

######  Theorem 1 (Harary [8]).

The connectedness of a graph $G$ is determined by the multiset of all vertex-
deleted subgraphs $G\setminus v$ of $G$.

Recall that the subgraphs $G\setminus v$ are considered up to isomorphism,
which means that the full structural information about them is available. It
is natural to ask which partial information suffices to decide the
connectedness of $G$. For example, the aforementioned “little theorem” of
Harary shows that it would be just enough to know the connectedness of each
card in the deck. The question addressed in this note is whether connectedness
can still be decided when the subgraphs $G\setminus v$ are considered up to an
important relaxation of the isomorphism relation, namely up to equivalence
under the color refinement test (also known as one-dimensional Weisfeiler-
Leman equivalence).

_Color refinement_ iteratively computes an isomorphism-invariant coloring of
the vertex set of an input graph $G$. Let $\operatorname{\mathsf{CR}}(G)$
denote the mutisets of colors assigned by the algorithm to the vertices of $G$
(see Section 2 for definitions). Color refinement _distinguishes_ two graphs
$G$ and $H$ if
$\operatorname{\mathsf{CR}}(G)\neq\operatorname{\mathsf{CR}}(H)$. If
$\operatorname{\mathsf{CR}}(G)=\operatorname{\mathsf{CR}}(H)$, the graphs are
called _CR-equivalent_. Since $\operatorname{\mathsf{CR}}(G)$ is a graph
invariant, CR-equivalence is a coarser equivalence relation than isomorphism
of graphs.

CR-equivalence has natural characterizations discovered in different research
contexts. Two graphs are CR-equivalent exactly when they are equivalent in the
two-variable first-order logic with counting quantifiers [10], when they have
equally many vertices and share a common covering [1], when they are
fractionally isomorphic [17], and when the number of homomorphisms from every
tree to these graphs is the same [5].

Color refinement is a versatile tool with applications in diverse fields [1,
7, 16]. Most prominently, it is a practical method of isomorphism testing. The
graph invariant $\operatorname{\mathsf{CR}}(G)$ can be computed efficiently
[10] and identifies the isomorphism type of $G$ for almost all graphs [2]. In
machine learning, color refinement is used for analysis and comparison of
graph-structured data [18] and for understanding the expressiveness of graph
neural networks (GNNs) [13].

Our main result is the following strengthening of Theorem 1.

######  Theorem 2.

The connectedness of a graph $G$ is determined by the multiset of the
invariants $\operatorname{\mathsf{CR}}(G\setminus v)$ for all vertex-deleted
subgraphs of $G$.

Note that the invariant $\operatorname{\mathsf{CR}}(G)$ is not powerful enough
to decide whether or not $G$ is connected. Indeed,
$\operatorname{\mathsf{CR}}(C_{6})=\operatorname{\mathsf{CR}}(2C_{3})$, where
$C_{n}$ denotes the cycle graph of length $n$, and $2C_{3}$ stands for the
vertex-disjoint union of two copies of $C_{3}$. This cuts off the easiest
route to Theorem 2 through the potential use of the “little theorem” of Harary
and shows that running color refinement on the deck of a graph can be more
beneficial than running it on the graph alone. Theorem 2 is interesting in
several respects, which we discuss now.

### Consequences for GNNs.

The expressive power of a representative class of GNNs is characterizable in
terms of color refinement [14]. In the search for a more expressive GNN
architecture, Cotta, Morris, Ribeiro [4] introduce the model of _$k$
-Reconstruction GNNs_, whose conceptual novelty consists in enhancing the
standard GNNs by running them on the deck of $k$-vertex subgraphs of a graph
$G$ rather than on $G$ itself. The simplest instructive example, highlighted
in [4], is given by considering the 6-cycle graph $C_{6}$. As we just
mentioned, $C_{6}$ is not identified by $\operatorname{\mathsf{CR}}(C_{6})$.
However, $C_{6}$ is reconstructible from its deck, which consists of six
copies of the 5-path graph $P_{5}$. Since $P_{5}$ is identified by
$\operatorname{\mathsf{CR}}(P_{5})$, the cycle graph $C_{6}$ is identified by
color refinement applied to its deck. We state an immediate consequence of
Theorem 2, providing further evidence that $k$-Reconstruction GNNs can often
be more powerful.

######  Corollary 3.

The connecteness of an $n$-vertex graph is recognizable by
$(n-1)$-Reconstruction GNNs.

### Graph reconstruction.

Theorem 2 shows that a prototypical reconstructible graph property can be
recognized even from an incomplete information about the deck. The research on
the reconstruction conjecture provides many results of this kind, and we here
mention two of them which concern connectedness. One can determine the
connectedness of an $n$-vertex graph from any $\lfloor n/2\rfloor+2$ of its
cards [3]. This can also be done using the $k$-deck of a graph whenever $k\geq
9n/10$ [6]. Note in this respect that the $k$-deck of a graph determines its
$k^{\prime}$-deck whenever $k^{\prime}\leq k$ and, therefore, the $k$-deck
with a smaller $k$ potentially provides less information. Turning back to
Theorem 2, it is natural to combine it with the logical characterization of
color refinement [10]. As a consequence, the connectedness of a graph can be
determined based on the information about its deck expressible in the two-
variable logic with counting quantifiers.

### Isomorphism invariants of graphs.

Let $\operatorname{\mathsf{D-CR}}(G)$ denote the multiset of
$\operatorname{\mathsf{CR}}(G\setminus v)$ over all vertices $v$ of $G$. This
graph invariant is interesting to compare not only with with
$\operatorname{\mathsf{CR}}(G)$ but as well with other related invariants.
Color refinement is also commonly known as the one-dimensional Weisfeiler-
Leman algorithm. Weisfeiler and Leman [20] actually introduced the two-
dimensional version, exploiting the same idea as color refinement but coloring
pairs of vertices instead of single vertices. Let
$\operatorname{\mathsf{WL}}(G)$ denote the graph invariant computed by this
algorithm on input $G$. A higher dimension results in a more powerful
isomorphism test and, in particular, $\operatorname{\mathsf{CR}}(G)$ is
determined by $\operatorname{\mathsf{WL}}(G)$. It is not hard to show that
$\operatorname{\mathsf{D-CR}}(G)$ is also determined by
$\operatorname{\mathsf{WL}}(G)$. Whether $\operatorname{\mathsf{CR}}(G)$,
$\operatorname{\mathsf{D-CR}}(G)$, and $\operatorname{\mathsf{WL}}(G)$ form a
hierarchy of graph invariants linearly ordered by their strength is an open
question. More precisely, it is unknown whether or not
$\operatorname{\mathsf{CR}}(G)$ is determined by
$\operatorname{\mathsf{D-CR}}(G)$. This question is discussed in [17, Section
6.8] and [4, Section 4.1]. As we mentioned above,
$\operatorname{\mathsf{CR}}(G)$ does not allow us to see whether $G$ is
connected. On the other hand, it is well known that the connectedness of a
graph $G$ can be determined from $\operatorname{\mathsf{WL}}(G)$. Theorem 2
provides a formal strengthening of the last fact by showing that the
connectedness of $G$ is determined even by $\operatorname{\mathsf{D-CR}}(G)$.

##  2 Color refinenemt

The vertex set of a graph $G$ is denoted by $V(G)$. The _neighborhood_
$N_{G}(x)$ of a vertex $x\in V(G)$ consists of all vertices adjacent to $x$.
We write $\deg_{G}x=|N_{G}(x)|$ to denote the degree of $x$.

On an input graph $G$, the _color-refinement_ algorithm (to be abbreviated as
_CR_) iteratively computes a sequence of colorings $C^{r}_{G}$ of $V(G)$. The
initial coloring $C^{0}_{G}$ is monochromatic, i.e.,
$C^{0}_{G}(x)=C^{0}_{G}(x^{\prime})$ for all $x,x^{\prime}\in V(G)$, where the
same initial color is used for all $G$. Each subsequent color of a vertex $x$
is obtained from the preceding coloring by the rule

$C^{r+1}_{G}(x)=\left\\{\\!\\!\left\\{C^{r}_{G}(y)\right\\}\\!\\!\right\\}_{y\in
N_{G}(x)},$

where $\left\\{\\!\\!\left\\{\right\\}\\!\\!\right\\}$ denotes a multiset.
Note that $C^{1}_{G}(x)=C^{1}_{G}(x^{\prime})$ exactly when
$\deg_{G}(x)=\deg_{G}(x^{\prime})$, that is, the color $C^{1}_{G}(x)$ can just
be seen as the degree of the vertex $x$.

We define $C_{G}(x)=(C^{0}_{G}(x),C^{1}_{G}(x),C^{2}_{G}(x),\ldots)$ as the
sequence of $C^{r}_{G}(x)$ for all $r$ and set
$\operatorname{\mathsf{CR}}(G)=\left\\{\\!\\!\left\\{C_{G}(x)\right\\}\\!\\!\right\\}_{x\in
V(G)}$. We say that CR _distinguishes_ two graphs $G$ and $H$ if
$\operatorname{\mathsf{CR}}(G)\neq\operatorname{\mathsf{CR}}(H)$. If $G$ and
$H$ are indistinguishable by CR, i.e.,
$\operatorname{\mathsf{CR}}(G)=\operatorname{\mathsf{CR}}(H)$, they are
referred to as _CR-equivalent_.

Defining $C_{G}(x)$ as an infinite sequence facilitates proving Theorem 2 in
the subsequent sections. For example, this definition is needed for the
equivalences in Lemma 8 below. Color refinement is, nevertheless, an
absolutely practical procedure. Indeed, if $n$-vertex graphs $G$ and $H$ are
distinguishable by CR, then they are distinguishable after the $n$-th
refinement round at latest in the sense that
$\left\\{\\!\\!\left\\{C^{n}_{G}(x)\right\\}\\!\\!\right\\}_{x\in
V(G)}\neq\left\\{\\!\\!\left\\{C^{n}_{H}(x)\right\\}\\!\\!\right\\}_{x\in
V(H)}$. Note in this respect that the color names of $C^{r}_{G}(x)$ increase
exponentially with $r$ if they are encoded in a straightforward way. This can
be avoided by renaming the colors synchronously on $G$ and $H$ after each
round.

##  3 Useful lemmas

As $C^{1}_{G}(x)=\deg_{G}x$, we have
$C^{2}_{G}(x)=\\{\\!\\!\\{\deg_{G}(y)\,:\,\begin{array}[]{@{}l@{}}y\in
N_{G}(x)\end{array}\\}\\!\\!\\}$, which can be referred to as the _iterated
degree_ of $x$. We begin with a particularly useful fact about the
reconstructibility of iterated degrees.

######  Lemma 4 (Nash-Williams [15, Lemma 3.3]).

Suppose that graphs $G$ and $H$ have the same deck up to CR-equivalence, that
is,

$\left\\{\\!\\!\left\\{\operatorname{\mathsf{CR}}(G\setminus
x)\right\\}\\!\\!\right\\}_{x\in
V(G)}=\left\\{\\!\\!\left\\{\operatorname{\mathsf{CR}}(H\setminus
y)\right\\}\\!\\!\right\\}_{y\in V(H)}.$

If $\operatorname{\mathsf{CR}}(G\setminus
u)=\operatorname{\mathsf{CR}}(H\setminus v)$ for a vertex $u\in V(G)$ and a
vertex $v\in V(H)$, then $C^{2}_{G}(u)=C^{2}_{H}(v)$.

If two graphs have different number of vertices, then they are
straightforwardly distinguished by CR. Nevertheless, such a pair of graphs
might still be similar in the sense that the CR colors occurring in them are
the same, and the difference is only in the color multiplicities. In the proof
of Theorem 2, this similarity concept will play an important role.

######  Definition 5.

We say that graphs $A$ and $B$ are _CR-similar_ and write $A\equiv B$ if
$\left\\{C_{A}(x)\right\\}_{x\in V(A)}=\left\\{C_{B}(x)\right\\}_{x\in V(B)}$.
That is, we require that the _sets_ of the colors produced by CR on $A$ and
$B$ be equal (although the _multisets_ of these colors might be unequal).

Note that two graphs with different numbers of vertices can be CR-similar, for
example, $C_{3}\equiv C_{4}$.

######  Lemma 6.

  1. 1.

If $A^{\prime}$ is a proper subgraph of a connected graph $A$, then
$A\not\equiv A^{\prime}$.

  2. 2.

Let $A^{\prime}$ be a proper subgraph of $A$ and $B^{\prime}$ be a proper
subgraph of $B$. Suppose that both $A$ and $B$ are connected. Then it is
impossible that simultaneously $A\equiv B^{\prime}$ and $B\equiv A^{\prime}$.

The proof is based on a characterization of the CR-color $C_{G}(x)$ in terms
of the universal cover of the graph $G$ rooted at the vertex $x$. We first
give the relevant definitions.

Let $G$ and $K$ be connected graphs, where $G$ is finite graph, while $K$ can
be either finite or infinite. Let $\alpha$ be a surjective homomorphism from
$K$ onto $G$. We call $\alpha$ a _covering map_ if its restriction to the
neighborhood of each vertex in $K$ is bijective. If such an $\alpha$ exists,
we say that $K$ _covers_ $G$. A graph $U$ is called a _universal cover_ of $G$
if $U$ covers every graph covering $G$. A universal cover $U=U_{G}$ of $G$ is
unique up to isomorphism. It can be seen as an unfolding of $G$ into a
(possibly infinite) tree starting from an arbitrarily chosen vertex $x$ of
$G$. Denote this unfolding by $U_{G,x}$ and note that $U_{G,x}\cong U_{G}$
whatever $x$. If $G$ is a tree, then $U_{G,x}\cong G$. If $G$ contains a
cycle, then $U_{G,x}$ is an infinite tree.

The unfolding $U_{G,x}$ can formally be described as follows. A sequence of
vertices $x_{0}x_{1}\ldots x_{k}$ in $G$ is called a _walk starting at
$x_{0}$_ if every two successive vertices $x_{i}$ and $x_{i+1}$ are adjacent.
If $x_{i+1}\neq x_{i-1}$ for all $0<i<k$, then the walk is _non-backtracking_.
Now, $U_{x}$ is the graph whose vertices are all non-backtracking walks in $G$
starting at $x_{0}=x$, with any two walks of the form $x_{0}x_{1}\ldots x_{k}$
and $x_{0}x_{1}\ldots x_{k}x_{k+1}$ being adjacent.

A straightforward inductive argument shows that a covering map $\alpha$
preserves the coloring produced by CR, that is,
$C_{K}^{r}(u)=C_{G}^{r}(\alpha(u))$ for all $r$. This has the following
consequence. Let $x$ be a vertex of a connected graph $G$ and $y$ be a vertex
of a connected graph $H$. If $U_{G,x}$ and $U_{H,y}$ are isomorphic as rooted
trees, then $C_{G}(x)=C_{H}(y)$. The converse implication is also true; see,
e.g., Lemmas 2.3 and 2.4 in [11]. We state this equivalence in a more precise
form. Let $U_{G,x}^{r}$ denote the rooted tree $U_{G,x}$ truncated at level
$r$.

######  Lemma 7.

$C_{G}^{r}(x)=C_{H}^{r}(y)$ if and only if $U_{G,x}^{r}\cong U_{H,y}^{r}$,
where $\cong$ denotes the isomorphism of rooted trees.

The characterization of CR-colors in terms of universal covers readily implies
the following equivalences.

######  Lemma 8.

The following conditions are equivalent for connected graphs $G$ and $H$:

  1. 1.

$G\equiv H$;

  2. 2.

$U_{G}\cong U_{H}$;

  3. 3.

$\left\\{C_{G}(x)\right\\}_{x\in V(G)}\cap\left\\{C_{H}(y)\right\\}_{y\in
V(H)}\neq\emptyset$.

Lemma 8 has a simple consequence about the CR-similarity of not necessarily
connected graphs. For graphs $G_{1},\ldots,G_{s}$, we write
$G_{1}+\cdots+G_{s}$ to denote the vertex-disjoint union of these graphs.

######  Lemma 9.

  1. 1.

If $G\equiv H$, then for every connected component $G^{\prime}$ of $G$ there
is a connected component $H^{\prime}$ of $H$ such that $G^{\prime}\equiv
H^{\prime}$ (and vice versa).

  2. 2.

If $G\equiv H$ where $G$ is connected and $H=H^{\prime}+H^{\prime\prime}$,
then $G\equiv H^{\prime}\equiv H^{\prime\prime}$.

###### Proof of Lemma 6.

Let us extend the notation $U_{G,x}^{r}$ to not necessarily connected graphs
by setting $U_{G,x}^{r}=U_{G^{\prime},x}^{r}$ where $G^{\prime}$ is the
connected component of $G$ containing the vertex $x$.

1. Let the positive integer $D$ be one more than the diameter of $A$. Fix a vertex $u$ of the subgraph $A^{\prime}$ maximizing the number of vertices in the truncated tree $U_{A^{\prime},u}^{D}$. The description of a universal cover in terms of non-backtracking walks readily implies that $U_{A,u}^{D}$ has strictly more vertices than $U_{A^{\prime},u}^{D}$ and, hence, than $U_{A^{\prime},w}^{D}$ for any vertex $w$ of $A^{\prime}$. Indeed, our choice of $D$ ensures that there is a walk of length at most $D$ starting from $u$ and passing through an edge of $A$ absent in $A^{\prime}$. Lemma 7, therefore, implies that $C_{A}(u)\neq C_{A^{\prime}}(w)$ for any $w\in V(A^{\prime})$.

2. Towards a contradiction, assume that $A\equiv B^{\prime}$ and $B\equiv A^{\prime}$. Let $D$ be the positive integer that is one more than the larger of the diameters of $A$ and $B$. Choose a vertex $u$ of the subgraph $A^{\prime}$ maximizing the number of vertices in the truncated tree $U_{A^{\prime},u}^{D}$. As in the proof of the first part, note that $U_{A,u}^{D}$ has strictly more vertices than $U_{A^{\prime},u}^{D}$ and, hence, than $U_{A^{\prime},w}^{D}$ for any vertex $w$ of $A^{\prime}$. The assumption $B\equiv A^{\prime}$ implies by Lemma 7 that every tree $U_{A^{\prime},w}^{D}$ for $w\in V(A^{\prime})$ has an isomorphic mate $U_{B,w^{\prime}}^{D}$ for some $w^{\prime}\in V(B)$ and vice versa. Therefore, $U_{A,u}^{D}$ has strictly more vertices than $U_{B,w}^{D}$ for any vertex $w$ of $B$. The same argument with the roles of $A$ and $B$ interchanged provides us with a vertex $v$ in $B$ which has strictly more vertices than $U_{A,w}^{D}$ for any vertex $w$ of $A$, in particular, than $U_{A,u}^{D}$. This contradiction completes the proof. ∎

##  4 Proof of Theorem 2

Let $G$ be a connected graph and $H$ any disconnected graph. We have to prove
that the multisets
$\left\\{\\!\\!\left\\{\operatorname{\mathsf{CR}}(G\setminus
x)\right\\}\\!\\!\right\\}_{x\in V(G)}$ and
$\left\\{\\!\\!\left\\{\operatorname{\mathsf{CR}}(H\setminus
y)\right\\}\\!\\!\right\\}_{y\in V(H)}$ are unequal. We will assume that these
multisets are equal and derive a contradiction.

Let $H_{1},\ldots,H_{s}$, where $s\geq 2$, be the connected components of $H$.
Let $V=V(H)$ and $V_{i}=V(H_{i})$ for $i\leq s$. Identify $V(G)$ with $V$ so
that

$\operatorname{\mathsf{CR}}(G\setminus
v)=\operatorname{\mathsf{CR}}(H\setminus v)$ (1)

for all $v\in V$. By Lemma 4, this yields

$C^{2}_{G}(v)=C^{2}_{H}(v)$

for all $v\in V$. Taking into account the decomposition of $H$ into connected
components, we have

$C^{2}_{G}(v)=C^{2}_{H_{i}}(v)\text{ if }v\in V_{i}.$ (2)

Denote the set of articulation (or cut) vertices of $G$ by $A$, and set
$B=V\setminus A$. Let $B_{i}=B\cap V_{i}$. If $v\in B_{i}$, then $G\setminus
v$ is connected and $H\setminus v=H_{1}+\ldots+H_{i}\setminus v+\ldots+H_{s}$.
By Part 2 of Lemma 9, Equality (1) implies that

$G\setminus v\equiv H_{i}\setminus v\equiv H_{j}\text{ if }v\in B_{i}\text{
and }j\neq i.$ (3)

Applying the similarity $G\setminus v\equiv H_{j}$ to any two vertices in
$B_{i}$, we derive the similarity $G\setminus v\equiv G\setminus v^{\prime}$
whenever $v,v^{\prime}\in B_{i}$. Applying Lemma 4 in the particular case of
$G=H$, we conclude that

$C^{2}_{G}(v)=C^{2}_{G}(v^{\prime})\text{ for all }v,v^{\prime}\in B_{i}.$ (4)

Let $i\neq j$ and assume that $B_{i}\neq\emptyset$ and $B_{j}\neq\emptyset$.
Let $x\in B_{i}$ and $y\in B_{j}$. By (3), we obtain $H_{i}\setminus x\equiv
H_{j}$ and $H_{j}\setminus y\equiv H_{i}$, which contradicts Part 2 of Lemma
6. Therefore, only one of the sets $B_{i}$ is not empty, that is, $B=B_{i}$
for some $i$. Without loss of generality, we suppose that $B=B_{1}$.

Claim A.  $B_{1}=V_{1}$.

Proof of Claim A.  Towards a contradiction, assume that $V_{1}$ contains a cut
vertex $v$ of $G$. By Part 1 of Lemma 9, Equality (1) implies that the
disconnected graph $G\setminus v$ has a connected component $G^{\prime}$ such
that $G^{\prime}\equiv H_{2}$. Choose a vertex $x\in B$ in a connected
component of $G\setminus v$ different from $G^{\prime}$. Since $x\in B_{1}$,
the similarity (3) implies that $G\setminus x\equiv H_{2}$ and, therefore,
$G\setminus x\equiv G^{\prime}$. However, $G^{\prime}$ is a subgraph of the
connected graph $G\setminus x$, contrary to Part 1 of Lemma 6.
$\,\triangleleft$

By Claim 4 and Equality (4), we have

$C^{2}_{G}(x)=C^{2}_{G}(x^{\prime})\text{ for all }x,x^{\prime}\in V_{1}.$

Taking into account Equality (2), this implies that

$C^{2}_{H_{1}}(x)=C^{2}_{H_{1}}(x^{\prime})\text{ for all }x,x^{\prime}\in
V_{i}.$ (5)

In particular, $H_{1}$ is a regular graph, say, of degree $d$.

Consider the block-cut tree of $G$. Recall that it describes the decomposition
of $G$ into biconnected components, that is, maximal biconnected subgraphs,
which are also called blocks. The vertex set of the tree consists of the
blocks and the cut vertices of $G$. A block $M$ is adjacent to a cut vertex
$c$ if $c\in M$.

Let $L$ be a block that is a leaf in the block-cut tree of $G$. Let $a$ be the
single cut vertex of $G$ belonging to $L$. By Claim 4, $a$ belongs to
$A=V\setminus B=V_{2}\cup\ldots\cup V_{s}$, and the other vertices of $L$
belong to $B=B_{1}=V_{1}$. Without loss of generality, suppose that $a\in
V_{2}$. Let $x$ be a vertex in $L$ adjacent to $a$. Since $x\in B_{1}$, the
similarity (3) implies that

$G\setminus x\equiv H_{1}\setminus x\equiv H_{2}.$ (6)

Since $H_{1}$ is a $d$-regular graph, every vertex of $H_{1}\setminus x$ has
degree $d-1$ or $d$. By the similarity (6), the neighbors of $a$ in $H_{2}$
are of degree $d-1$ or $d$. By Equality (2), the same holds true for $a$ in
$G$. Therefore, every connected component $G^{\prime}$ of $G\setminus a$ has a
vertex of degree $d-2$ or $d-1$ (namely a neighbor of $a$ in $G$). As a
consequence, $G^{\prime}\not\equiv H_{1}$ for every $G^{\prime}$. On the other
hand, Equality (1) implies that

$G\setminus a\equiv H_{1}+H_{2}\setminus a+\cdots+H_{s}.$

By Part 1 of Lemma 9, we conclude that $H_{1}\equiv G^{\prime}$ for some
component $G^{\prime}$ of $G\setminus a$. This contradiction completes the
proof of Theorem 2.

## References

  * [1] D. Angluin.  Local and global properties in networks of processors (extended abstract).  In Proc. of the 12th Annual ACM Symposium on Theory of Computing (STOC’80), pages 82–93. ACM, 1980.
  * [2] L. Babai, P. Erdös, and S. M. Selkow.  Random graph isomorphism.  SIAM J. Comput., 9(3):628–635, 1980.
  * [3] A. Bowler, P. Brown, T. Fenner, and W. Myrvold.  Recognizing connectedness from vertex-deleted subgraphs.  Journal of Graph Theory, 67(4):285–299, 2011.
  * [4] L. Cotta, C. Morris, and B. Ribeiro.  Reconstruction for powerful graph representations.  In M. Ranzato, A. Beygelzimer, Y. Dauphin, P. Liang, and J. W. Vaughan, editors, Advances in Neural Information Processing Systems (NeurIPS 2021), volume 34, pages 1713–1726. Curran Associates, Inc., 2021.
  * [5] Z. Dvořák.  On recognizing graphs by numbers of homomorphisms.  Journal of Graph Theory, 64(4):330–342, 2010.
  * [6] C. Groenland, T. Johnston, A. Scott, and J. Tan.  Reconstruction from smaller cards.  Technical report, arxiv.org/abs/2103.13359, 2023.
  * [7] M. Grohe, K. Kersting, M. Mladenov, and P. Schweitzer.  Color refinement and its applications.  In An Introduction to Lifted Probabilistic Inference. The MIT Press, 2021.
  * [8] F. Harary.  On the reconstruction of a graph from a collection of subgraphs.  In Theory of Graphs and its Applications, pages 47–52. Academic Press: New York, 1964.
  * [9] F. Harary.  A survey of the reconstruction conjecture.  In Graphs and Combinatorics, volume 406 of Lect. Notes Math, pages 18–28. Springer: Berlin, Heidelberg, 1974.
  * [10] N. Immerman and E. Lander.  Describing graphs: A first-order approach to graph canonization.  In A. Selman, editor, Complexity Theory Retrospective, pages 59–81. Springer-Verlag, 1990.
  * [11] A. Krebs and O. Verbitsky.  Universal covers, color refinement, and two-variable counting logic: Lower bounds for the depth.  In 30th Annual ACM/IEEE Symposium on Logic in Computer Science (LICS’15), pages 689–700. IEEE Computer Society, 2015.
  * [12] R. D. Mauldin, editor.  The Scottish Book.  Birkhäuser, Boston, MA, 1981.  Mathematics from the Scottish Café, Including selected papers presented at the Scottish Book Conference held at North Texas State University, Denton, Tex., May 1979.
  * [13] C. Morris, Y. Lipman, H. Maron, B. Rieck, N. M. Kriege, M. Grohe, M. Fey, and K. Borgwardt.  Weisfeiler and Leman go machine learning: The story so far.  Journal of Machine Learning Research, 24(333):1–59, 2023.
  * [14] C. Morris, M. Ritzert, M. Fey, W. L. Hamilton, J. E. Lenssen, G. Rattan, and M. Grohe.  Weisfeiler and Leman go neural: Higher-order graph neural networks.  In The 33-rd Conference on Artificial Intelligence (AAAI’19), pages 4602–4609. AAAI Press, 2019.
  * [15] C. S. J. A. Nash-Williams.  The reconstruction problem.  In Selected topics in graph theory, pages 205–236. Academic Press: New York, 1978.
  * [16] C. Riveros, B. Scheidt, and N. Schweikardt.  Using color refinement to boost enumeration and counting for acyclic CQs of binary schemas.  Technical report, arxiv.org/abs/2405.12358, 2024.
  * [17] E. R. Scheinerman and D. H. Ullman.  Fractional graph theory. A rational approach to the theory of graphs.  Wiley: John Wiley & Sons, 1997.
  * [18] N. Shervashidze, P. Schweitzer, E. J. van Leeuwen, K. Mehlhorn, and K. M. Borgwardt.  Weisfeiler-Lehman graph kernels.  Journal of Machine Learning Research, 12:2539–2561, 2011.
  * [19] S. M. Ulam.  A collection of mathematical problems, volume 8 of Intersci. Tracts Pure Appl. Math. Interscience Publishers, New York, NY, 1960.
  * [20] B. Weisfeiler and A. Leman.  The reduction of a graph to canonical form and the algebra which appears therein.  NTI, Ser. 2, 9:12–16, 1968.  English translation is available at https://www.iti.zcu.cz/wl2018/pdf/wl_paper_translation.pdf.


#  PerSign: Personalized Bangladeshi Sign Letters Synthesis

Mohammad Imrul Jubair  University of Colorado Boulder, United States
<EMAIL_ADDRESS>,  Ali Ahnaf  Ahsanullah University of Science
and Technology, Bangladesh<EMAIL_ADDRESS>,  Tashfiq Nahiyan Khan
Ahsanullah University of Science and Technology, Bangladesh
<EMAIL_ADDRESS>,  Ullash Bhattacharjee  Ahsanullah University of
Science and Technology, Bangladesh<EMAIL_ADDRESS>and  Tanjila Joti
Ahsanullah University of Science and Technology, Bangladesh
<EMAIL_ADDRESS>
(2022)

###### Abstract.

Bangladeshi Sign Language (BdSL)—like other sign languages— is tough to learn
for general people, especially when it comes to expressing letters. In this
poster, we propose PerSign, a system that can reproduce a person’s image by
introducing sign gestures in it. We make this operation “personalized”, which
means the generated image keeps the person’s initial image profile–face, skin
tone, attire, background—unchanged while altering the hand, palm, and finger
positions appropriately. We use an image-to-image translation technique and
build a corresponding unique dataset to accomplish the task. We believe the
translated image can reduce the communication gap between signers111person who
uses sign language. and non-signers without having prior knowledge of BdSL.

Bangladeshi Sign Language (BdSL), Image-to-Image Translation, Sign Letters
Synthesis.

††journalyear: 2022††copyright: rightsretained††conference: The Adjunct
Publication of the 35th Annual ACM Symposium on User Interface Software and
Technology; October 29-November 2, 2022; Bend, OR, USA††booktitle: The Adjunct
Publication of the 35th Annual ACM Symposium on User Interface Software and
Technology (UIST ’22 Adjunct), October 29-November 2, 2022, Bend, OR,
USA††doi: 10.1145/3526114.3558712††isbn: 978-1-4503-9321-8/22/10††ccs: Human-
centered computing Text input

Figure 1. An overview of our PerSign. Here,  shows the novelty of the work.
Top row of  indicates the Bangladeshi sign language (BdSL) classification
where an image with gesture is classified as the Bengali letter . Middle
exhibits the vice versa of the previous one, where symbols of signs are
generated using generative adversarial networks (GANs). In contrast, ours
(last) generate personalized signs—it reproduces a person’s image with the
desired sign for a given letter while keeping the profile—face, dress,
etc.—unchanged.  shows an abstraction of our working pipeline (see section-1).
shows image translation module that takes letter and profile image, and GAN
translates the image by learning from our dataset.

An illustration of PerSign showing the process of the application and also
illustrating the overall architecture of the system

##  1\. Introduction

About $13$ million people in Bangladesh are suffering from different degrees
of hearing loss, of which $3$ million have hearing disability (Alauddin and
Joarder, 2004). There are around $1$ million using Bangladeshi Sign Language
(BdSL) in their everyday life(Syed Tauhid Ahmed, 2016). While communicating
with a signer , there are two major tasks for a non-signer: (i) understanding
the signs and (ii) expressing the signs. Researchers made impressive
contributions to task (i) by developing sign letters222signs that represent
letters only. recognition techniques from images (Fig. 1 ). Several works has
been proposed for BdSL letters classification via machine learning techniques
(Islam et al., 2022; Rahim et al., 2022; Miah et al., 2022; Hasan et al.,
2021; Khatun et al., 2021; Talukder et al., 2021; Hoque et al., 2020, 2018).
Task (ii) still has less research attention since it is a difficult process
for non-signers. A naive and tiresome approach to expressing signs is to use
flashcards with signs and symbols. Being inspired by that, (Shishir et al.,
2020) proposed a system that generates symbols of signs using generative
adversarial networks (GANs) (Goodfellow et al., 2014). However, their work
only produces symbols of signs—which may raise questions regarding the
necessity of such a system. Another way is by using animated avatars of signs,
i.e. (Kipp et al., 2011); but there is no such system for BdSL. Moreover, an
avatar-based system does not provide a realistic environment for
communication.

All the above scenarios inspired us to make task (ii) more realistic yet
effortless. Hence, we introduced PerSign: Personalized Bangladeshi Sign
Letters Synthesis which converts the image of a user into an image showing
signs while keeping the person’s profile unchanged. Fig. 1 explains the
working pipeline of our prototype. A user first uploads their profile photo
($I_{P}$) only once to our system. The $I_{P}$ must follow a specific rule of
showing hand and palm (as shown in a). After that, the user inserts the
desired letter ($L$) to be expressed (e.g.  in b). Our system converts $I_{P}$
into $I_{L}$ by considering $L$. This can be seen as $I_{L}\leftarrow
I_{P}+L$, where $I_{L}$ contains the same person in $I_{P}$ with unchanged
face, skin tone, attire, and background (c). In that case, the person does not
need any expertise in sign language. We believe, a signer will feel a natural
environment if $I_{L}$ is shown, thus, making the communication more realistic
and affectionate.

Do we really need such a system?—in order to address this question, we
performed a survey on a group of $6$ guardians and $11$ teachers of deaf
children—who were also signers—regarding the necessity of our system. We let
participants upload profile images to PerSign and asked them to rate the
results on a scale of \small1⃝ to \small5⃝ according to Likert rating method
(Likert, 1932), with \small1⃝ being not necessary at all and \small5⃝ being
very necessary. Out of total $17$ participants, $13$ and $3$ rated PerSign
with \small5⃝ and \small4⃝ respectively with an average rating of $4.705$.
Most of the sign language teachers commented that personalized signs are very
helpful for general people to get closer to signers, especially when it comes
to children.

##  2\. Implementation

We employed the Generative Adversarial Network (Fig. 1 ), an unsupervised deep
learning technique that automatically learns the patterns from datasets in
order for the model to produce a new output (Goodfellow et al., 2014). Our
problem lies under a sub-domain of GAN which image-to-image translation (Isola
et al., 2017) and we adopted GestureGAN (Tang et al., 2018)—a gesture-to-
gesture translation method—to implement a prototype system. For this purpose,
we built a dataset of images with hand gestures of arbitrary poses, sizes, and
backgrounds. Since we needed a paired dataset $\\{I_{P},I_{L}\\}$ we could not
reuse any of the existing unpaired ones. For localizing gestures and
appearances, we exploited OpenPose (Cao et al., 2019) to make the skeleton of
the hands and face of the images and store the pair of input and output
images. We then trained and tested the GestureGAN model with our dataset, and
constructed our system’s working prototype. Fig. 2 presents more results of
PerSign.

Figure 2. Result analysis. (a) input profile image. (b) generated image with
signs. (c) zoomed view of faces from input (left) and output (right). We can
see, that the face is retained.

A picture of the result obtained exclusively showing the face and hand-
gestures.

##  3\. Conclusion and Future work

In this poster, we proposed a framework—PerSign—for synthesizing Bangladeshi
Sign Letters that can be personalized. Through our method, anyone will be able
to communicate as a signer without having any expertise in BdSL. We built our
own dataset and exploited GestureGAN method to accomplish the task. Our work
is still in progress and we gathered comments from the participants during the
survey to find areas for improvement. Most of the users recommended applying
the technique for a sequence of images to render video for a stream of input
letters. We can achieve better results by increasing the number of diverse
examples in our dataset. Though our dataset has samples for all BdSL letters,
some specific letters need to be treated carefully because of their
similarities in patterns. Our final aim is to merge task (i) and (ii) into a
single system to provide an one-stop solution for two-way communication. Some
of the users also suggested extending our work for full gestures, rather than
letters only. We are also intent on implementing better GUI with voice input.
Last but not the least, we plan to conduct a thorough evaluation from experts
and signers. We believe, this poster opens new avenues in sign language for
further research.

## References

  * (1)
  * Alauddin and Joarder (2004) Mohammad Alauddin and Abul Hasnat Joarder. 2004.  Deafness in bangladesh.  In _Hearing Impairment_. Springer, 64–69.
  * Cao et al. (2019) Z. Cao, G. Hidalgo Martinez, T. Simon, S. Wei, and Y. A. Sheikh. 2019.  OpenPose: Realtime Multi-Person 2D Pose Estimation using Part Affinity Fields.  _IEEE Transactions on Pattern Analysis and Machine Intelligence_ (2019).
  * Goodfellow et al. (2014) Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio. 2014\.  Generative adversarial nets.  _Advances in neural information processing systems_ 27 (2014).
  * Hasan et al. (2021) SK Nahid Hasan, Md Jahid Hasan, and Kazi Saeed Alam. 2021\.  Shongket: A Comprehensive and Multipurpose Dataset for Bangla Sign Language Detection. In _2021 International Conference on Electronics, Communications and Information Technology (ICECIT)_. IEEE, 1–4.
  * Hoque et al. (2020) Oishee Bintey Hoque, Mohammad Imrul Jubair, Al-Farabi Akash, and Saiful Islam. 2020.  BdSL36: A Dataset for Bangladeshi Sign Letters Recognition. In _Proceedings of the Asian Conference on Computer Vision (ACCV) Workshops_.
  * Hoque et al. (2018) Oishee Bintey Hoque, Mohammad Imrul Jubair, Md. Saiful Islam, Al-Farabi Akash, and Alvin Sachie Paulson. 2018.  Real Time Bangladeshi Sign Language Detection using Faster R-CNN. In _2018 International Conference on Innovation in Engineering and Technology (ICIET)_. 1–6.  https://doi.org/10.1109/CIET.2018.8660780
  * Islam et al. (2022) Sohal Islam, Showni Rudra Titli, Kazi Arham Kabir, Md Abdullah Al Hossain, and Md Altaf Hossain. 2022\.  Improving Real-time Hand Gesture Recognition System for Translation: Sensor Development. In _2022 17th Annual System of Systems Engineering Conference (SOSE)_. IEEE, 254–259.
  * Isola et al. (2017) Phillip Isola, Jun-Yan Zhu, Tinghui Zhou, and Alexei A Efros. 2017.  Image-to-image translation with conditional adversarial networks. In _Proceedings of the IEEE conference on computer vision and pattern recognition_. 1125–1134.
  * Khatun et al. (2021) Ayesha Khatun, Mohammad Sajid Shahriar, Md Hasibul Hasan, Krishna Das, Sabbir Ahmed, and Md Sakibul Islam. 2021.  A Systematic Review on the Chronological Development of Bangla Sign Language Recognition Systems. In _2021 Joint 10th International Conference on Informatics, Electronics & Vision (ICIEV) and 2021 5th International Conference on Imaging, Vision & Pattern Recognition (icIVPR)_. IEEE, 1–9.
  * Kipp et al. (2011) Michael Kipp, Alexis Héloir, and Quan Nguyen. 2011\.  Sign Language Avatars: Animation and Comprehensibility. 113–126.  https://doi.org/10.1007/978-3-642-23974-8_13
  * Likert (1932) Rensis Likert. 1932\.  A technique for the measurement of attitudes.  _Archives of psychology_ (1932).
  * Miah et al. (2022) Abu Saleh Musa Miah, Jungpil Shin, Md Al Mehedi Hasan, and Md Abdur Rahim. 2022.  BenSignNet: Bengali Sign Language Alphabet Recognition Using Concatenated Segmentation and Convolutional Neural Network.  _Applied Sciences_ 12, 8 (2022), 3933.
  * Rahim et al. (2022) Md Abdur Rahim, Jungpil Shin, and Keun Soo Yun. 2022\.  Soft Voting-based Ensemble Model for Bengali Sign Gesture Recognition.  _Annals of Emerging Technologies in Computing (AETiC)_ 6, 2 (2022).
  * Shishir et al. (2020) Fairuz Shadmani Shishir, Tonmoy Hossain, and Faisal Muhammad Shah. 2020.  Esharagan: An approach to generate disentangle representation of sign language using infogan. In _2020 IEEE Region 10 Symposium (TENSYMP)_. IEEE, 1383–1386.
  * Syed Tauhid Ahmed (2016) M. A. H. Akhand Syed Tauhid Ahmed. 2016\.  Bangladeshi Sign Language Recognition using Fingertip Position. In _International Conference on Medical Engineering, Health Informatics and Technology (MediTec)_.
  * Talukder et al. (2021) Dipon Talukder, Fatima Jahara, Suvadra Barua, and Md Mokammel Haque. 2021.  OkkhorNama: BdSL Image Dataset For Real Time Object Detection Algorithms. In _2021 IEEE Region 10 Symposium (TENSYMP)_. IEEE, 1–6.
  * Tang et al. (2018) Hao Tang, Wei Wang, Dan Xu, Yan Yan, and Nicu Sebe. 2018.  GestureGAN for Hand Gesture-to-Gesture Translation in the Wild. In _ACM MM_.