# Nonlinear approximation spaces for inverse problems

Albert Cohen, Matthieu Dolbeault11footnotemark: 1, Olga Mula, Agustin
Somacal11footnotemark: 1 Laboratoire Jacques-Louis Lions, Sorbonne Université,
4 place Jussieu, 75005 Paris, France<EMAIL_ADDRESS><EMAIL_ADDRESS>and agustin.somacal@sorbonne-
universite.fr)Technical University of Eindhoven, Netherlands<EMAIL_ADDRESS>
###### Abstract

This paper is concerned with the ubiquitous inverse problem of recovering an
unknown function $u$ from finitely many measurements possibly affected by
noise. In recent years, inversion methods based on linear approximation spaces
were introduced in [7, 27] with certified recovery bounds. It is however known
that linear spaces become ineffective for approximating simple and relevant
families of functions, such as piecewise smooth functions that typically occur
in hyperbolic PDEs (shocks) or images (edges). For such families, nonlinear
spaces [18] are known to significantly improve the approximation performance.
The first contribution of this paper is to provide with certified recovery
bounds for inversion procedures based on nonlinear approximation spaces. The
second contribution is the application of this framework to the recovery of
general bidimensional shapes from cell-average data. We also discuss how the
application of our results to $n$-term approximation relates to classical
results in compressed sensing.

##  1 Introduction

###  1.1 The recovery problem

In this paper, we treat the following state estimation problem in a general
Banach space $V$. We want to recover an approximation to an unknown function
$u\in V$ from data given by $m$ observations

$z_{i}\coloneqq\ell_{i}(u)+\eta_{i},\quad i=1,\dots,m,$ (1.1)

where $\ell_{i}:V\mapsto\mathbb{R}$ are known measurement functionals, and
$\eta_{i}$ is additive noise. The functionals $\ell_{i}$ often correspond to
the response of a physical measurement device but they can have a different
interpretation depending on the application. Their behavior can be linear (in
which case the $\ell_{i}$ are linear functionals from $V^{\prime}$, the dual
of $V$) or nonlinear. This type of recovery problem is clearly ill-posed when
the dimension of $V$ exceeds $m$. It arises ubiquitously in sampling and
inverse problem applications where $V$ is infinite dimensional (to name a few,
see [1, 3, 22, 26]).

One natural strategy to address this difficulty is to search for a recovery of
$u$ by an element of a low-dimensional reconstruction space $V_{n}\subset V$.
The space $V_{n}$ could be either an $n$-dimensional linear subspace, or more
generally a nonlinear approximation space parametrized by $n$ degrees of
freedom, with $n\leqslant m$.

In order to obtain quantitative results for such recovery procedures, it is
necessary to possess additional information about $u$, usually as an
assumption that $u$ belongs to a certain model class ${\cal K}$ contained in
$V$. The approximation space $V_{n}$ is chosen in order to collectively
approximate the elements of ${\cal K}$ as well as possible, in the sense that

$\text{dist}({\cal K},V_{n})_{V}\coloneqq\max_{u\in{\cal K}}\min_{v\in
V_{n}}||u-v||_{V}$

is as small as possible for moderate values of $n$.

Numerous theoretical results and numerical algorithms have been proposed in
several fields to study and solve the above recovery problem (below we recall
some relevant results). However, to the best of our knowledge, they all
involve at least one or several of the following assumptions:

  * •

The $\ell_{i}$ are _linear_ functionals,

  * •

$V_{n}$ is a _linear_ (or affine) subspace of $V$,

  * •

$V$ is a _Hilbert_ space,

  * •

The model class ${\cal K}$ is a _ball in a smoothness space_ , e.g., a unit
ball in Lipschitz, Sobolev, or Besov spaces. Results involving this type of
model classes have been intensively studied in the field of optimal recovery
(see [9, 30, 32]).

The goal of this paper is to develop and analyze inversion procedures that do
not require any of the above assumptions. Our analysis and numerical
algorithms can thus be applied to virtually any recovery problem. The starting
point of our development is based on algorithms introduced for inverse state
estimation using reduced order models of parametrized Partial Differential
Equations (PDEs). We next recall the specific framework. The presentation will
also serve to explain more in depth the motivations leading to propose the
present generalization.

###  1.2 State estimation with reduced models for parametrized PDE’s

A relevant scenario in inverse state estimation is when the model class ${\cal
K}$ is given by the set of solutions to some parameter-dependent PDE of the
general form

${\cal P}(u,y)=0,$ (1.2)

where ${\cal P}$ is a differential operator, $y$ a vector of parameters
ranging in some domain $Y$ in $\mathbb{R}^{d}$, and $u$ is the solution. If
well-posedness holds in some Banach space $V$ for each $y\in Y$, we denote by
$u(y)\in V$ the corresponding solution for the given parameter value $y$ and
by

${\cal M}\coloneqq\\{u(y)\;:\;y\in Y\\},$

the solution manifold.

In inverse state estimation, we take ${\cal K}={\cal M}$ for the model class
so the unknown $u$ to recover belongs to ${\cal M}$. However, the parameter
$y$ that satisfies $u=u(y)$ is unknown, so we cannot solve the forward problem
(1.2) to approximate $u$. Instead, we must approximate $u$ from the partial
observational data (1.1), and the knowledge of the model class ${\cal K}={\cal
M}$.

For the manifold ${\cal M}$, efficient approximation spaces $V_{n}$ are
usually obtained by reduced modelling techniques. In their most simple format,
reduced models consist into linear spaces $(V_{n})_{n\geqslant 0}$ with
$\dim(V_{n})=n$. The ideal benchmark in this linear approximation setting is
provided by the Kolmogorov $n$-width

$d_{n}({\cal M})_{V}\coloneqq\inf_{\dim(V_{n})\leqslant n}{\rm dist}({\cal
M},V_{n})_{V},$

which describes the optimal approximation performance achievable by an
$n$-dimensional space over the set ${\cal M}$.

Apart from very simplified cases, the space $V_{n}$ achieving the above
infimum is usually out of reach. Practical model reduction techniques such as
polynomial approximation in the parametrized domain [15, 16, 39] or reduced
bases [19, 25, 28, 37, 41] construct spaces $V_{n}$ that are “suboptimal yet
good”. In particular, the reduced basis method, which generates $V_{n}$ by a
specific selection of particular solution instances $u^{1},\dots,u^{n}\in{\cal
M}$, has been proved to have approximation error ${\rm dist}({\cal
M},V_{n})_{V}$ that decays with the same polynomial or exponential rates as
$d_{n}({\cal M})_{V}$, and in that sense are close to optimal [17].

###  1.3 The PBDW method

We take the Parametrized Background Data Weak (PBDW) method as a starting
point for our analysis. The PBDW method, first introduced in [27], as well as
several extensions, has been the object of a series of works [7, 8, 13, 14] on
its optimality properties as a recovery algorithm. It has also been used for
different practical applications, see [3, 22, 24]. We refer to [31] for an
overview of the state of the art on this approach, and its connections with
different fields. For our current purposes, it will suffice to recall the
first version of the algorithm, which is the goal of this section.

The PBDW method uses a linear approximation space $V_{n}$ of dimension
$n\leqslant m$. Usually this space is a reduced model in applications. It is
assumed that the $\ell_{i}$ are continuous linear functionals, that is
$\ell_{i}\in V^{\prime}$, and that $V$ is a Hilbert space. Then, introducing
the Riesz representers $\omega_{i}\in V$ such that
$\ell_{i}(v)=\langle\omega_{i},v\rangle_{V}$, the data of the noise-free
observation

$\ell(u)\coloneqq(\ell_{1}(u),\dots,\ell_{m}(u)),$

is equivalent to that of the orthogonal projection $w=P_{W}u$ on the Riesz
measurement space

$W\coloneqq{\rm span}\\{\omega_{1},\dots,\omega_{m}\\}.$

Assuming linear independence of the $\ell_{i}$, this space has dimension $m$.
A critical quantity is the number

$\mu=\mu(V_{n},W)\coloneqq\max_{v\in V_{n}}\frac{\|v\|_{V}}{\|P_{W}v\|_{V}},$
(1.3)

that describes the “stability” of the description of an element of $V_{n}$ by
its projection onto $W$, and may be thought of as the inverse cosine of the
angle between $W$ and $V_{n}$. In particular, this quantity is finite only
when $n\leqslant m$. It can be explicitly computed as the inverse of the
smallest singular value of a cross-grammian matrix between orthonormal bases
of $V_{n}$ and $W$ (see [7, 31]).

The PBDW method consists in solving the minimization problem

$\min_{v\in V_{w}}\min_{\tilde{v}\in V_{n}} \|v-\tilde{v}\|_{V},$

where $V_{w}\coloneqq w+W^{\perp}$ is the set of all states $v$ such that
$P_{W}v=w$. We denote by $(u^{*},\tilde{u})\in V_{w}\times V_{n}$ the
minimizing pair, which is unique when $\mu<\infty$, and can be computed by
solving an $n\times n$ linear system. The function $\tilde{u}$ may be seen as
a particular best-fit estimator of $u$ on $V_{n}$, since it is also defined by

$\tilde{u}\coloneqq{\rm argmin}\\{\|P_{W}v-w\|_{V}\;:\;v\in V_{n}\\}.$

The function $u^{*}$ can be derived from $\tilde{u}$ by the correction
procedure

$u^{*}\coloneqq\tilde{u}+(w-P_{W}\tilde{u}),$

which shows that $u^{*}\in V_{n}+W$. It may be thought of as a generalized
interpolation estimator, since it agrees with the observed data
($P_{W}u^{*}=P_{W}u$). In the case of noise-free data, it is proved in [7, 27]
that these estimators satisfy the recovery bounds

$\|u-\tilde{u}\|_{V}\leqslant\mu\min_{v\in V_{n}}\|u-v\|_{V}\quad{\rm and}
\quad\|u-u^{*}\|_{V}\leqslant\mu\min_{v\in V_{n}\oplus(V_{n}^{\perp}\cap
W)}\|u-v\|_{V}.$

These bounds reflect a typical trade-off in the choice of the reduced basis
space, since making $n$ larger has both effect of decreasing the approximation
error $\min_{v\in V_{n}}\|u-v\|_{V}$ and increasing the stability constant
$\mu=\mu(V_{n},W)$.

When the PBDW method is applied to noisy data, amounting in observing a
perturbed version $\overline{w}$ of $w=P_{W}u$, the recovery bounds remain
valid up to the additional term $\mu\|w-\overline{w}\|_{V}$. In summary, one
has for both estimators

$\max\\{\|u-\tilde{u}\|_{V},\|u-u^{*}\|_{V}\\}\leqslant\mu(e_{n}(u)+\kappa),$
(1.4)

where

$e_{n}(u)\coloneqq\min_{v\in V_{n}}\|u-v\|_{V}$

is the reduced model approximation error and
$\kappa\coloneqq\|w-\overline{w}\|_{V}$ is the noise error measured in the
space $W$. Note that since the additive perturbations $\eta_{i}$ are applied
to the data $\ell_{i}(u)$, a natural model for the measurement noise is to
assume a bound of the form

$\|\eta\|_{p}\leqslant\varepsilon,$ (1.5)

for the vector $\eta=(\eta_{1},\dots,\eta_{m})$, typically in the max norm
$p=\infty$ or euclidean norm $p=2$. Therefore, one has
$\kappa\leqslant\beta\varepsilon$, where

$\beta\coloneqq\max_{v\in W}\frac{\|v\|_{V}}{\|\ell(v)\|_{p}},$

resulting in a bound of the form $\mu e_{n}(u)+\mu\beta\varepsilon$ for both
estimators.

###  1.4 Towards nonlinear approximation spaces

The simplicity of the PBDW method and its variants comes together with a
fundamental limitation on its performance: it is by essence a linear
reconstruction method with recovery bounds tied to the approximation error
$e_{n}(u)$. When the only prior information is that the unknown function $u$
belongs to a class ${\cal K}$, with ${\cal K}={\cal M}$ the solution manifold
in the case of parametric PDEs, its best performance over ${\cal K}$ is thus
limited by the $n$-width $d_{n}({\cal K})_{V}$ and in turn by $d_{m}({\cal
K})_{V}$ since $n\leqslant m$.

In several simple yet relevant settings, it is known that $n$-widths have poor
decay with $n$. One instance is when the class ${\cal K}$ contains piecewise
smooth states, with a state-dependent location of jump discontinuities. As an
elementary example, one can easily check that if $V=L^{2}([0,1])$ and ${\cal
K}$ is the set all indicator functions $u=\raise
1.29167pt\hbox{\large$\chi$}_{[a,b]}$ with $a,b\in[0,1]$, one has $d_{n}({\cal
K})_{V}\sim n^{-1/2}$. This decay is of course even slower for more general
classes of piecewise smooth function in higher dimension, see in particular
[5, Chapter 3, equation (3.76)]. Such functions are typical in parametrized
hyperbolic PDEs, due to the presence of shocks with positions that differ when
parameters entering the velocity vary. We refer to [4, 8, 20, 23, 33, 40] for
other examples of parametric PDEs whose solution manifold has slow Kolmogorov
$n$-width decay.

For such classes of functions, nonlinear approximation methods are well known
to perform significantly better than their linear counterparts. Typical
representatives of such methods include approximation by rational fractions,
free knot splines or adaptive finite elements, best $n$-term approximation in
a basis or dictionary, neural network or various tensor formats. In these
instances the space $V_{n}$ still depends on $n$ or ${\cal O}(n)$ parameters
but is not anymore a linear space. We refer to [18] for a general introduction
on the topic of nonlinear approximation.

###  1.5 Objective and outline

The objective of this paper is to study the natural extensions of the PBDW
method to such nonlinear approximation spaces and identify the basic
structural properties that lead to near optimal recovery estimates similar to
(1.4).

We begin in §2 by considering the most general setting where $V$ is a Banach
space, $V_{n}$ a nonlinear approximation family, and the $\ell_{i}$ are
functionals defined on $V$ that are not necessarily linear, but Lipschitz
continuous, that is

$\|\ell(v)-\ell(\tilde{v})\|_{Z}\leqslant\alpha_{Z}\|v-\tilde{v}\|_{V},\quad
v,\tilde{v}\in V.$ (1.6)

Here $\|\cdot\|_{Z}$ can be any given norm defined over $\mathbb{R}^{m}$ with
the constant $\alpha_{Z}$ depending on this choice of norm. In this framework,
we discuss the best-fit estimation procedure that consists in minimizing the
distance to the observed data in a given norm $\|\cdot\|_{Z}$.

Our main structural assumption on $V_{n}$ is the following inverse stability
property: the reduced model is stable with respect to the measurement
functionals if there exists a finite constant $\mu_{Z}$ such that

$\|v-\tilde{v}\|_{V}\leqslant\mu_{Z}\|\ell(v)-\ell(\tilde{v})\|_{Z},\quad
v,\tilde{v}\in V_{n}.$ (1.7)

The stability constant $\mu_{Z}$ depends on the $Z$ norm and plays a role
similar to that of $\mu$ in the linear case. In particular, we show that this
constant is finite only if $n\leqslant m$. The resulting estimator $\tilde{u}$
is then proved to satisfy a general recovery bound of the form

$\|u-\tilde{u}\|_{V}\leqslant C_{1}e_{n}(u)+C_{2}\|\eta\|_{p},$

where $e_{n}(u)\coloneqq\min_{v\in V_{n}}\|u-v\|_{V}$ is the nonlinear reduced
model approximation error, $\|\eta\|_{p}$ the level of measurement noise in
$\ell^{p}$ norm, and the constants $C_{1}$ and $C_{2}$ depend on $\alpha_{Z}$
and $\mu_{Z}$.

In §3, we consider the more particular setting where the $\ell_{i}$ are linear
functionals. Then, we show that constants $C_{1}$ and $C_{2}$ are each
minimized by a different choice of norm $\|\cdot\|_{Z}$, resulting in two
different best fit estimators $\tilde{u}$, as already observed in [6] in the
case of linear reduced models. This particular setting also allows us to
introduce a generalized interpolation estimator $u^{*}$ and establish similar
recovery estimates for $\|u-u^{*}\|_{V}$.

We next apply our framework to the inverse problem that consists in recovering
a general shape $\Omega$, identified to its characteristic function $\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}$, based on cell average data

$a_{T}(\Omega)\coloneqq\frac{1}{|T|}\int_{T}\raise
1.29167pt\hbox{\large$\chi$}_{\Omega},\quad T\in{\cal T},$

where ${\cal T}$ is a fixed cartesian mesh. One motivation for this problem is
the design of finite volume schemes for the computation of solutions to
transport PDEs on such meshes.

We first discuss in §4 the best estimation rate in terms of the mesh size $h$
that can be achieved by standard linear reconstructions, and which is
essentially that of piecewise constant approximations, that is ${\cal
O}(h^{1/q})$ regardless of the smoothness of the boundary $\partial\Omega$.
This intrinsic limitation is due to the presence of the jump discontinuity
that is not well resolved by the mesh.

We then discuss in §5 a local recovery strategy based on a nonlinear
approximation space $V_{n}$ that consists of characteristic functions of half-
planes which can fit the boundary of $\Omega$ at a subcell resolution level,
as already proposed in [2, 34, 35, 36]. One main result, whose proof is given
in an appendix, is that this approximation space is stable in the sense of
(1.7) with respect to cell average measurements on a stencil of $3\times 3$
squares. In turn, if $\Omega$ has a $C^{2}$ boundary, the recovered shape
$\tilde{\Omega}$ is proved to satisfy an estimate of the form

$\|\raise 1.29167pt\hbox{\large$\chi$}_{\Omega}-\raise
1.29167pt\hbox{\large$\chi$}_{\tilde{\Omega}}\|_{L^{q}}\leqslant Ch^{2/q},$

where $h$ is the mesh size, which cannot be achieved by any linear
reconstruction. This paves the way to higher order reconstruction methods for
smoother boundaries by using local nonlinear approximation spaces with curved
boundaries and larger stencils.

Finally, we discuss in §6 the application of our results to the recovery of
large vectors of size $N$ from $m<N$ linear measurements, up to the error of
best $n$-term approximation. This problem is well-known in compressed sensing
[10, 21], and was in particular studied in [12] which discusses the importance
of the recovery norm $\|\cdot\|_{V}$ to understand if near-optimal recovery
bounds can be achieved with $m$ not much larger than $n$. We show that the
structural assumptions identified in our general setting are naturally related
to the so-called null space property introduced in [12].

##  2 Nonlinear reduction of inverse problems

###  2.1 A general framework

In full generality we are interested in recovering functions $u$ in a general
Banach space $V$ with norm $\|\cdot\|_{V}$, from the measurement vector
$z=(z_{1},\dots,z_{m})\in\mathbb{R}^{m}$ given by (1.1). A recovery (or
inversion) map

$z\to R(z),$

takes this vector to an approximation $R(z)$ of $u$. We are interested in
controlling the recovery error $\|u-R(z)\|_{V}$.

To build the recovery map $R$, we use a nonlinear approximation space of
dimension $n$ is a family of functions that can be described by $n$
parameters. Loosely speaking, this means that there exists a set
$S\subset\mathbb{R}^{n}$ and a continuous map $\varphi:S\to V$ such that

$V_{n}\coloneqq\\{\varphi(x)\;:\;x\in S\\}.$

Note that this definition covers the case of an $n$ dimensional linear
subspace since we can choose $S=\mathbb{R}^{n}$ and $\varphi$ a linear map.

Our main assumptions are the Lipschitz stability of the functionals $\ell_{i}$
over the whole space $V$ and their inverse Lipschitz stability over the
nonlinear approximation space $V_{n}$, expressed by (1.6) and (1.7),
respectively. Note that since $\mathbb{R}^{m}$ is finite dimensional, the norm
$\|\cdot\|_{Z}$ that is chosen in $\mathbb{R}^{m}$ to express these properties
could be arbitrary up to a modification of the stability constants
$\alpha_{Z},\mu_{Z}$. These constants can be optimally defined as

$\alpha_{Z}=\sup_{v_{1},v_{2}\in
V}\frac{\|\ell(v_{1})-\ell(v_{2})\|_{Z}}{\|v_{1}-v_{2}\|_{V}},$

and

$\mu_{Z}=\sup_{v_{1},v_{2}\in
V_{n}}\frac{\|v_{1}-v_{2}\|_{V}}{\|\ell(v_{1})-\ell(v_{2})\|_{Z}}.$

Note that one always has $\alpha_{Z}\mu_{Z}\geqslant 1$.

######  Remark 2.1.

Note that when $V_{n}$ is an $n$-dimensional space and the $\ell_{i}$ are
linear functionals, the quantity $\mu_{Z}$ may be rewritten as

$\mu_{Z}=\max_{v\in V_{n}}\frac{\|v\|_{V}}{\|\ell(v)\|_{Z}}.$

As discussed further, the quantity $\mu$ defined in (1.3) for the analysis of
the PBDW method is an instance of $\mu_{Z}$ corresponding to a particular
choice of norm $\|\cdot\|_{Z}$. Assuming the $\ell_{i}$ are independent
functionals, one easily checks that finiteness of this quantity imposes that
$n\leqslant m$. Indeed, if $n>m$, there exists a non-trivial $v\in
V_{n}\cap{\cal N}$, where

${\cal N}\coloneqq\\{v\;:\;\ell(v)=0\\}$

is the null space of the measurement map that has codimension $m$, and
therefore $\mu_{Z}$ is infinite.

######  Remark 2.2.

The restriction $n\leqslant m$ is also needed for nonlinear spaces $V_{n}$ and
measurement $\ell$, under assumptions expressing that $m$ and $n$ are local
dimensions. More precisely, assume that the map $\varphi$ defining $V_{n}$ is
differentiable at some $x_{0}$ in the interior of $S$, that $\ell$ is
differentiable at $v_{0}=\varphi(x_{0})$, and that both tangent maps have full
rank at these points, that is,

${\rm dim}(d\varphi_{x_{0}}(\mathbb{R}^{n}))=n\quad{\rm and}\quad{\rm
dim}(d\ell_{v_{0}}(V))=m.$

Then, by taking $v_{1}=v_{0}$ and $v_{2}=\varphi(x_{0}+tx)$ in the quotient
that defines $\mu_{Z}$, and letting $t\to 0$ for arbitrary
$x\in\mathbb{R}^{n}$, one finds that

$\mu_{Z}\geqslant\max_{v\in
d\varphi_{x_{0}}(\mathbb{R}^{n})}\frac{\|v\|_{V}}{\|d\ell_{v_{0}}(v)\|_{Z}},$

and therefore it is infinite if $n\leqslant m$, by the same argument as in the
previous remark.

###  2.2 The best fit estimator

We define a first recovery map $z\mapsto\tilde{u}=R(z)$ as the best fit
estimator in the $Z$ norm

$\tilde{u}\coloneqq{\rm argmin}\\{\|z-\ell(v)\|_{Z}\;:\;v\in V_{n}\\}.$ (2.8)

The existence of such a minimizer is trivial if the space $V_{n}$ and the
measurement map $\ell$ are linear. It can also be ensured in the nonlinear
case under additional assumptions, for example compactness of the set $S$
defining the nonlinear space $V_{n}$, which will be the case in the
application to shape recovery discussed in §5. If the minimizer does not
exist, we may consider a near minimizer, that is $\tilde{u}\in V_{n}$
satisfying

$\|z-\ell(\tilde{u})\|_{Z}\leqslant C\|z-\ell(v)\|_{Z},\quad v\in V_{n},$

for some fixed $C>1$. Inspection of the proofs of our main results below
reveals that similar recovery bounds can be obtained for such a near
minimizer, up to the multiplicative constant $C$.

Recall that our assumption (1.5) on the noise model is a control on
$\|\eta\|_{p}$ for some $1\leqslant p\leqslant\infty$. For this value of $p$,
we introduce the quantity

$\beta_{Z}\coloneqq\max_{z\in\mathbb{R}^{m}}\frac{\|z\|_{Z}}{\|z\|_{p}}$

We are now in position to state a recovery bound in this general framework.

######  Theorem 2.3.

The best fit estimator $\tilde{u}$ from (2.8) satisfies the estimate

$\|u-\tilde{u}\|_{V}\leqslant C_{1}e_{n}(u)+C_{2}\|\eta\|_{p},$ (2.9)

where $C_{1}\coloneqq 1+2\alpha_{Z}\mu_{Z}$ and $C_{2}\coloneqq
2\beta_{Z}\mu_{Z}$.

Proof: Consider any $v\in V_{n}$ and write

$\|u-\tilde{u}\|_{V}\leqslant\|u-v\|_{V}+\|v-\tilde{u}\|_{V}\leqslant\|u-v\|_{V}+\mu_{Z}\|\ell(v)-\ell(\tilde{u})\|_{Z},$

where we have used (1.7). On the other hand, the minimizing property of
$\tilde{u}$ ensures that

$\|\ell(v)-\ell(\tilde{u})\|_{Z}\leqslant\|z-\ell(v)\|_{Z}+\|z-\ell(\tilde{u})\|_{Z}\leqslant
2\|z-\ell(v)\|_{Z}.$

Furthermore, using the stability (1.6) of $\ell$ and the definition of
$\beta_{Z}$, we have

$\|z-\ell(v)\|_{Z}\leqslant\|\ell(v)-\ell(u)\|_{Z}+\|\eta\|_{Z}\leqslant\alpha_{Z}\|u-v\|+\beta_{Z}\|\eta\|_{p}.$

Combining the three estimates, we reach

$\|u-\tilde{u}\|_{V}\leqslant(1+2\alpha_{Z}\mu_{Z})\|u-v\|_{V}+2\beta_{Z}\mu_{Z}\|\eta\|_{p},$

which gives (2.9) by optimizing over $v\in V_{n}$. $\Box$

The constants $C_{1}$ and $C_{2}$ in the above recovery estimate depend on the
choice of norm $\|\cdot\|_{Z}$. Note that they are invariant when this norm is
scaled by a factor $t>0$, since this has the effect of multiplying
$\alpha_{Z}$ and $\beta_{Z}$ by $t$ and dividing $\mu_{Z}$ by $t$, which is
consistant with the fact that the resulting estimator $\tilde{u}$ is left
unchanged by such a scaling. In the next section we show, in the particular
setting of linear measurements, that specific choices of $\|\cdot\|_{Z}$ can
be used to minimize $C_{1}$ or $C_{2}$. This setting also allows us to
introduce and study a generalized interpolation estimator, which is not
relevant to the present section since the nonlinear measurement map $\ell$ is
not assumed to be surjective: in the presence of noise, there might exist no
$v\in V$ that agrees with the data, in the sense that $z=\ell(u)+\eta$ does
not belong to the range of $\ell$.

##  3 Linear observations

In this section, we assume that the $\ell_{i}\in V^{\prime}$ are independent
linear functionals, still allowing $V_{n}$ to be a general nonlinear space. In
this framework, which contains the example of shape recovery discussed in §5,
one has

$\alpha_{Z}=\max_{v\in V}\frac{\|\ell(v)\|_{Z}}{\|v\|_{V}}$

and

$\mu_{Z}=\max_{v\in V_{n}^{\text{diff}}}\frac{\|v\|_{V}}{\|\ell(v)\|_{Z}},$

where

$V_{n}^{\text{diff}}=V_{n}-V_{n}\coloneqq\\{v_{1}-v_{2}\;:\;v_{1},v_{2}\in
V_{n}\\}.$

In this particular setting, we can identify the norms $\|\cdot\|_{Z}$ that
minimize the constants $C_{1}\coloneqq 1+2\alpha_{Z}\mu_{Z}$ and
$C_{2}\coloneqq 2\beta_{Z}\mu_{Z}$, respectively.

###  3.1 Optimal norms

As $\ell:V\to\mathbb{R}^{m}$ is continuous and surjective, we can define a
norm on $\mathbb{R}^{m}$ through

$\|z\|_{W}=\min\\{\|v\|_{V}\;:\;\ell(v)=z\\}.$ (3.10)

######  Remark 3.1.

If $V$ is a Hilbert space, the minimizer is unique by strict convexity of
$\|\cdot\|_{V}$, and the $m$-dimensional space

$W:=\\{\operatorname{argmin}_{\ell(v)=z}\|v\|_{V},\;z\in\mathbb{R}^{m}\\}$

is exactly the span of the Riesz representers of the observation functionals
$\ell_{i}\in V^{\prime}$. Moreover, denoting $P_{W}$ the orthogonal projection
on $W$, we have

$\|\ell(v)\|_{W}=\|P_{W}v\|_{V},\quad v\in V.$

For this reason, we sometimes refer to $\|\cdot\|_{W}$ as the Riesz norm even
in the case of a more general Banach space.

The following result shows that the choice
$\|\cdot\|_{Z}\coloneqq\|\cdot\|_{W}$ is the one that minimizes the constant
$C_{1}$, while $C_{2}$ is minimized by simply taking the $\ell^{p}$ norm
$\|\cdot\|_{Z}=\|\cdot\|_{p}$.

######  Theorem 3.2.

For any norm $\|\cdot\|_{Z}$, one has

$\alpha_{W}\mu_{W}=\mu_{W}\leqslant\alpha_{Z}\mu_{Z},$

and

$\beta_{p}\mu_{p}=\mu_{p}\leqslant\beta_{Z}\mu_{Z},$

where $(\alpha_{W},\beta_{W},\mu_{W})$ and $(\alpha_{p},\beta_{p},\mu_{p})$
are the triplets $(\alpha_{Z},\beta_{Z},\mu_{Z})$ when
$\|\cdot\|_{Z}\coloneqq\|\cdot\|_{W}$ and $\|\cdot\|_{Z}=\|\cdot\|_{p}$,
respectively.

Proof: One has

$\alpha_{W}=\max_{v\in
V}\frac{\|\ell(v)\|_{W}}{\|v\|_{V}}=\max_{z\in\mathbb{R}^{m}}\max_{\ell(v)=z}\frac{\|z\|_{W}}{\|v\|_{V}}=1,$

and so

$\alpha_{W}\mu_{W}=\mu_{W}=\max_{v\in
V_{n}^{\text{diff}}}\frac{\|v\|_{V}}{\|\ell(v)\|_{W}}\leqslant\max_{v\in
V_{n}^{\text{diff}}}\frac{\|\ell(v)\|_{Z}}{\|\ell(v)\|_{W}}\max_{v\in
V_{n}^{\text{diff}}}\frac{\|v\|_{V}}{\|\ell(v)\|_{Z}}=\max_{v\in
V_{n}^{\text{diff}}}\frac{\|\ell(v)\|_{Z}}{\|\ell(v)\|_{W}}\mu_{Z}.$

We now observe that from the definition of $W$, one has

$\max_{v\in
V_{n}^{\text{diff}}}\frac{\|\ell(v)\|_{Z}}{\|\ell(v)\|_{W}}\leqslant\max_{z\in\mathbb{R}^{m}}\frac{\|z\|_{Z}}{\|z\|_{W}}=\max_{z\in\mathbb{R}^{m}}\max_{\ell(v)=z}\frac{\|z\|_{Z}}{\|v\|_{V}}=\alpha_{Z}.$

We have thus obtained the first claim
$\alpha_{W}\mu_{W}=\mu_{W}\leqslant\alpha_{Z}\mu_{Z}$. For the second claim,
note that we trivially have $\beta_{p}=1$, and so

$\beta_{p}\mu_{p}=\mu_{p}=\max_{v\in
V_{n}^{\text{diff}}}\frac{\|v\|_{V}}{\|\ell(v)\|_{p}}\leqslant\max_{v\in
V_{n}^{\text{diff}}}\frac{\|\ell(v)\|_{Z}}{\|\ell(v)\|_{p}}\max_{v\in
V_{n}^{\text{diff}}}\frac{\|v\|_{V}}{\|\ell(v)\|_{Z}}\leqslant\beta_{Z}\mu_{Z}.$

$\Box$

######  Remark 3.3.

In the particular case where $V$ is a Hilbert space, $V_{n}$ a linear subspace
and $p=2$, it was already observed in [6] that the reconstruction operators
based on the choice $\|\cdot\|_{Z}=\|\cdot\|_{W}$ or
$\|\cdot\|_{Z}=\|\cdot\|_{2}$ are the most stable with respect to the
approximation error and the noise error, respectively. The above result may
thus be seen as a generalization of this state of affairs to the case of
nonlinear subspaces of Banach spaces, and $\ell^{p}$ noise.

###  3.2 The generalized interpolation estimator

Thanks to the surjectivity of $\ell$, we may introduce the space

$V_{z}:=\\{v\in V\;:\;\ell(v)=z\\},$

and consider the minimization problem

$\min_{v\in V_{z}}\min_{\tilde{v}\in V_{n}} \|v-\tilde{v}\|_{V}.$

If $(u^{*},\tilde{u})\in V_{z}\times V_{n}$ is a minimizing pair, the function
$u^{*}$ is given by

$u^{*}=u^{*}(z)\in{\rm argmin}\\{{\rm dist}(v,V_{n})_{V}\;:\;\ell(v)=z\\},$

and is called the generalized interpolation estimator, since it exactly
matches the data.

######  Remark 3.4.

The best fit and generalized interpolation estimation may be thought of as the
two extreme cases, $t\to\infty$ and $t\to 0$, of the penalized estimator

$u_{t}\coloneqq{\rm argmin}\\{\|z-\ell(v)\|_{Z}+t\,{\rm
dist}(v,V_{n})_{V}\\}.$

As explained earlier, the generalized interpolation operator may not be well
defined in the general case where the $\ell_{i}$ are nonlinear. As opposed to
the best fit, or the above penalized estimator $u_{t}$ when $t>0$, the
generalized interpolation estimator does not involve the choice of a
particular norm $Z$.

On the other hand, we see that $\tilde{u}$ is the solution to the problem

$\min_{\tilde{v}\in V_{n}}{\rm dist}(\tilde{v},V_{z})_{V}.$

Observing that

${\rm
dist}(\tilde{v},V_{z})_{V}=\min_{\ell(v)=z}\|\tilde{v}-v\|_{V}=\min_{\ell(v^{\prime})=\ell(\tilde{v})-z}\|v^{\prime}\|_{V}=\|\ell(\tilde{v})-z\|_{W},$

we thus find that $\tilde{u}$ is precisely the best fit estimator for the
Riesz norm $\|\cdot\|_{Z}\coloneqq\|\cdot\|_{W}$.

In the Hilbert space setting, the generalized interpolation estimator $u^{*}$
is therefore the orthogonal projection of this particular best fit estimator
$\tilde{u}$ onto the affine space $V_{z}$. It may thus also be derived from
$\tilde{u}$ by the correction procedure

$u^{*}=\tilde{u}+w-P_{W}\tilde{u},$

where $w=\operatorname{argmin}_{\ell(v)=z}\|v\|_{V}\in W$ is the preimage by
$\ell$ of the measurements $z$. In the noiseless case when $w=P_{W}u$, this
correction can only improve the approximation since it reduces the component
of $u-\tilde{u}$ in the $W$ direction while leaving unchanged the orthogonal
component, and so, in view of Theorems 2.3 and 3.2, we are ensured that

$\|u-u^{*}\|_{V}\leqslant C_{1}e_{n}(u),$

where $C_{1}\coloneqq 1+2\mu_{W}$.

More generally, in the noisy case, and without the assumption that $V$ is a
Hilbert space, there is no guarantee that $u^{*}$ performs better than
$\tilde{u}$, but we still obtain an error estimate on $u^{*}$ that is similar
in nature to that satisfied by $\tilde{u}$.

######  Theorem 3.5.

The generalized interpolation estimator $u^{*}$ satisfies the estimate

$\|u-u^{*}\|_{V}\leqslant C_{1}e_{n}(u)+C_{2}\|\eta\|_{p},$ (3.11)

where $C_{1}\coloneqq 2+2\mu_{W}$ and $C_{2}\coloneqq(1+2\mu_{W})\beta_{W}$.

Proof: Take $\delta\in\operatorname{argmin}_{\ell(v)=\eta}\|v\|_{V}$, so that
$\ell(\delta)=\eta$ and $\|\eta\|_{W}=\|\delta\|_{V}$. For $v$ and $v^{*}$ in
$V_{n}$, decompose

$\|u-u^{*}\|_{V}\leqslant\|u-v\|_{V}+\|v-v^{*}\|_{V}+\|v^{*}-u^{*}\|_{V}.$
(3.12)

For the middle term, using (1.7), we write

$\displaystyle\|v-v^{*}\|_{V}$
$\displaystyle\leqslant\mu_{W}\|\ell(v-v^{*})\|_{W}$
$\displaystyle\leqslant\mu_{W}(\|\ell(v-u)\|_{W}+\|\ell(u-u^{*})\|_{W}+\|\ell(u^{*}-v^{*})\|_{W})$
$\displaystyle\leqslant\mu_{W}(\|v-u\|_{V}+\|\eta\|_{W}+\|u^{*}-v^{*}\|_{V})$

since $\alpha_{W}=1$, so the decomposition (3.12) becomes

$\|u-u^{*}\|_{V}\leqslant(1+\mu_{W})\|u-v\|_{V}+\mu_{W}\|\eta\|_{W}+(1+\mu_{W})\|v^{*}-u^{*}\|_{V}.$

To bound the last term, we optimize over the choice of $v^{*}\in V_{n}$ and
use the definition of $u^{*}$ to obtain

$\inf_{v^{*}\in V_{n}}\|v^{*}-u^{*}\|_{V}={\rm dist}(u^{*},V_{n})\leqslant{\rm
dist}(u+\delta,V_{n})\leqslant{\rm
dist}(u,V_{n})+\|\delta\|_{V}=e_{n}(u)+\|\eta\|_{W}$

since $\ell(u+\delta)=\ell(u)+\eta=z$. Combining the last two estimates and
optimizing over $v\in V_{n}$ gives

$\|u-u^{*}\|_{V}\leqslant(2+2\mu_{W})e_{n}(u)+(1+2\mu_{W})\|\eta\|_{W},$

and the result follows from the definition of $\beta_{W}$. $\Box$

##  4 Shape recovery from cell averages

###  4.1 The shape recovery problem

The problem of reconstructing a function $u$ from its cell averages

$a_{T}(u)\coloneqq\frac{1}{|T|}\int_{T}u,\quad T\in{\cal T},$

where ${\cal T}$ is a partition of the domain $D\subset\mathbb{R}^{d}$ in
which $u$ is defined, appears naturally in two areas:

  * •

In $2d$ or $3d$ image processing, it corresponds to the so-called super-
resolution problem, that is, reconstructing a high resolution image from its
low resolution version defined on the coarse grid ${\cal T}$ of pixels or
voxels.

  * •

In numerical simulation of hyperbolic conservation laws, it plays a central
role when developing finite volume schemes on the computation mesh ${\cal T}$.

Standard reconstruction methods are challenged when the function $u$ exhibits
jump discontinuities which are not well resolved by the partition ${\cal T}$.
Such discontinuities correspond to edges in image processing or shocks in
conservation laws. Here we may focus on the very simple case of characteristic
functions of sets

$u=\raise 1.29167pt\hbox{\large$\chi$}_{\Omega},$

that already carry the main difficulty. Therefore we are facing a problem of
reconstructing a shape $\Omega$ from local averages of $\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}$.

As a simple example we work in the domain $D=[0,1]^{2}$ with a uniform grid
based on square cells of sidelength $h=\frac{1}{L}$ for some $L>1$, therefore
of the form

${\cal T}={\cal
T}_{h}\coloneqq\\{T_{i,j}=[(i-1)h,ih]\times[(j-1)h,jh]\;:\;i,j=1,\dots,L\\}.$

The cardinality of the grid is therefore

$n:=\\#({\cal T})=L^{2}=h^{-2}.$

We consider classes of characteristic functions $\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}$ of sets $\Omega\subset D$ with boundary
of a prescribed Hölder smoothness. The definition of these classes requires
some precision.

######  Definition 4.1.

For $s\geqslant 1$, $0<R<1/2$ and $M>0$, we define the class ${\cal
F}_{s,R,M}$ as consisting of all characteristic functions $\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}$ of domains
$\Omega\subset[R,1-R]^{2}\subset D$ with the following property: for all $x\in
D$ there exists an orthonormal system $(e_{1},e_{2})$ and a function
$\psi\in{\cal C}^{s}$ with $\|\psi\|_{{\cal C}^{s}}\leqslant M$, such that

$y\in\Omega\iff z_{2}\leqslant\psi(z_{1}),$

for any $y=x+z_{1}e_{1}+z_{2}e_{2}$ with $|z_{1}|,|z_{2}|\leqslant R$.

Here, we have used the usual definition

$\|\psi\|_{{\cal C}^{s}}=\sup_{0\leqslant k\leqslant\lfloor
s\rfloor}\|\psi^{(k)}\|_{L^{\infty}([-R,R])}+\sup_{s,t\in[-R,R]}|s-t|^{\lfloor
s\rfloor-s}\Big{|}\psi^{(\lfloor s\rfloor)}(s)-\psi^{(\lfloor
s\rfloor)}(t)\Big{|},$

for the Hölder norm. In the case of integer smoothness, we use the convention
that ${\cal C}^{s}$ denotes functions with Lipschitz derivatives up to order
$s-1$, so that in particular the case $s=1$ corresponds to domains with
Lipschitz boundaries.

######  Remark 4.2.

The condition $\Omega\subset[R,1-R]^{2}$ imposing that $\Omega$ remains away
from the boundary $\partial D$ might be quite restrictive in some
applications; instead, one can assume that the domains $\Omega$ and $D$ are
periodic, or symmetrize $\Omega$ with respect to $\partial D$.

###  4.2 The failure of linear reconstruction methods

The most trivial linear reconstruction method consists in the piecewise
constant approximation

$\tilde{u}=\sum_{T\in{\cal T}}a_{T}(u)\raise
1.29167pt\hbox{\large$\chi$}_{T}.$ (4.13)

The approximation rate of this reconstruction over the class ${\cal
F}_{s,R,M}$ is as follows.

######  Proposition 4.3.

Let $u=\raise 1.29167pt\hbox{\large$\chi$}_{\Omega}\in{\cal F}_{s,R,M}$, its
piecewise constant approximation $\tilde{u}$ by average values on each cell,
defined in (4.13), satisfies

$\|\raise 1.29167pt\hbox{\large$\chi$}_{\Omega}-\tilde{u}\|_{L^{q}}\leqslant
Ch^{\frac{1}{q}}=Cn^{-\frac{1}{2q}},$

where the constant $C$ depends on $R$ and $M$.

Proof: Let $N=\lceil(\sqrt{2}R)^{-1}\rceil$, and partition the domain
$D=[0,1]^{2}$ into $N^{2}$ squares of side $1/N$. Then each subsquare $Q$ is
contained in the set $\\{x+z_{1}e_{1}+z_{2}e_{2},|z_{1}|,|z_{2}|\leqslant
R\\}$ from Definition 4.1, where $x$ is the center of $Q$. Thus
$\partial\Omega$ is the restriction of the graph of an $M$-Lipschitz function
on $Q$, so its arc length is bounded by

$|\partial\Omega\cap Q|\leqslant{\rm diam}(Q)\sqrt{1+M^{2}}\leqslant
2R\sqrt{1+M^{2}}.$

As any curve of arclength $h$ intersects at most four cells from ${\cal T}$,
$\partial\Omega\cap Q$ intersects at most $4\lceil 2R\sqrt{1+M^{2}}/h\rceil$
cells, and summing over all subsquares, $\partial\Omega$ intersects at most
$4N^{2}\lceil 2R\sqrt{1+M^{2}}/h\rceil$ cells. Denoting ${\cal
T}_{\partial\Omega}$ the set of these cells, and observing that $u|_{T}\equiv
a_{T}(u)\in\\{0,1\\}$ for $T\notin{\cal T}_{\partial\Omega}$, we get

$\|\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}-\tilde{u}\|_{L^{q}}^{q}=\sum_{T\in{\cal
T}}\int_{T}|u-a_{T}(u)|^{q}\leqslant\sum_{T\in{\cal
T}_{\partial\Omega}}|T|=h^{2}|{\cal T}_{\partial\Omega}|\leqslant
24\frac{\sqrt{1+M^{2}}}{R}h$

for $h\leqslant R$, and this bound also holds for $h>R$ since $\|\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}-\tilde{u}\|_{L^{q}}^{q}\leqslant 1$.
$\Box$

The next result shows, for the particular case $q=2$, that no better rate can
actually be achieved by any linear method, regardless of the smoothness $s$ of
the boundary. We conjecture that a similar result holds for $1\leqslant
q\leqslant\infty$. This motivates the use of nonlinear recovery methods, which
are the object of the next section.

We recall that the Kolmogorov $n$-width of a compact set $S$ from some Banach
space $V$ is defined by

$d_{n}(S)_{V}:={\rm inf}_{\dim(E)\leqslant n}{\rm dist}(S,E)_{V},$

where ${\rm dist}(S,E)_{V}:=\max_{u\in S}\min_{v\in E}\|u-v\|_{V}$ and the
infimum is taken over all finite dimensional spaces $E$ of dimension at most
$n$.

######  Proposition 4.4.

Let $s\geqslant 1$ be arbitrary. Then for $R$ sufficiently small, and $M$
sufficiently large, there exists $c>0$ such that the Kolmogorov $n$-widths of
the class ${\cal F}_{s,R,M}$ satisfy

$d_{n}({\cal F}_{s,R,M})_{L^{2}}\geqslant cn^{-\frac{1}{4}},\quad n\geqslant
1.$

Proof: The proof of this result relies on similar lower bounds for
dictionaries of $d$-dimensional ridge functions

$\mathbb{P}_{k}^{d}:=\\{x\mapsto\sigma_{k}(\omega\cdot
x+b)\;:\;\|\omega\|_{2}=1,\;c_{1}\leqslant b\leqslant c_{2}\\}$

where $\sigma_{k}(t):=\max\\{0,t\\}^{k}$ is the so-called RELU-$k$ function.
Here, we work in the space $L^{2}(B)$ where $B$ is an arbitrary ball of
$\mathbb{R}^{d}$, and the constants $(c_{1},c_{2})$ are taken as the inf and
sup of $\omega\cdot x$ as $x\in B$ and $\|\omega\|_{2}=1$, respectively, that
is we take all $b$ such that the line discontinuity of the $k$-th derivative
of $\sigma_{k}(\omega\cdot x+b)$ crosses the ball $B$. Theorem 9 from [38],
which improves on earlier results from [29], shows that if

$B_{1}(\mathbb{P}_{k}^{d}):=\overline{\Big{\\{}\sum_{j=1}^{n}a_{j}g_{j}\;:\;n\in\mathbb{N},\;g_{j}\in\mathbb{P}_{k}^{d},\;\sum_{j=1}^{n}|a_{j}|\leqslant
1\Big{\\}}}$

denotes the symmetrized convex hull of this dictionary (the closure being
taken in $L^{2}(B)$), then

$d_{n}(B_{1}(\mathbb{P}_{k}^{d}))_{L^{2}(B)}\geqslant
cn^{-\frac{2k+1}{2d}},\quad n\geqslant 1,$

where $c$ depends on $k$, $d$, and the diameter of $B$.

In our case of interest we work with the value $d=2$ and $k=0$, so that the
ridge functions are simply the characteristic functions of half-planes. By
convexity, we have

$d_{n}(\mathbb{P}_{0}^{2}))_{L^{2}(B)}=d_{n}(B_{1}(\mathbb{P}_{0}^{2}))_{L^{2}(B)}\geqslant
cn^{-\frac{1}{4}}.$

We take for $B$ the ball of center $(1/2,1/2)$ and radius $1/4$, which is
inside our domain $D=[0,1]^{2}$. It is then readily seen that for $R$ small
enough and $M$ large enough, we can extend any ridge function
$g\in\mathbb{P}_{0}^{2}$ into a characteristic function $\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}$ from ${\cal F}_{s,R,M}$, as illustrated
in Figure 1.

$B$$\Omega$$D$

Figure 1: Example of extension of the indicator of a half-plane on $B$ to the
indicator of a smooth domain $\Omega$ on $D$

Observing that if $E_{D}$ is a linear subspace of $L^{2}(D)$ of dimension at
most $n$, its restriction $E_{B}$ to $B$ is a linear subspace of $L^{2}(B)$ of
dimension at most $n$, and one has

${\rm dist}(\raise
1.29167pt\hbox{\large$\chi$}_{\Omega},E_{B})_{L^{2}(B)}\leqslant{\rm
dist}(\raise 1.29167pt\hbox{\large$\chi$}_{\Omega},E_{D})_{L^{2}(D)}.$

By infimizing, it follows that

$d_{n}({\cal F}_{s,R,M})_{L^{2}(D)}\geqslant
d_{n}(\mathbb{P}_{0}^{2})_{L^{2}(B)}\geqslant cn^{-\frac{1}{4}},$

which concludes the proof. $\Box$

######  Remark 4.5.

The fact that we impose conditions on $R$ and $M$ in the above statement is
natural since the class ${\cal F}_{s,R,M}$ becomes empty if $R$ is not small
enough and $M$ not large enough, due to the fact that the sets $\Omega$ are
assumed to be contained in the interior of $D$.

######  Remark 4.6.

The above results are easily extended to higher dimension $d\geqslant 2$, with
a similar definition for the class ${\cal F}_{s,R,M}$. The rate of
approximation in $L^{q}$ norm by piecewise constant functions on uniform
partitions is then $n^{-\frac{1}{dq}}$, which in the case $q=2$ is proved by a
similar argument to be the best achievable by any linear reconstruction
method. We conjecture that the same holds for more general $1\leqslant
q\leqslant\infty$.

##  5 Shape recovery by nonlinear least-squares

###  5.1 Nonlinear reconstruction on a stencil

We now discuss a nonlinear reconstruction method for $u\in{\cal F}_{s,R,M}$,
whose output $\tilde{u}$ is the indicator of a domain $\tilde{\Omega}$ with
polygonal boundary : on each cell $T$, the domain $\tilde{\Omega}$ coincides
with a certain half plane. In order to define the delimiting line we only use
the average values of $u$ on a $3\times 3$ stencil of cells centered at $T$.

We assume that $h<R$, so that $\Omega$ does not intersect the boundary cells
$T_{i,j}$ with $i$ or $j$ in $\\{1,L\\}$, and fix indices $1<i,j<L$. For the
cell $T=T_{i,j}$, denote $\overline{x}=((i-\frac{1}{2})h,(j-\frac{1}{2})h)$
its center, and

$S=[(i-2)h,(i+1)h]\times[(j-2)h,(j+1)h]=\bigcup_{i-1\leqslant
i^{\prime}\leqslant i+1,\,j-1\leqslant j^{\prime}\leqslant
j+1}T_{i^{\prime}j^{\prime}}$

the stencil composed of $T$ and its $8$ neighboring cells. We define the
nonlinear approximation space

$V_{2}:=\left\\{\raise
1.29167pt\hbox{\large$\chi$}_{\vec{n}\cdot(x-\overline{x})\geqslant
c}\;:\;\>\vec{n}\in\mathbb{S}^{1},c\in\mathbb{R}\right\\},$ (5.14)

which is a two-parameter family as each function is determined by
$(\arg\vec{n},c)\in[0,2\pi)\times\mathbb{R}$, where $\arg\vec{n}$ is the angle
of $\vec{n}$ with respect to the horizontal axis.

Here, our measurements are the average values of $u$ on the cells contained in
$S$

$\ell(u)=(a_{T^{\prime}}(u))_{T^{\prime}\subset S}\in\mathbb{R}^{9}.$

In order to find a reconstruction of $u$ in $V_{2}$ based on these
measurements, we need an inverse stability property of the form (1.7). This is
not possible here, since $\ell$ cancels on all functions $\raise
1.29167pt\hbox{\large$\chi$}_{\Omega}\in V_{2}$ with $\Omega\cap S=\emptyset$.
We therefore restrict the nonlinear family $V_{2}$, and consider only
indicators of half-planes whose boundary passes through the central cell $T$:

$V_{2,T}:=\Big{\\{}\raise 1.29167pt\hbox{\large$\chi$}_{\Omega}\in
V_{2},\partial\Omega\cap T\neq\emptyset\Big{\\}}=\Big{\\{}\raise
1.29167pt\hbox{\large$\chi$}_{\vec{n}\cdot(x-\overline{x})\geqslant
c},\>\vec{n}\in\mathbb{S}^{1},|c|\leqslant\frac{h}{2}|\vec{n}|_{1}\Big{\\}}.$
(5.15)

In this setting, we prove the existence of the following stability constants
for $V=L^{1}(S)$ and $Z=\ell^{1}$, which is the best norm on $\mathbb{R}^{m}$
in view of Theorem 3.2. For notational simplicity, we omit the reference to
$Z$ in these constants.

######  Proposition 5.1.

One has

$\|\ell(u)\|_{1}\leqslant\alpha\|u\|_{L^{1}(S)},\quad u\in L^{1}(D),$ (5.16)

and

$\|u-v\|_{L^{1}(S)}\leqslant\mu\|\ell(u-v)\|_{1},\quad u,v\in V_{2,T},$ (5.17)

where $\alpha=h^{-2}$ and $\mu=C_{0}h^{2}$ with $C_{0}$ an absolute constant.

The proof of the stability property (5.16) is trivial since on each cell

$|a_{T^{\prime}}(u)|\leqslant|T^{\prime}|^{-1}\|u\|_{L^{1}(T^{\prime})}=h^{-2}\|u\|_{L^{1}(T^{\prime})}.$

The proof of the inverse stability (5.17) is quite technical and left to the
appendix. The value of the constant $C_{0}$ obtained in this proof is quite
pessimistic. On the other hand, one can estimate numerically the best value of
this constant by taking a large random collection $(v_{j})_{j=1,\dots,J}$ of
elements in $V_{2,T}$, and computing the empirical value

$\hat{C}_{0}:=\max_{i,j}\frac{\|\ell(v_{i}-v_{j})\|_{1}}{\|v_{i}-v_{j}\|_{L^{1}(S)}}.$

In such numerical computations, we systematically found that
$\hat{C}_{0}\leqslant\frac{3}{2}$, and we observe that the value $\frac{3}{2}$
is attained when taking $v_{1}$ and $v_{2}$ to have value $1$ above and below,
respectively, a horizontal discontinuity line that crosses $T$ at its centers.
We thus conjecture that the exact value of $C_{0}$ is $\frac{3}{2}$.

Given the noisy observation

$z=\ell(u)+\eta\in\mathbb{R}^{9},$

we define the estimator of $u$ on the cell $T$ by

$\tilde{u}_{T}\in\operatorname{argmin}_{v\in V_{2}}\|z-\ell(v)\|_{1}.$ (5.18)

Here we minimize over all $V_{2}$, that is on all indicators of half planes,
but we note that we may restrict to half-planes whose boundary passes through
the stencil $S$.

The following result, which uses Proposition 5.1, shows that its distance to
$u$ in $L^{1}(T)$ is comparable to the error between $u$ and its best
approximation in the $L^{1}(S)$ norm

$\overline{u}_{S}:=\operatorname{argmin}_{v\in V_{2}}\|u-v\|_{L^{1}(S)}.$

######  Lemma 5.2.

For all $u\in{\cal F}_{s,R,M}$, one has

$\|u-\tilde{u}_{T}\|_{L^{1}(T)}\leqslant
C_{1}\|u-\overline{u}_{S}\|_{L^{1}(S)}+2\beta\mu\|\eta\|_{p},$

where $C_{1}=1+2\alpha\mu=1+2C_{0}$ and
$C_{2}=2\beta\mu=2C_{0}9^{1-\frac{1}{p}}h^{2}$, with $\alpha,\mu$ as in
Proposition 5.1, and $\beta=9^{1-\frac{1}{p}}$ the maximal ratio between
$\ell^{p}$ and $\ell^{1}$ norm in $\mathbb{R}^{9}$.

###### Proof.

We distinguish two cases:

  * •

If $\tilde{u}_{T}\in V_{2,T}$ and $\overline{u}_{S}\in V_{2,T}$, that is, both
boundaries pass through the central cell $T$, we apply Theorem 2.9 together
with Proposition 5.1

$\displaystyle\|u-\tilde{u}_{T}\|_{L^{1}(T)}\leqslant\|u-\tilde{u}_{T}\|_{L^{1}(S)}$
$\displaystyle\leqslant C_{1}\min_{v\in
V_{2,T}}\|u-v\|_{L^{1}(S)}+C_{2}\|\eta\|_{p}$
$\displaystyle=C_{1}\|u-\overline{u}_{S}\|_{L^{1}(S)}+C_{2}\|\eta\|_{p}.$

with $C_{1}=1+2\alpha\mu$, $C_{2}=2\beta\mu$.

  * •

Otherwise, either $\tilde{u}_{T}$ or $\overline{u}_{S}$ has constant value $0$
or $1$ on $T$, so $\tilde{u}_{T}-\overline{u}_{S}$ has constant sign on $T$,
and thus

$\displaystyle\|\overline{u}_{S}-\tilde{u}_{T}\|_{L^{1}(T)}$
$\displaystyle=h^{2}|a_{T}(\tilde{u}_{T}-\overline{u}_{S})|\leqslant
h^{2}\|\ell(\tilde{u}_{T}-\overline{u}_{S})\|_{1}\leqslant
h^{2}(\|\ell(\overline{u}_{S})-z\|_{1}+\|\ell(\tilde{u}_{T})-z\|_{1})$
$\displaystyle\leqslant 2h^{2}\|\ell(\overline{u}_{S})-z\|_{1}\leqslant
2h^{2}\|\ell(\overline{u}_{S}-u)\|_{1}+2h^{2}\|\eta\|_{1}\leqslant
2\|u-\overline{u}_{S}\|_{L^{1}(S)}+2h^{2}\beta\|\eta\|_{p}.$

By triangle inequality, it follows that

$\|u-\tilde{u}_{T}\|_{L^{1}(T)}\leqslant
3\|u-\overline{u}_{S}\|_{L^{1}(S)}+2h^{2}\beta\|\eta\|_{p},$

which has better constants than in the estimate obtained in the first case,
since the constant $C_{0}$ is larger than $1$. ∎

The order of the best local approximation error
$\|u-\overline{u}_{S}\|_{L^{1}(S)}$ that appears as a bound for the
reconstruction error $\|u-\tilde{u}_{T}\|_{L^{1}(T)}$ depends on the
smoothness of the boundary, as expressed in the following lemma.

######  Lemma 5.3.

For all $u\in{\cal F}_{s,R,M}$, with $R\geqslant\frac{3}{\sqrt{2}}h$, one has

$\|u-\overline{u}_{S}\|_{L^{1}(S)}\leqslant M(3\sqrt{2}h)^{\min(s,2)+1}.$

###### Proof.

We apply the definition of ${\cal F}_{s,R,M}$ at point $\overline{x}$: as
$R\geqslant\frac{3}{\sqrt{2}}h$, the stencil $S$ is contained in the domain

$\\{\overline{x}+z_{1}e_{1}+z_{2}e_{2},\;|z_{1}|,|z_{2}|\leqslant R\\},$

so $u|_{S}$ is the indicator of a domain delimited by a ${\cal C}^{s}$
function $\psi$, with $\|\psi\|_{{\cal C}^{s}}\leqslant M$. From the
definition of ${\cal C}^{s}$, there exists an affine function $\xi$ such that

$|\psi(z_{1})-\xi(z_{1})|\leqslant
M(3\sqrt{2}h)^{\min(s,2)},\quad|z_{1}|\leqslant\frac{3}{\sqrt{2}}h.$

Then the function $v:\overline{x}+z_{1}e_{1}+z_{2}e_{2}\mapsto\raise
1.29167pt\hbox{\large$\chi$}_{z_{2}\leqslant\xi(z_{1})}$ belongs to $V_{2}$,
and we have

$\|u-\overline{u}_{S}\|_{L^{1}(S)}\leqslant\|u-v\|_{L^{1}(S)}\leqslant
M(3\sqrt{2}h)^{\min(s,2)+1}.$

∎

###  5.2 Global nonlinear reconstruction

We now consider the process of recovering $u\in{\cal F}_{s,R,M}$ globally from
its data

$z=\ell(u)+\eta,$

where now $\ell(u):=(a_{T}(u))_{T\in{\cal T}}\in\mathbb{R}^{n}$ and
$\eta\in\mathbb{R}^{n}$ is the noise vector. Applying to each inner cell
$T\in{\cal T}$ the previous reconstruction procedure based on the $3\times 3$
stencil $S$ centered at $T$, we obtain a global recovery
$\tilde{u}=\tilde{u}(z)$ such that

$\tilde{u}|_{T}=\tilde{u}_{T}|_{T},\quad T=T_{i,j}\in{\cal T},\quad 1<i,j<L,$

where $\tilde{u}_{T}$ is the local estimator from (5.18). On the boundary
cells $T=T_{i,j}$ with $i$ or $j$ in $\\{1,L\\}$, $u|_{T}$ is zero by
Definition 4.1 so we simply set $\tilde{u}|_{T}=0$. Note that $\tilde{u}$ is
of the form

$\tilde{u}=\raise 1.29167pt\hbox{\large$\chi$}_{\tilde{\Omega}},$

where $\tilde{\Omega}$ has piecewise linear boundary with respect to the mesh
${\cal T}$. The following result gives a global approximation bound, which
confirms the improvement over linear methods when $s>1$.

######  Theorem 5.4.

For all $u\in{\cal F}_{s,R,M}$, one has

$\|u-\tilde{u}\|_{L^{q}(D)}\leqslant
C_{1}n^{-\frac{\min(1,s/2)}{q}}+C_{2}n^{-\frac{1}{pq}}\|\eta\|_{p}^{\frac{1}{q}}.$

###### Proof.

First notice that if the result is proved for $p=q=1$, as $u-v$ has values in
$\\{-1,0,1\\}$,

$\|u-v\|_{L^{q}(D)}^{q}=\|u-v\|_{L^{1}(D)}\leqslant
C_{1}n^{-1}+C_{2}n^{-1}\|\eta\|_{1}\leqslant\left(C_{1}^{\frac{1}{q}}n^{-\frac{1}{q}}+C_{2}^{\frac{1}{q}}n^{-\frac{1}{pq}}\|\eta\|_{p}^{\frac{1}{q}}\right)^{q},$

so it suffices treat the case $p=q=1$.

By an argument similar to the proof of Proposition 4.3, $\partial\Omega$
intersects at most $16N^{2}\lceil 2R\sqrt{1+M^{2}}/h\rceil$ stencils of $9$
cells. Using the fact that $u=\overline{u}_{S}$ is a constant on any other
stencil, we get

$\displaystyle\|u-\tilde{u}\|_{L^{1}(D)}$ $\displaystyle=\sum_{T\text{ inner
cell}}\|u-\tilde{u}\|_{L^{1}(T)}\leqslant\sum_{T\text{ inner
cell}}(1+2\alpha\mu)\|u-\overline{u}\|_{L^{1}(S)}+2\beta\mu\|\eta\|_{\ell^{1}(S)}$
$\displaystyle\leqslant
16N^{2}\left\lceil\frac{2R\sqrt{1+M^{2}}}{h}\right\rceil
M(3\sqrt{2}h)^{\min(s,2)+1}+18\beta\mu\|\eta\|_{1}\leqslant
C_{1}h^{\min(s,2)}+C_{2}h^{2}\|\eta\|_{1}.$

We conclude by recalling that $n=h^{-2}$. ∎

######  Remark 5.5.

Here the convergence rate for the noiseless term $n^{-\frac{\min(1,s/2)}{q}}$
is limited due to the use of polygonal domains in the reconstruction. So the
best approximation rate $h^{\frac{2}{q}}=n^{-\frac{1}{q}}$ is already attained
for ${\cal C}^{2}$ boundaries. When the smoothness parameter $s$ is larger
than $2$, better rates $n^{-\frac{s}{2q}}$ should be reachable if we use non-
linear approximation spaces that are richer than the space $V_{2}$, for
example indicator functions of domains with boundary that have a higher order
polynomial description rather than straight lines. Of course, the stable
identification of these approximants in the sense of (1.7) might require
stencils that are of larger size than $3\times 3$.

######  Remark 5.6.

If $\|\eta\|_{\infty}\leqslant\frac{1}{9}$, then $\tilde{u}$ is exactly equal
to $u$ on any cell whose corresponding stencil does not intersect
$\partial\Omega$, so the error is concentrated on $\mathcal{O}(\sqrt{n})$
cells, leading to an improved rate $n^{-\frac{p+1}{2pq}}$ instead of
$n^{-\frac{1}{pq}}$ for the noise term.

###  5.3 Numerical illustration

We study the behavior of the above discussed linear and non-linear recovery
methods from cell averages for the particular target function
$u=\chi_{\Omega}$, with $\Omega$a slightly decentered disk of radius
$r=0.325$.

The linear method consists of the piecewise constant approximation (4.13),
referred to as PiecewiseConstant. As to the nonlinear method, for the local
best fit problem, we use the $\ell^{2}$ norm on $\mathbb{R}^{9}$ instead of
the $\ell^{1}$ norm. By norm equivalence on $\mathbb{R}^{9}$, the same
convergence results can be proved to hold with different constants. This
method, which we refer to as LinearInterface, does not ensure consistency of
the reconstruction in the sense that $a_{T}(\tilde{u})=a_{T}(u)$. One way to
approach this consistency property is to modify the $\ell^{2}$ norm by putting
a large weight on the central cell. We refer to this variant as
LinearInterfaceCC, here taking the weight $100$.

Figure 2: Convergence curves for the linear and nonlinear recovery methods

(a)

(b)

(c)

(d)

Figure 3: (a) The target function, (b) its recovery by PiecewiseConstant
showing the cell-average data, and the recovered boundaries by (c)
LinearInterface and (d) LinearInterfaceCC methods.

Figure 2 shows the convergence rates of the three methods in the $L^{1}$ norm.
The expected $h^{-2}$ decay is observed in both non-linear methods while the
linear method lays behind with a decay rate of $h^{-1}$. It is relevant to
note that although both non-linear methods benefit from the same rate, the
associated constants differ by an order of magnitude, showing the practical
improvement gained by imposing consistency. This improvement is also visible
on Figure 3 which shows that in the LinearInterface method, the interfaces
that minimize the $l_{2}$ error on the $9$ surrounding cells lay always inside
the circle as the curvature of the boundary pushes them towards the center. On
the contrary, LinearInterfaceCC seems to find the right compromise between
sticking to the cell average while capturing at the same time the curvature
trend hinted by the surrounding cell averages.

##  6 Relation to compressed sensing

###  6.1 Compressed sensing and best $n$-term approximation

In this section we discuss the application of our setting to the sparse
recovery of large vectors from a few linear observations. We thus take

$V=\mathbb{R}^{N},$

equipped with some given norm $\|\cdot\|_{V}$ of interest. The linear
measurements of $u=(u_{1},\dots,u_{N})^{\top}\in\mathbb{R}^{N}$ are given by

$(\ell_{1}(u),\dots,\ell_{m}(u))^{\top}=\Phi u,$

where $\Phi$ is an $m\times N$ measurement matrix, with typically $m\ll N$.

The topic of compressed sensing deals with sparse recovery of $u$ from such
measurements, that is, searching to recover an accurate approximation to $u$
by a vector with only a few non-zero components. We refer to [10] for some
first highly celebrated breakthrough results and to [21] for a general
treatment.

We define the nonlinear space of $n$-sparse vectors as

$V_{n}:=\Big{\\{}u\in\mathbb{R}^{N}\;:\;\|u\|_{0}:=\\#\\{i\;:\;u_{i}\neq
0\\}\leqslant n\Big{\\}},$

and the best $n$-term approximation error in the $V$ norm as

$e_{n}(u)_{V}:=\min_{v\in V_{n}}\|u-v\|_{V}.$

One natural question is to understand for which type of measurement matrices
$\Phi$ does the noise-free measurement $y=\Phi u$ contain enough information,
in order to recover any $u$ up to an error $e_{n}(u)_{V}$. In other words, one
asks if there exists a recovery map $R:\mathbb{R}^{m}\to\mathbb{R}^{N}$ such
that one has the instance optimality property at order $n$

$\|u-R(\Phi u)\|_{V}\leqslant C_{0}e_{n}(u)_{V},\quad u\in\mathbb{R}^{N},$
(6.19)

with $C_{0}$ a fixed constant, which we denote by $IOP(n,C_{0})$. This
question has been answered in [12] in terms of the null space ${\cal
N}:=\\{v\in \mathbb{R}^{N}\;:\;\Phi v=0\\}$. We say that $\Phi$ satisfies the
null space property at order $k$ with constant $C_{1}$, denoted by
$NSP(k,C_{1})$ if and only if

$\|v\|_{V}\leqslant C_{1}e_{k}(v)_{V},\quad v\in{\cal N}.$ (6.20)

This property quantifies how much vectors from the null space can be
concentrated on a few coordinates. One main result of [12] is the equivalence
between $IOP$ at order $n$ and $NSP$ at order $2n$ in the following sense.

######  Theorem 6.1.

One has $IOP(n,C_{0})\Rightarrow NSP(2n,C_{0})$ and conversely
$NSP(2n,C_{1})\Rightarrow IOP(n,2C_{1})$.

One natural question is whether matrices $\Phi$ with such properties can be
constructed with a number of rows/measurements $m$ barely larger than $n$. As
we recall further the answer to this question is strongly tied to the norm $V$
used on $\mathbb{R}^{N}$.

###  6.2 Stability and the null space property

The nonlinear estimation results that we have obtained in §2 and §3 can be
applied to the setting of sparse recovery, offering us a different vehicle
than the null space property to establish instance optimality.

In the present setting, for a given norm $\|\cdot\|_{Z}$, the stability
property (1.6) takes the form

$\|\Phi u\|_{Z}\leqslant\alpha_{Z}\|u\|_{V},\quad u\in \mathbb{R}^{N}$ (6.21)

and the inverse stability property (1.7) takes the form

$\|v\|_{V}\leqslant\mu_{Z}\|\Phi v\|_{Z},\quad v\in V_{2n},$ (6.22)

since for sparse vectors we have $V_{n}^{\rm diff}=V_{n}-V_{n}=V_{2n}$. We
refer to these properties as $S(\alpha_{Z})$ and $IS(2n,\mu_{Z})$,
respectively.

Application of Theorem 2.3 in the noiseless case immediately gives us that the
nonlinear best fit recovery $R(\Phi u)=\tilde{u}$ satisfies the instance
optimality bound (6.19) with constant $C_{0}=1+2\alpha_{Z}\mu_{Z}$. In other
words

$S(\alpha_{Z})\;{\rm and}\;IS(2n,\mu_{Z})\Rightarrow IOP(n,C_{0}),\quad\quad
C_{0}=1+2\alpha_{Z}\mu_{Z}.$ (6.23)

The following result shows that $(S,IS)$ is actually equivalent to $NSP$, and
thus to $IOS$, in the sense that a converse result holds when $\|\cdot\|_{Z}$
is chosen to be the Riesz norm (3.10).

######  Theorem 6.2.

For any norm $\|\cdot\|_{Z}$, one has

$S(\alpha_{Z})\;{\rm and}\;IS(2n,\mu_{Z})\Rightarrow NSP(2n,C_{1}),\quad\quad
C_{1}=1+\alpha_{Z}\mu_{Z}.$ (6.24)

Conversely, let $\|\cdot\|_{W}$ be the Riesz norm so that $\|\Phi
u\|_{W}=\min_{\Phi v=\Phi u}\|v\|_{V}$, then

$NSP(2n,C_{1})\Rightarrow S(\alpha_{W})\;{\rm
and}\;IS(2n,\mu_{W}),\quad\quad\alpha_{W}=1\;{\rm and}\;\mu_{W}=1+C_{1}.$
(6.25)

###### Proof.

Assume that $S(\alpha_{Z})$ and $IS(2n,\mu_{Z})$ hold. Let $v\in{\cal N}$ and
$\tilde{v}$ its best approximation in $V_{2n}$, then

$\displaystyle\|v\|_{V}$
$\displaystyle\leqslant\|v-\tilde{v}\|_{V}+\|\tilde{v}\|_{V}$
$\displaystyle\leqslant e_{2n}(v)_{V}+\mu_{Z}\|\Phi\tilde{v}\|_{W}$
$\displaystyle=e_{2n}(v)_{V}+\mu_{Z}\|\Phi(v-\tilde{v})\|_{W}\leqslant(1+\alpha_{Z}\mu_{Z})e_{2n}(x)_{V}.$

This shows that $NSP(2n,C_{1})$ holds with $C_{1}=1+\alpha_{Z}\mu_{Z}$.

Conversely, assume that $NSP(2n,C_{1})$ holds. From the definition of the
Riesz norm, it is immediate that $S(\alpha_{W})$ holds with $\alpha_{W}=1$.
For $v\in V_{2n}$, let $\tilde{v}$ be the minimizer of
$\min_{\Phi\tilde{v}=\Phi v}\|\tilde{v}\|_{V}$. Then, one has

$\|v\|_{V}\leqslant\|\tilde{v}\|_{V}+\|v-\tilde{v}\|_{V}\leqslant\|\tilde{v}\|_{V}+C_{1}\sigma_{2n}(v-\tilde{v})_{V}\leqslant(1+C_{1})\|\tilde{v}\|_{V},$

by using $v$ as a sparse approximation to $v-\tilde{v}$. Since
$\|\tilde{v}\|_{V}=\|\Phi v\|_{W}$ , this shows that $S(2n,\mu_{W})$ holds
with $\mu_{W}=1+C_{1}$. ∎

###  6.3 The case of $\ell^{p}$ norms

The range of $m$ allowing the properties to be fulfilled is best understood in
the case of the $\ell^{p}$ norms, that is $\|\cdot\|_{V}=\|\cdot\|_{p}$, as
discussed in [12] which points out a striking difference between the case
$p=2$ and $p=1$:

  1. 1.

In the case $p=2$, it is proved that $NSP(2,C_{1})$ cannot hold unless
$N\leqslant C_{1}^{2}m$. In other words, instance optimality in $\ell^{2}$
even at order $n=1$ requires a number of measurements that is proportional to
the full space dimension.

  2. 2.

In the more favorable case $p=1$, it is proved that for matrices which satisfy
the $\ell^{2}$-RIP property of order $3n$

$(1-\delta)\|v\|_{2}^{2}\leqslant\|\Phi
v\|_{2}^{2}\leqslant(1+\delta)\|v\|^{2}_{2},\quad v\in V_{3n},$

with parameter $0<\delta<\frac{(\sqrt{2}-1)^{2}}{3}$, the $NSP(2n,C_{1})$
holds with $C_{1}$ depending on $\delta$. Such matrices are known to exists
with $m\sim n\log(N/n)$ rows.

Our setting based on the stability properties $S$ and $IS$ applies more
naturally to a different class of matrices built from graphs, which is also
known to be well adapted for sparse recovery in the $\ell^{1}$ norm. A
bipartite graph with $(N,m)$ left and right vertices, and of left degree $d$,
is an $(l,\varepsilon)$-graph expander if

$|X|\leqslant l\Rightarrow|N(X)|\geqslant d(1-\varepsilon)|X|,\quad
X\subset\\{1,\dots,N\\},$

where $N(X)\subset\\{1,\dots,m\\}$ is the set of vertices connected to $X$. We
necessarily have $|N(X)|\leqslant d|X|$, and $(1-\varepsilon)dl\geqslant m$.
From [11], it is known that there exists a $(2n,\frac{1}{2})$-graph expander
with $d\sim\log\frac{N}{n}$ and $m\sim nd\sim n\log(N/n)$.

Now denote $\Phi\in\\{0,1\\}^{m\times N}$ the adjacency matrix of this graph,
so that each column of $\Phi$ has $d$ nonzero entries. Then

$\|\Phi x\|_{1}\leqslant d\|x\|_{1},\quad x\in\mathbb{R}^{N},$

and

$\|\Phi x\|_{1}\geqslant d(1-\varepsilon)\|x\|_{1},\quad x\in V_{2n}.$

Therefore $S(\alpha_{1})$ and $IS(2n,\mu_{1})$, hold with $\alpha_{1}=d$ and
$\mu_{1}=\frac{1}{d(1-\varepsilon)}=\frac{2}{d}$, which by (6.24) and (6.23)
gives $NSP(2n,C_{1})$ with $C_{1}=3$ and $IOP(n,C_{0})$ with $C_{0}=5$.

##  7 Appendix: Proof of Proposition 5.1

$T_{i-1,j-1}$$T_{i-1,j}$$T_{i-1,j+1}$$T_{i,j-1}$$T=T_{i,j}$$T_{i,j+1}$$T_{i+1,j-1}$$T_{i+1,j}$$T_{i+1,j+1}$

$c_{u}$${\vec{n}_{u}}$$\overline{x}$

Figure 4: Partition of the stencil $S$ and representation of a function $u\in
V_{2}$ with the shaded region corresponding to $u=1$, and the dotted line to
$\Delta_{u}$

Take $w=u-v\in V_{2,T}^{\rm diff}$, with $u,v\in V_{2,T}$, and denote
$\vec{n}_{u},\vec{n}_{v}$ and $c_{u},c_{v}$ the corresponding unit vectors and
offsets from the definition (5.14) of $V_{2,T}$. Define also
$\Delta_{u}=\\{x,\vec{n}_{u}\cdot(x-\overline{x})=c_{u}\\}$ as the delimiting
line between the regions $\\{u=0\\}$ and $\\{u=1\\}$, and similarly for
$\Delta_{v}$. We distinguish two cases.

Case 1: $\vec{n}_{u}\cdot\vec{n}_{v}\leqslant 0$. Up to exchanging $u$ and
$v$, we can assume that $c_{u}\leqslant c_{v}$, i.e that $\overline{x}$ is in
the half-plane

$P=\\{x\in\mathbb{R}^{2},\,(\vec{n}_{u}-\vec{n}_{v})\cdot(x-\overline{x})\geqslant
c_{u}-c_{v}\\}$

of points closer to the cone $\\{w=1\\}$ than to $\\{w=-1\\}$.

${\vec{n}_{u}}$${\vec{n}_{v}}$$\overline{x}$$P$$P^{c}$$B$

Figure 5: Top blue region: $w=1$, bottom red region: $w=-1$. The dotted line
is $\partial P$, the dashed ones are $\Delta_{u}$ and $\Delta_{v}$, and $y$ is
picked in the tiny ball $B$

Let $\varepsilon=\left(\frac{1}{2}-\frac{\sqrt{2}}{3}\right)h>0$, and take
$y\in
B:=B\left(\overline{x}+(\sqrt{2}h+2\varepsilon)\frac{\vec{n}_{u}-\vec{n}_{v}}{|\vec{n}_{u}-\vec{n}_{v}|},\varepsilon\right)$.
Then

$|y-\overline{x}|<\sqrt{2}h+3\varepsilon\leqslant\frac{3}{2}h,$

so $y\in S$. Moreover

$d(y,\\{w=-1\\})\geqslant
d(y,P^{c})=(z-\overline{x})\cdot\frac{\vec{n}_{u}-\vec{n}_{v}}{|\vec{n}_{u}-\vec{n}_{v}|}+c_{v}-c_{u}>\sqrt{2}h+\varepsilon>\sqrt{2}h,$

so $w\geqslant 0$ on the cell containing $y$. Finally, as
$\vec{n}_{u}\cdot\vec{n}_{v}\leqslant 0$,

$\vec{n}_{u}\cdot(y-\overline{x})>(\sqrt{2}h+2\varepsilon)\frac{1-\vec{n}_{v}\cdot\vec{n}_{u}}{|\vec{n}_{u}-\vec{n}_{v}|}-\varepsilon\geqslant\frac{\sqrt{2}h+2\varepsilon}{2}-\varepsilon=\frac{h}{\sqrt{2}}\geqslant
c_{u}$

and

$\vec{n}_{v}\cdot(y-\overline{x})<(\sqrt{2}h+2\varepsilon)\frac{\vec{n}_{u}\cdot\vec{n}_{v}-1}{|\vec{n}_{u}-\vec{n}_{v}|}+\varepsilon\leqslant-\frac{\sqrt{2}h+2\varepsilon}{2}+\varepsilon=-\frac{h}{\sqrt{2}}\leqslant
c_{v},$

where the last inequalities come from the bound
$|c|\leqslant\frac{h}{2}|n|_{1}$ in (5.14). Therefore $u(y)=1$ and $v(y)=0$,
so $w(y)=1$. In the end, $w$ is non-negative on any cell intersecting $B$, and
$w=1$ on $B$, so

$\|\ell(w)\|_{1}\geqslant\sum_{T\subset S,\,T\cap
B\neq\emptyset}\frac{1}{|T|}\int_{T}w\geqslant h^{-2}\int_{B}w=h^{-2}|B|>0.$

As $\|w\|_{L^{1}(S)}\leqslant|S|=9h^{2}$, bound (5.17) is proved for
$\vec{n}_{u}\cdot\vec{n}_{v}\leqslant 0$ and $r=1$, with

$\mu\leqslant\frac{|S|}{h^{-2}|B|}=\frac{9}{\pi}\left(\frac{1}{2}-\frac{\sqrt{2}}{3}\right)^{-2}h^{2}.$

Case 2. Now assume that $\vec{n}_{u}\cdot\vec{n}_{v}>0$. If the borders
$\Delta_{u}$ and $\Delta_{v}$ do not intersect inside $S$ then $w$ is of
constant sign on $S$, so $\|\ell(w)\|_{1}=h^{-2}\|w\|_{L^{1}(S)}$. Otherwise,
denote $I$ their point of intersection, and define $\theta_{u},\theta_{v}$
such that

$\vec{n}_{u}=\begin{pmatrix}\cos\theta_{u}\\\
\sin\theta_{u}\end{pmatrix}\text{ and
}\vec{n}_{v}=\begin{pmatrix}\cos\theta_{v}\\\ \sin\theta_{v}\end{pmatrix}.$

$I$$x_{u}$$x_{v}$

$I$$x_{u}$$x_{v}$

Figure 6: Examples of the two possible cases when
$\vec{n}_{u}\cdot\vec{n}_{v}>0$. On the left $x_{u}$ and $x_{v}$ are on the
same side of $I$, so all the lines $\Delta_{\theta}$ of the shaded cone pass
through $T$. On the right, $x_{u}$ and $x_{v}$ are on opposite sides with
respect to $I$, but $I$ is inside the shape delimited by four semicircles.
$\overline{S}$ is the band between the two dotted lines

Up to exchanging $u$ and $v$, we can assume that
$\theta_{v}\in[\theta_{u},\theta_{u}+\frac{\pi}{2}[$. We next consider the
lines

$\Delta_{\theta}=\\{I+t\vec{n}_{\theta}^{\perp},t\in\mathbb{R}\\},\quad\vec{n}_{\theta}^{\perp}=\begin{pmatrix}-\sin\theta\\\
\cos\theta\end{pmatrix},$

for $\theta\in[\theta_{u},\theta_{v}]$, and show that they intersect $S$ on a
line segment $S_{\theta}:=\Delta_{\theta}\cap S$ of length at least
$\sqrt{6}h$. To do so, we take two points
$x_{u}=I+t_{u}\vec{n}_{\theta_{u}}^{\perp}\in\Delta_{\theta_{u}}\cap T$ and
$x_{v}=I+t_{v}\vec{n}_{\theta_{v}}^{\perp}\in\Delta_{\theta_{v}}\cap T$, which
is possible since $\Delta_{\theta_{u}}=\Delta_{u}$ and
$\Delta_{\theta_{v}}=\Delta_{v}$ intersect $T$, and consider two cases:

  * •

If $t_{u}\,t_{v}\geqslant 0$, then the line segment $[x_{u},x_{v}]$ is
contained in $T$ and intersects all the $\Delta_{\theta}$ for
$\theta\in[\theta_{u},\theta_{v}]$. Thus all lines $\Delta_{\theta}$ pass
through $T$, which implies that $|S_{\theta}|\geqslant
2\sqrt{2}h\geqslant\sqrt{6}h$.

  * •

Otherwise, the triangle $(I,x_{u},x_{v})$ is obtuse at point $I$, which
implies that $I$ is inside the disc of diameter $[x_{u},x_{v}]$. As $x_{u}$
and $x_{v}$ belong to $T$, point $I$ is in one of the four discs
$B(y,\frac{h}{2})$, $y=\overline{x}+(\pm\frac{h}{2},\pm\frac{h}{2})$. A
minimization over all possible values of $\theta$ ensures that

$|S_{\theta}|\geqslant\min\left(2\sqrt{2},\frac{3-|\\!\cos\theta|}{2\,|\\!\sin\theta|}+\frac{2-|\\!\sin\theta|}{2\,|\\!\cos\theta|}\right)h\geqslant\sqrt{6}h.$

Writing
$S_{\theta}=[I+s_{\theta}\vec{n}_{\theta}^{\perp},I+t_{\theta}\vec{n}_{\theta}^{\perp}]$
with $s_{\theta}\leqslant 0\leqslant t_{\theta}$, we now compute

$\displaystyle\int_{S}|w|$
$\displaystyle=\int_{\theta_{u}}^{\theta_{v}}\int_{s_{\theta}}^{t_{\theta}}|t|dtd\theta=\int_{\theta_{u}}^{\theta_{v}}\frac{t_{\theta}^{2}+s_{\theta}^{2}}{2}d\theta$
$\displaystyle\geqslant\int_{\theta_{u}}^{\theta_{v}}\left(\frac{t_{\theta}-s_{\theta}}{2}\right)^{2}d\theta=\int_{\theta_{u}}^{\theta_{v}}\frac{|S_{\theta}|^{2}}{4}d\theta\geqslant\frac{3}{2}(\theta_{v}-\theta_{u})h^{2}.$

Similarly, we compute the integral of $|w|$ on

$\overline{S}=\left\\{x\in
S,|(\vec{n}_{u}-\vec{n}_{v})\cdot(x-\overline{x})-c_{u}+c_{v}|\leqslant\frac{h}{\sqrt{2}}\right\\}$

the band of width $\sqrt{2}h$ around $\partial P$. As
$\Delta_{\theta}\cap\overline{S}\subset[I-(\sqrt{2}\cos(\theta-\overline{\theta}))^{-1}h\vec{n}_{\theta}^{\perp},I+(\sqrt{2}\cos(\theta-\overline{\theta}))^{-1}h\vec{n}_{\theta}^{\perp}]$,
where $\overline{\theta}=\frac{\theta_{u}+\theta_{v}}{2}$, we have

$\displaystyle h^{-2}\int_{\overline{S}}|w|$
$\displaystyle\leqslant\int_{\theta_{u}}^{\theta_{v}}\int_{-(\sqrt{2}\cos(\theta-\overline{\theta}))^{-1}}^{(\sqrt{2}\cos(\theta-\overline{\theta}))^{-1}}|t|dtd\theta=\int_{\theta_{u}}^{\theta_{v}}\frac{d\theta}{2\cos^{2}(\theta-\overline{\theta})}$
$\displaystyle=\frac{1}{2}\tan(\theta_{v}-\overline{\theta})-\frac{1}{2}\tan(\theta_{u}-\overline{\theta})=\tan\left(\frac{\theta_{v}-\theta_{u}}{2}\right)\leqslant\frac{2}{\pi}(\theta_{v}-\theta_{u})$

since $0\leqslant\theta_{v}-\theta_{u}\leqslant\frac{\pi}{2}$. In the end, as
two points $x,y\in S$ such that $w(x)=1$ and $w(y)=-1$ can only be in the same
sub-square if one of them is in $\overline{S}$,

$\|\ell(w)\|_{1}\geqslant
h^{-2}\left(\int_{S}|w|-2\int_{\overline{S}}|w|\right)\geqslant\left(\frac{3}{2}-\frac{4}{\pi}\right)(\theta_{v}-\theta_{u}),$

whereas

$\|w\|_{L^{1}(S)}\leqslant\int_{\theta_{u}}^{\theta_{v}}\frac{t_{\theta}^{2}+s_{\theta}^{2}}{2}d\theta\leqslant
18(\theta_{v}-\theta_{u})h^{2}$

since $|s_{\theta}|,|t_{\theta}|\leqslant 3\sqrt{2}h$. Therefore, the
Proposition 5.1 is proved for $\vec{n}_{u}\cdot\vec{n}_{v}>0$ and $r=1$, with

$\mu\leqslant 18\left(\frac{3}{2}-\frac{4}{\pi}\right)^{-1}h^{2}.$

## References

  * [1] B. Adcock, A. C. Hansen, and C. Poon.  Beyond consistent reconstructions: optimality and sharp bounds for generalized sampling, and application to the uniform resampling problem.  SIAM Journal on Mathematical Analysis, 45(5):3132–3167, 2013.
  * [2] F. Arandiga, A. Cohen, R. Donat, and N. Dyn.  Interpolation and approximation of piecewise smooth functions.  SIAM Journal on Numerical Analysis, 43(1):41–57, 2005.
  * [3] J.-P. Argaud, B. Bouriquet, F. de Caso, H. Gong, Y. Maday, and O. Mula.  Sensor placement in nuclear reactors based on the generalized empirical interpolation method.  Journal of Computational Physics, 363:354 – 370, 2018.
  * [4] B. Battisti, T. Blickhan, G. Enchery, V. Ehrlacher, D. Lombardi, and O. Mula.  Wasserstein model reduction approach for parametrized flow problems in porous media.  working paper or preprint, May 2022.
  * [5] P. Benner, A. Cohen, M. Ohlberger, and K. Willcox.  Model Reduction and Approximation: Theory and Algorithms, volume 15.  SIAM, 2017.
  * [6] P. Berger, K. Gröchenig, and G. Matz.  Sampling and reconstruction in distinct subspaces using oblique projections.  arXiv preprint arXiv:1706.06444, 2017.
  * [7] P. Binev, A. Cohen, W. Dahmen, R. DeVore, G. Petrova, and P. Wojtaszczyk.  Data assimilation in reduced modeling.  SIAM/ASA Journal on Uncertainty Quantification, 5(1):1–29, 2017\.
  * [8] P. Binev, A. Cohen, O. Mula, and J. Nichols.  Greedy algorithms for optimal measurements selection in state estimation using reduced models.  SIAM/ASA Journal on Uncertainty Quantification, 6(3):1101–1126, 2018.
  * [9] B. Bojanov.  Optimal recovery of functions and integrals.  In First European Congress of Mathematics, pages 371–390. Springer, 1994.
  * [10] Emmanuel J. Candès, Justin K. Romberg, and Terence Tao.  Stable signal recovery from incomplete and inaccurate measurements.  Commun. Pure Appl. Math., 59(8):1207–1223, 2006.
  * [11] Michael Capalbo, Omer Reingold, Salil Vadhan, and Avi Wigderson.  Randomness conductors and constant-degree lossless expanders.  In Proceedings of the thiry-fourth annual ACM symposium on Theory of computing, pages 659–668, 2002.
  * [12] A. Cohen, W. Dahmen, and R. DeVore.  Compressed sensing and best k-term approximation.  Journal of the American mathematical society, 22(1):211–231, 2009\.
  * [13] A. Cohen, W. Dahmen, R. DeVore, J. Fadili, O. Mula, and J. Nichols.  Optimal reduced model algorithms for data-based state estimation.  SIAM Journal on Numerical Analysis, 58(6):3355–3381, 2020.
  * [14] A. Cohen, W. Dahmen, O. Mula, and J. Nichols.  Nonlinear reduced models for state and parameter estimation.  SIAM/ASA Journal on Uncertainty Quantification, 10(1):227–267, 2022\.
  * [15] A. Cohen and R. DeVore.  Approximation of high-dimensional parametric pdes.  Acta Numerica, 24:1–159, 2015.
  * [16] A. Cohen, R. DeVore, and C. Schwab.  Analytic regularity and polynomial approximation of parametric and stochastic elliptic PDE’s.  Analysis and Applications, 09(01):11–47, 2011.
  * [17] R. DeVore, G. Petrova, and P. Wojtaszczyk.  Greedy algorithms for reduced bases in Banach spaces.  Constructive Approximation, 37(3):455–466, 2013.
  * [18] R.A. DeVore.  Nonlinear approximation.  Acta numerica, 7:51–150, 1998.
  * [19] J. L. Eftang, A. T. Patera, and E. M. Rönquist.  An ”hp” certified reduced basis method for parametrized elliptic partial differential equations.  SIAM Journal on Scientific Computing, 32(6):3170–3200, 2010.
  * [20] V. Ehrlacher, D. Lombardi, O. Mula, and F.-X. Vialard.  Nonlinear model reduction on metric spaces. application to one-dimensional conservative pdes in wasserstein spaces.  ESAIM M2AN, 54(6):2159–2197, 2020.
  * [21] Simon Foucart and Holger Rauhut.  An invitation to compressive sensing.  In A mathematical introduction to compressive sensing, pages 1–39. Springer, 2013.
  * [22] F. Galarce, D. Lombardi, and O. Mula.  State estimation with model reduction and shape variability. application to biomedical problems.  SIAM Journal on Scientific Computing, 44(3):B805–B833, 2022.
  * [23] C. Greif and K. Urban.  Decay of the kolmogorov n-width for wave problems.  Applied Mathematics Letters, 96:216–222, 2019.
  * [24] J. K. Hammond, R. Chakir, F. Bourquin, and Y. Maday.  PBDW: A non-intrusive reduced basis data assimilation method and its application to an urban dispersion modeling framework.  Applied Mathematical Modelling, 76:1–25, 2019.
  * [25] J. S. Hesthaven, G. Rozza, and B. Stamm.  Certified reduced basis methods for parametrized partial differential equations.  SpringerBriefs in Mathematics, 2015.
  * [26] T. Hrycak and K. Gröchenig.  Pseudospectral fourier reconstruction with the modified inverse polynomial reconstruction method.  Journal of Computational Physics, 229(3):933–946, 2010.
  * [27] Y. Maday, A. T. Patera, J. D. Penn, and M. Yano.  A parameterized-background data-weak approach to variational data assimilation: formulation, analysis, and application to acoustics.  International Journal for Numerical Methods in Engineering, 102(5):933–965, 2015.
  * [28] Y. Maday and B. Stamm.  Locally adaptive greedy approximations for anisotropic parameter reduced basis spaces.  SIAM J. Scientific Computing, 35(6):A2417–A2441, 2013.
  * [29] Yuly Makovoz.  Random approximants and neural networks.  Journal of Approximation Theory, 85(1):98–109, 1996.
  * [30] C. A. Micchelli and Th. J. Rivlin.  A survey of optimal recovery.  Springer, 1977.
  * [31] O. Mula.  Inverse problems: A deterministic approach using physics-based reduced models.  Submitted, 2022.
  * [32] E. Novak and H. Wozniakowski.  Tractability of multivariate problems, volume i: Linear information, european math.  Soc., Zürich, 2(3), 2008.
  * [33] M. Ohlberger and S. Rave.  Reduced basis methods: Success, limitations and future challenges.  In Proceedings of the Conference Algoritmy, pages 1–12, 2016.
  * [34] James Edward Pilliod.  An analysis of piecewise linear interface reconstruction algorithms for volume-of-fluid methods.  U. of Calif., Davis, 1992.
  * [35] James Edward Pilliod Jr and Elbridge Gerry Puckett.  Second-order accurate volume-of-fluid algorithms for tracking material interfaces.  Journal of Computational Physics, 199(2):465–502, 2004.
  * [36] Elbridge Gerry Puckett.  A volume-of-fluid interface tracking algorithm with applications to computing shock wave refraction.  In proceedings of the fourth international symposium on Computational Fluid Dynamics, pages 933–938, 1991.
  * [37] G. Rozza, D. B. P. Huynh, and A. T. Patera.  Reduced basis approximation and a posteriori error estimation for affinely parametrized elliptic coercive partial differential equations.  Archives of Computational Methods in Engineering, 15(3):1, 9 2007\.
  * [38] Jonathan W Siegel and Jinchao Xu.  Sharp bounds on the approximation rates, metric entropy, and $n$-widths of shallow neural networks.  arXiv preprint arXiv:2101.12365, 2021.
  * [39] H. Tran, C. G. Webster, and G. Zhang.  Analysis of quasi-optimal polynomial approximations for parameterized pdes with deterministic and stochastic coefficients.  Numerische Mathematik, 137(2):451–493, 2017.
  * [40] G. Welper.  Transformed snapshot interpolation.  arXiv preprint arXiv:1505.01227, 2015.
  * [41] Z. Zou, D. Kouri, and W. Aquino.  An adaptive local reduced basis method for solving pdes with uncertain inputs and evaluating risk.  Computer Methods in Applied Mechanics and Engineering, 345:302–322, 2019.