<figure><img src="image_1.png" /></figure>

HAL
open science

# Optimal pointwise sampling for $L^2$ approximation

Albert Cohen, Matthieu Dolbeault

► To cite this version:

Albert Cohen, Matthieu Dolbeault. Optimal pointwise sampling for $L^2$ approximation. Journal of Complexity, 2022, 68, pp.101602. hal-03221181

HAL Id: hal-03221181

https://hal.science/hal-03221181

Submitted on 7 May 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Optimal pointwise sampling for $L^2$ approximation

Albert Cohen and Matthieu Dolbeault

May 7, 2021

## Abstract

Given a function $u \in L^2 = L^2(D, \mu)$, where $D \subset \mathbb{R}^d$ and $\mu$ is a measure on $D$, and a linear subspace $V_n \subset L^2$ of dimension $n$, we show that near-best approximation of $u$ in $V_n$ can be computed from a near-optimal budget of $Cn$ pointwise evaluations of $u$, with $C > 1$ a universal constant. The sampling points are drawn according to some random distribution, the approximation is computed by a weighted least-squares method, and the error is assessed in expected $L^2$ norm. This result improves on the results in [6, 8] which require a sampling budget that is sub-optimal by a logarithmic factor, thanks to a sparsification strategy introduced in [17, 18]. As a consequence, we obtain for any compact class $\mathcal{K} \subset L^2$ that the sampling number $\rho_{Cn}^{\text{rand}}(\mathcal{K})_{L^2}$ in the randomized setting is dominated by the Kolmogorov $n$-width $d_n(\mathcal{K})_{L^2}$. While our result shows the existence of a randomized sampling with such near-optimal properties, we discuss remaining issues concerning its generation by a computationally efficient algorithm.

MSC 2020: 41A65, 41A81, 93E24, 62E17, 94A20

# 1 Introduction

We study the approximation of a function $u \in L^2(D, \mu)$, where $D$ is a domain in $\mathbb{R}^d$ and $\mu$ a measure on $D$, by an element $\tilde{u}$ of $V_n$, a subspace of $L^2(D, \mu)$ of finite dimension $n$, based on pointwise data of $u$. Therefore, to construct $\tilde{u}$, we are allowed to evaluate $u$ on a sample of $m$ points $X = \{x^1, \dots, x^m\} \in D^m$. In addition, we consider randomized sampling and reconstruction, in the sense that $X$ will be drawn according to a distribution $\sigma$ over $D^m$, so the error $u - \tilde{u}$ should
be evaluated in some probabilistic sense. For the sake of notational simplicity, having fixed $D$ and $\mu$, we write throughout the paper

$$L^2 := L^2(D, \mu) \quad \text{and} \quad \|v\| := \|v\|_{L^2} = \left( \int_D |v|^2 d\mu \right)^{1/2}, \qquad (1)$$

as well as

$$e_n(u) := \min_{v \in V_n} \|u - v\|. \qquad (2)$$

Our main result is the following:

**Theorem 1.** For some universal constants $C, K \ge 1$, and for any $n$-dimensional space $V_n \subset L^2$, there exists a random sampling $X = \{x^1, \dots, x^m\}$ with $m \le Cn$ and a reconstruction map $R: D^m \times \mathbb{C}^m \to V_n$, such that for any $u \in L^2$,

$$\mathbb{E}_X (\|u - \tilde{u}\|^2) \le K e_n(u)^2 \qquad (3)$$

where $\tilde{u} := R(x^1, \dots, x^m, u(x^1), \dots, u(x^m))$.

The reconstruction map $R$ is obtained through a weighted least-squares method introduced in [6], which has already been discussed in several papers, see [1, 5, 8, 9, 7, 15, 16]. The weights involved are given by the expression

$$w: x \in D \mapsto n \min_{v \in V_n} \frac{\|v\|^2}{|v(x)|^2} = \frac{n}{\sum_{j=1}^n |L_j(x)|^2}, \qquad (4)$$

where the last formula holds for any $L^2$-orthonormal basis $(L_1, \dots, L_n)$ of $V_n$, which, up to the factor $n$, is the *Christoffel function* associated to the space $V_n$ and the space $L^2(D, \mu)$. The weighted least-squares solution is then simply defined as

$$\tilde{u} := \underset{v \in V_n}{\operatorname{argmin}} \sum_{i=1}^{m} w(x^i) |u(x^i) - v(x^i)|^2. \qquad (5)$$

Introducing the discrete $\ell^2$ norm

$$\|v\|_{X}^{2} := \frac{1}{m} \sum_{i=1}^{m} w(x^{i}) |v(x^{i})|^{2} \qquad (6)$$

and its associated scalar product $\langle \cdot, \cdot \rangle_X$, we get a computable formula for $\tilde{u}$:

$$\tilde{u} = \underset{v \in V_n}{\operatorname{argmin}} \|u - v\|_{X}^{2} = P_{V_n}^{X} u, \qquad (7)$$
where $P_{V_n}^X$ denotes the orthogonal projection on $V_n$ with respect to $\langle \cdot, \cdot \rangle_X$. Note that strictly speaking $\| \cdot \|_X$ is not a norm over $L^2$, however the existence and uniqueness of $P_{V_n}^X$ will be ensured by the second condition in Lemma 1 below, see Remark 2.

Therefore our main achievements lie in the particular choice of the random sample $X$ for ensuring the near-optimal approximation and sampling budget in Theorem 1.

Now, the proof of Theorem 1 relies on two conditions: first, the expectation of $\| \cdot \|_X^2$ has to be bounded by $\| \cdot \|^2$ up to a constant. Second, an inverse bound should hold almost surely, instead of just in expectation, for functions $v$ in $V_n$. More precisely, one has:

**Lemma 1.** Assume that $m$ and the law $\sigma$ of $X = \{x^1, \dots, x^m\}$ are such that

$$ \mathbb{E}(\|v\|_X^2) \le \alpha \|v\|^2, \quad v \in V, \qquad (8) $$

and

$$ \|v\|^2 \le \beta \|v\|_X^2 \quad \text{a.s.,} \quad v \in V_n. \qquad (9) $$

Then

$$ \mathbb{E}(\|u - \tilde{u}\|^2) \le (1 + \alpha\beta)e_n(u)^2. \qquad (10) $$

*Proof.* Denote $u_n$ the orthogonal projection of $u$ on $V_n$ with respect to the $L^2(D, \mu)$ norm. Applying Pythagoras theorem both for $\| \cdot \|$ and $\| \cdot \|_X$, one obtains

$$ \begin{aligned} \mathbb{E}(\|u - \tilde{u}\|^2) &= \|u - u_n\|^2 + \mathbb{E}(\|u_n - \tilde{u}\|^2) \\ &\le \|u - u_n\|^2 + \beta \mathbb{E}(\|u_n - \tilde{u}\|_X^2) \\ &= \|u - u_n\|^2 + \beta \mathbb{E}(\|u_n - u\|_X^2 - \|u - \tilde{u}\|_X^2) \\ &\le \|u - u_n\|^2 + \beta \mathbb{E}(\|u_n - u\|_X^2) \\ &\le (1 + \alpha\beta)\|u - u_n\|^2, \end{aligned} $$

which proves (10) since $\|u - u_n\| = e_n(u)$. $\square$

In section 2, we recall how both conditions (8) and (9) can be obtained with $m$ quasi-linear in $n$, that is, of order $n \log n$. We reduce this budget to $m$ of order $n$ in section 3, by randomly subsampling the set of evaluation points, based on results from [17, 18]. The proof of Theorem 1 follows. We compare it to the recent results [13, 14, 19] in section 4, in particular regarding the domination of sampling numbers by $n$-widths. We conclude in section 5 by a discussion on the offline computational cost for practically generating the sample $X$.
## 2 Weighted least-squares

A first approach consists in drawing the $x^i$ independently according to the same distribution $\rho$, that is, taking $\sigma = \rho^{\otimes m}$. The natural choice for $\rho$ is $d\rho = \frac{1}{w}d\mu$, which is a probability measure since

$$ \int_D \frac{1}{w} d\mu = \frac{1}{n} \sum_{j=1}^{n} \int_D |L_j(x)|^2 d\mu(x) = \frac{1}{n} \sum_{j=1}^{n} \|L_j\|^2 = 1. \quad (11) $$

We denote by $Z = \{x^1, \dots, x^m\}$ this first random sample and by $\| \cdot \|_Z$ the corresponding discrete $\ell^2$ norm. With this sampling measure,

$$ \mathbb{E}(\|v\|_Z^2) = \frac{1}{m} \sum_{i=1}^{m} \int_D w(x) |v(x)|^2 d\rho = \int_D |v|^2 d\mu = \|v\|^2, \quad (12) $$

so condition (8) is ensured for $X = Z$ with $\alpha = 1$. To study the second condition, we introduce the Hermitian positive semi-definite Gram matrix

$$ G_Z := (\langle L_j, L_k \rangle_Z)_{j,k=1,\dots,n} \quad (13) $$

and notice that (9) is equivalent to

$$ |\nu|^2 = \left\| \sum_{j=1}^{n} \nu_j L_j \right\|^2 \le \beta \left\| \sum_{j=1}^{n} \nu_j L_j \right\|_Z^2 = \beta \nu^* G_Z \nu, \quad \nu \in \mathbb{C}^n, \quad (14) $$

which in turn rewrites as $\lambda_{\min}(G_Z) \ge \beta^{-1}$.

By the central limit theorem, as $m$ tends to infinity, the scalar products $\langle L_j, L_k \rangle_Z$ converge almost surely to $\langle L_j, L_k \rangle = \delta_{j,k}$, so $G_Z$ converges to the identity matrix, and we expect that $\lambda_{\min}(G_Z) \ge \beta^{-1}$ holds for $\beta > 1$ with high probability as $m$ gets large. A quantitative formulation can be obtained by studying the concentration of $G_Z$ around $I$ in the matrix spectral norm

$$ \|M\|_2 := \max\{|Mx| : |x| = 1\}. $$

This is based on the matrix Chernoff bound, see [2, 22] for the original inequality and [5], Lemma 2.1, for its application to our problem:

**Lemma 2.** For $m \ge 10n \ln(\frac{2n}{\epsilon})$, if $X \sim \rho^{\otimes m}$, then

$$ \mathbb{P} \left( \|G_Z - I\|_2 \le \frac{1}{2} \right) \ge 1 - \epsilon. \quad (15) $$

In particular, $\mathbb{P}(\lambda_{\min}(G_Z) \ge \frac{1}{2}) \ge 1 - \epsilon$.
Thus assumption (9) is satisfied with $\beta = 2$, but only with probability $1 - \varepsilon$. As we would like it to hold almost surely, we condition the sampling to the event

$$E := \left\{ \|G_Z - I\|_2 \le \frac{1}{2} \right\} \qquad (16)$$

which defines a new sample

$$Y = Z|E. \qquad (17)$$

In practice, $Y$ can be obtained through a rejection method, which consists in drawing successively set of points $Z^1, Z^2, \dots$ according to $\rho^{\otimes m}$, and defining $Y = Z^k$ for the first value $k$ such that $E$ holds. We then define $\tilde{u}$ as the weighted least-square estimator based on this conditioned sample, that is

$$\tilde{u} := P_{V_n}^Y u. \qquad (18)$$

This approach was introduced and analyzed in [8], see in particular Theorem 3.6 therein. A simpler version of this result, sufficient for our purposes, is the following:

**Lemma 3.** For $m \ge 10n \ln(4n)$, if $Z \sim \rho^{\otimes m}$ and $Y = Z|E$, then

$$\|G_Y - I\|_2 \le \frac{1}{2}, \qquad (19)$$

and

$$\mathbb{E}_Y(\|u - \tilde{u}\|^2) \le 5e_n(u)^2. \qquad (20)$$

*Proof.* The first part immediately results from the definition of $Y$ and $E$, and implies condition (9) with $\beta = 2$. Moreover, $\mathbb{P}(E) \ge \frac{1}{2}$ by Lemma 2 with $\varepsilon = \frac{1}{2}$, so for any $v \in L^2(D, \mu)$,

$$\mathbb{E}_Y(\|v\|_Y^2) = \mathbb{E}_Z(\|v\|_Z^2|E) = \frac{\mathbb{E}_Z(\|v\|_Z^2 \chi_E)}{\mathbb{P}(E)} \le \frac{\mathbb{E}_Z(\|v\|_Z^2)}{\mathbb{P}(E)} \le 2\|v\|^2, \qquad (21)$$

so condition (8) holds with $\alpha = 2$. The conclusion follows from Lemma 1. $\square$

**Remark 1.** The number of redraws $k$ for reaching $Y$ follows a geometric law of expectation $\mathbb{E}(k) = \mathbb{P}(E)^{-1} = (1 - \varepsilon)^{-1}$, that is $\mathbb{E}(k) \le 2$ for the particular choice of $m$ in the above lemma. It should be well noted that $u$ is not evaluated at the intermediately generated samples $Z^1, \dots, Z^{k-1}$, which thus enter the offline cost of the samping algorithm.
**Remark 2.** The fact that the Grammian $G_Y$ is non-singular implies that for any $u$ with given values $y^i$ at the points $x^i$, we can uniquely define

$$ \tilde{u} = P_{V_n}^Y u = \sum_{j=1}^{n} a_j L_j, \qquad (22) $$

since $a = (a_1, \dots, a_n)^*$ solves the system of normal equations

$$ G_Y a = b, \qquad (23) $$

where the right-side vector has coordinates

$$ b_j = \langle L_j, u \rangle_Y = \frac{1}{m} \sum_{i=1}^{m} w(x^i) L_j(x^i) \overline{y^i}. \qquad (24) $$

If $u$ is in $L^2$, the $y^i$ are only defined up to a representer, however since two representers $u^1$ and $u^2$ coincide $\mu$-almost surely, we find that $P_{V_n}^Y u$ is well defined almost surely over the draw of $Y$.

### 3 Random subsampling

With Lemma 3, we already have an error bound similar to that of Theorem 1. However, the sampling budget is larger than $n$ by a logarithmic factor, which we seek to remove in this section. To do so, we partition the sample $Y$ into subsets of size comparable to $n$, and randomly pick one of these subsets to define the new sample. An appropriate choice of the partitioning is needed to circumvent the main obstacle, namely the preservation of condition (9). It relies on the following lemma, taken from Corollary B of [18], itself a consequence of Corollary 1.5 in [17]. The relevance of these two results to sampling problems were exploited in [19] and noticed in [10], respectively.

**Lemma 4.** Let $a_1, \dots, a_m \in \mathbb{C}^n$ be vectors of norm $|a_i|^2 \le \delta$ for $i = 1, \dots, m$, and satisfying

$$ \alpha I \le \sum_{i=1}^{m} a_i a_i^* \le \beta I \qquad (25) $$

for some constants $\delta < \alpha \le \beta$. Then there exists a partition of $\{1, \dots, m\}$ into two sets $S_1$ and $S_2$ such that

$$ \frac{1 - 5\sqrt{\delta/\alpha}}{2} \alpha I \le \sum_{i \in S_j} a_i a_i^* \le \frac{1 + 5\sqrt{\delta/\alpha}}{2} \beta I, \quad j = 1, 2. \qquad (26) $$
In Lemma 2 of [18] this result is applied inductively in order to find a smaller set $J \subset \{1, \dots, k\}$ of cardinality $|J| \le cn$ such that

$$ C^{-1}I \le \frac{m}{n} \sum_{i \in J} a_i a_i^* \le CI, \qquad (27) $$

for some universal constants $c, C > 1$. We adapt this approach in order to obtain a complete partition of $\{1, \dots, m\}$ by sets having such properties.

**Lemma 5.** Let $a_1, \dots, a_m \in \mathbb{C}^n$ be vectors of norm $|a_i|^2 = \frac{n}{m}$ for $i = 1, \dots, m$, and satisfying

$$ \frac{1}{2} I \le \sum_{i=1}^{m} a_i a_i^* \le \frac{3}{2} I. \qquad (28) $$

Then there exists an integer $L$ and a partition of $\{1, \dots, m\}$ into $2^L$ sets $J_1, \dots, J_{2^L}$ such that

$$ c_0 I \le \frac{m}{n} \sum_{i \in J_k} a_i a_i^* \le C_0 I, \quad 1 \le k \le 2^L, \qquad (29) $$

with universal constants $c_0$ and $C_0$. In addition, each set $J_k$ satisfies

$$ |J_k| \le C_0 n. \qquad (30) $$

*Proof.* The cardinality estimate (30) follows from the upper inequality in (29) by taking the trace

$$ nC_0 = \operatorname{tr}(C_0 I) \ge \frac{m}{n} \sum_{i \in J_k} \operatorname{tr}(a_i a_i^*) = \frac{m}{n} |J_k| \frac{n}{m} = |J_k|. \qquad (31) $$

For the proof of (29), if $n/m \ge 1/200$, then the result holds with $L = 0$, $J_1 = \{1, \dots, m\}$ and $c_0 = C_0 = 200$. Now assuming $\delta := n/m < 1/200$, define by induction $\alpha_0 = \frac{1}{2}, \beta_0 = \frac{3}{2}$, and

$$ \alpha_{\ell+1} := \alpha_{\ell} \frac{1 - 5\sqrt{\delta/\alpha_{\ell}}}{2}, \quad \beta_{\ell+1} := \beta_{\ell} \frac{1 + 5\sqrt{\delta/\alpha_{\ell}}}{2}, \quad \ell \ge 0. \qquad (32) $$

As $\alpha_{\ell+1} \le \frac{\alpha_{\ell}}{2}$, the minimal integer $L$ such that $\alpha_L \le 100\delta$ is well defined, and satisfies

$$ \alpha_L = \alpha_{L-1} \frac{1 - 5\sqrt{\delta/\alpha_{L-1}}}{2} > 100\delta \frac{1 - 5\sqrt{1/100}}{2} = 25\delta. \qquad (33) $$
Moreover $\alpha_{\ell} \ge 2^{L-\ell-1}\alpha_{L-1} \ge 2^{L-\ell-1} 100\delta$ for $\ell = 0, \dots, L-1$, so

$$ \beta_L = 3\alpha_L \prod_{\ell=0}^{L-1} \frac{1 + 5\sqrt{\delta/\alpha_\ell}}{1 + 5\sqrt{\delta/\alpha_\ell}} \le C\delta, \quad (34) $$

with $C := 300 \prod_{\ell \ge 2} \frac{1+\sqrt{2}^{-\ell}}{1-\sqrt{2}^{-\ell}}$.

Finally, we inductively define partitions $\{S_1^\ell, \dots, S_{2^\ell}^\ell\}$ for $0 \le \ell \le L$: start with $S_1^0 = \{1, \dots, m\}$ and for any $\ell, j$, noticing that

$$ \alpha_{\ell} I \le \sum_{i \in S_{j}^{\ell}} a_{i} a_{i}^{*} \le \beta_{\ell} I, \quad (35) $$

apply Lemma 4 to split $S_j^\ell$ into subsets $S_{2j-1}^{\ell+1}$ and $S_{2j}^{\ell+1}$ satisfying the same property. At the last step, we define

$$ J_k = S_k^L. \quad (36) $$

The framing (29) thus holds with $c_0 = \alpha_L/\delta \ge 25$ and $C_0 = \beta_L/\delta \le 11000$. $\square$

**Proof of Theorem 1:** Define

$$ a_i = \left( \sqrt{\frac{w(x^i)}{m}} L_j(x^i) \right)_{j=1,\dots,n} \quad (37) $$

the normalised random vectors corresponding to the sample $Y = \{x^1, \dots, x^m\}$ introduced in the previous section. As

$$ \frac{1}{2} I \le G_Y = \sum_{i=1}^{m} a_i a_i^* \le \frac{3}{2} I \quad (38) $$

and

$$ |a_i|^2 = \frac{1}{m} w(x^i) \sum_{j=1}^{n} |L_j(x^i)|^2 = \frac{n}{m} \quad (39) $$

thanks to the choice of weights (4), the assumptions of Lemma 5 are satisfied. Applying this lemma, we obtain sets $J_1, \dots, J_{2^L}$ partitioning $\{1, \dots, m\}$. Let $\kappa$ be a random variable taking value $k \in \{1, \dots, 2^L\}$ with probability $p_k = |J_k|/m$, and create a random subsampling $X$ of $Y$ through

$$ X = \{x^i \in Y : i \in J_{\kappa}\}. \quad (40) $$
Then the budget condition $|X| = |J_κ| ≤ C_0n$ is satisfied according to (30). Here, we define the discrete norm as

$$ \|v\|_{X}^{2} := \frac{1}{|X|} \sum_{i \in J_{\kappa}} w(x^{i})|v(x^{i})|^{2}, \qquad (41) $$

and the associated Gram matrix

$$ G_X := (\langle L_j, L_k \rangle_X)_{j,k=1,...,n} = \frac{m}{|J_κ|} \sum_{i \in J_κ} a_i a_i^*. \qquad (42) $$

The weighted least-squares estimate is now defined as

$$ \tilde{u} := \underset{v \in V_n}{\operatorname{argmin}} \frac{1}{|X|} \sum_{i \in J_κ} w(x^i) |u(x^i) - v(x^i)|^2, \qquad (43) $$

and it thus depends on the random draws of both $Y$ and $κ$. Condition (9) follows from the lower inequality in (29) with $β = \frac{C_0}{c_0}$ since

$$ G_X \geq \frac{m}{|J_κ|} \frac{n}{m} c_0 I \geq \frac{c_0}{C_0} I. \qquad (44) $$

Finally, we have for any $v \in L^2(D, \mu)$

$$ \mathbb{E}_X (\|v\|_X^2) = \sum_{k=1}^{2^L} \frac{p_k}{|J_k|} \sum_{i \in J_k} w(x^i) |v(x^i)|^2 = \mathbb{E}_Y (\|v\|_Y^2) \leq 2e_n(u)^2, \qquad (45) $$

so condition (8) holds with $α = 2$. Applying Lemma 1, we conclude that (10) holds with $C = C_0$ and $K = 1 + 2\frac{C_0}{c_0}$. $\square$

# 4 Comparison with related results

In order to compare Theorem 1 with several recent results [10, 12, 19, 23], we consider its implication when the target function $u$ belongs to a certain class of functions $K$ that describes some prior information on $u$, such as smoothness.

Recall that if $V$ is a Banach space of functions defined on $D$ and $K \subset V$ is a compact set, its *Kolmogorov n-width* is defined by

$$ d_n(K)_V := \inf_{\dim V_n = n} \sup_{u \in K} \inf_{v \in V_n} \|u - v\|_V, \qquad (46) $$
where the first infimum is taken over all linear spaces $V_n \subset V$ of dimension $n$. This quantity thus describes the best approximation error that can be achieved uniformly over the class $\mathcal{K}$ by an $n$-dimensional linear space.

On the other hand, building a best approximation of $u$ requires in principle full knowledge on $u$, and we want to consider the situation where we only have access to a limited number of point evaluations. This leads one to consider the *sampling numbers*, also called *optimal recovery numbers*, both in the deterministic and randomized settings.

For deterministic samplings, we define the (linear) sampling numbers

$$ \rho_m^{\det}(\mathcal{K})_V := \inf_{X,R_X} \max_{u \in \mathcal{K}} \|u - R_X(u(x^1), \dots, u(x^m))\|_V, \quad (47) $$

where the infimum is taken over all samples $X = \{x^1, \dots, x^m\} \in D^m$ and linear reconstruction maps $R_X : \mathbb{C}^m \to V$. For random samplings, we may define similar quantities by

$$ \rho_m^{\text{rand}}(\mathcal{K})_V^2 := \inf_{X,R_X} \max_{u \in \mathcal{K}} \mathbb{E} \left( \|u - R_X(u(x^1), \dots, u(x^m))\|_V^2 \right), \quad (48) $$

where the infimum is taken over all random variables $X = \{x^1, \dots, x^m\} \in D^m$ and linear reconstruction maps $R_X : \mathbb{C}^m \to V$. Note that a deterministic sample can be viewed as a particular choice of random sample following a Dirac distribution in $D^m$, and therefore

$$ \rho_m^{\text{rand}}(\mathcal{K})_V \le \rho_m^{\det}(\mathcal{K})_V. \quad (49) $$

Sampling numbers may also be defined without imposing the linearity of $R_X$, leading to smaller quantities. In what follows, we shall establish upper bound on the linear sampling numbers, which in turns are upper bounds for the nonlinear ones. We refer to [20] for an introduction and study of sampling numbers in the context of general linear measurements, and to [21] that focuses on point evaluation, also termed as *standard information*.

By optimizing the choice of the space $V_n$ used in Theorem 1, we obtain as a consequence that, for $V = L^2 = L^2(D, \mu)$, the sampling numbers in the randomized setting are dominated by the Kolmogorov $n$-widths.

**Corollary 1.** For any compact set $\mathcal{K} \subset L^2$, one has

$$ \rho_{Cn}^{\text{rand}}(\mathcal{K})_{L^2} \le K d_n(\mathcal{K})_{L^2}, \quad (50) $$

where $C$ and $K$ are the same constants as in Theorem 1.
**Remark 3.** The bound (50) cannot be attained with independent and identically distributed sampling points $x^1, \dots, x^m$. Indeed, consider the simple example, already evoked in [22], where $D = [0, 1]$, $\mu$ is the Lebesgue measure,

$$V_n = \left\{ \sum_{i=1}^{n} a_i \chi_{[\frac{i}{n}, \frac{i+1}{n}]} , (a_1, \dots, a_n) \in \mathbb{C}^n \right\} \quad (51)$$

is a space of piecewise constant functions, and $K = \{u \in V_n, \|u\|_{L^\infty} \le 1\}$. Then $K \subset V_n$ so $d_n(K)_{L^2} = 0$, and an exact reconstruction $R_X u = u$ is possible if and only if $X$ contains at least one point in each interval $[\frac{i}{n}, \frac{i+1}{n}]$. Thus $\rho_n^{\det}(K)_{L^2} = 0$, but in the case of i.i.d measurements, $m$ has to grow like $n \log n$ to ensure this constraint, due to the coupon collector's problem.

**Remark 4.** In [11], a result similar to Theorem 1 is obtained under the extra assumption of a uniform bound on $e_n(u)/e_{2n}(u)$, yielding the validity of (50) assuming a uniform bound on $d_n(K)_{L^2}/d_{2n}(K)_{L^2}$. The recovery method used in [11] is not of least-square type, but rather an elaboration of the pseudo-spectral approach that would simply approximate the inner products $\langle u, L_j \rangle = \int_D uL_jd\mu$ by a quadrature, using a hierarchical approach introduced in [24].

Ideally, one would like a “worst case” or “uniform” version of Theorem 1, in the form

$$\rho_{C_n}^{\det}(K)_{L^2} \le K d_n(K)_{L^2}, \quad (52)$$

but it is easily seen that such an estimate cannot be expected for general compacts sets of $L^2$, due to the fact that pointwise evaluations are not continuous in $L^2$ norm.

It is however possible to recover such uniform estimates by mitigating the non-achievable estimate (52) in various ways. One first approach, developed in [14, 23], gives an inequality similar to (52), with $d_n(K)_{L^2}$ replaced by $d_n(K)_{L^\infty}$. It is based on the following lemma, see Theorem 2.1 in [23], which we recall for comparison with our Lemma 1:

**Lemma 6.** Assume that $\mu$ is a finite measure of mass $\mu(D) = M < \infty$, that the constant functions belong to $V_n$, and that there exists a sample $X = \{x^1, \dots, x^m\}$ and weights $w^i$ such that the discrete norm

$$\|v\|_{X}^{2} = \frac{1}{|X|} \sum_{i=1}^{m} w_{i} |v(x^{i})|^{2} \quad (53)$$
satisfies a framing

$$ \beta^{-1} \|v\|^2 \le \|v\|_X^2 \le \alpha \|v\|^2, \quad v \in V_n. \quad (54) $$

Then

$$ \|u - P_{V_n}^X u\| \le \sqrt{M} (1 + \sqrt{\alpha\beta}) e_n(u)_{L^\infty}, \quad (55) $$

where $e_n(u)_{L^\infty} = \min_{v \in V_n} \|u - v\|_{L^\infty}$.

*Proof.* For any $v \in L^2$, we have $\|v\|^2 \le M \|v\|_{L^\infty}^2$, and as $1 \in V_n$,

$$ \|v\|_X^2 \le \|1\|_X^2 \|v\|_{L^\infty}^2 \le \alpha \|1\|^2 \|v\|_{L^\infty}^2 = \alpha M \|v\|_{L^\infty}^2. \quad (56) $$

Hence

$$ \begin{aligned} \|u - P_{V_n}^X u\| &\le \|u - v\| + \|v - P_{V_n}^X u\| \\ &\le \|u - v\| + \sqrt{\beta} \|v - P_{V_n}^X u\|_X \\ &\le \|u - v\| + \sqrt{\beta} \|v - u\|_X \\ &\le (\sqrt{M} + \sqrt{\alpha\beta M}) \|u - v\|_{L^\infty}, \end{aligned} $$

and we conclude by optimizing over $v \in V_n$. $\square$

<figure><img src="image_3.png" /></figure>

Here, in contrast to the derivation of (10) in Lemma 1, one only uses the framing property (54), and does not need the condition $\mathbb{E}(\|v\|_X^2) \le \alpha \|v\|^2$. For this reason, one may achieve the above objective with a simpler sparsification approach proposed in [4] and adapted in [14], which performs a greedy selection of the points $x^i$ within the sample Y, together with the definition of weights $w_i$ associated with these points. If the initial sample Y satisfies

$$ \frac{1}{2} I \le G_Y \le \frac{3}{2} I, \quad (57) $$

then, for any $c > 1$ the selection algorithm produces a sample X of at most $cn$ points such that (54) holds with $\alpha = \frac{3}{2} \left(1 + \frac{1}{\sqrt{c}}\right)^2$ and $\beta^{-1} = \frac{1}{2} \left(1 - \frac{1}{\sqrt{c}}\right)^2$.

Optimizing the choice of $V_n$ (but imposing that constant functions are contained in this space), this leads to the following comparison result between deterministic optimal recovery numbers in $L^2$ and n-widths in $L^\infty$: for any compact set $K \in C(D)$, one has

$$ \rho_{cn}^{\det}(K)_{L^2} \le C \sqrt{M} d_{n-1}(K)_{L^\infty}, \quad (58) $$
where C depends on c > 1. For c = 2, one can take C = 11. We refer to [14, 23] where this type of result is established.

Another approach consists in making pointwise evaluations continuous by restriction to the case where $\mathcal{K} = B_H$ is the unit ball of a reproducing kernel Hilbert space $H \subset L^2$, and assuming that the sequence $(d_n(B_H)_{L^2})_{n \ge 1}$ is $\ell^2$-summable. The following result from [19], also based on the sparsification techniques from [17], improves on a bound found in [12]

$$ \rho_{C_n}^{\det}(B_H)^2 \le K \frac{\log n}{n} \sum_{k \ge n} d_k(B_H)_{L^2}^2, \quad (59) $$

In [13], a similar inequality is established for more general compact classes $\mathcal{K}$ of $L^2$ such that point evaluation is well defined on functions of $\mathcal{K}$, however now with $\ell^p$ sums with $p < 2$ replacing $\ell^2$ sums.

In the above inequality, the logarithmic factor appears as a residual of the result obtained before sparsification, contrarily to the bounds (50) and (58), which do not explicitely depend on the size of the initial sample $Y$. This results in a gap of a factor $\log n$ between (59) and known lower bounds for $\rho_{C_n}^{\det}(B_H)^2$, see [19].

# 5 Computational aspects

The various results (50), (58), (59) ensure the existence of good sampling and reconstructions algorithms in various settings. We end by a discussion on the computational cost of these strategies.

For the weighted least-squares methods corresponding to samples $Z$ and $Y$, the most expensive step consists in assembling the matrix $G_Z$ as a sum of $m$ matrices of size $n$, so the algorithmic complexity is of order $O(mn^2) = O(n^3 \log n)$. Besides, to obtain $G_Y$, this step may need to be repeated a few times, as explained in Remark 1, but this only affects the offline complexity by a small random factor.

Note that we assumed that an orthogonal basis $(L_1, \dots, L_n)$ of $V_n$ is explicitly known, which might not be the case for irregular domains $D$. However, under reasonable assumptions on $D$ or $V_n$, one can compute an approximately orthogonal basis $(\tilde{L}_1, \dots, \tilde{L}_n)$, either by performing a first discretization of $D$ with a large number of points, or by using a hierarchical method on a sequence of nested spaces $V_1 \subset \dots \subset V_n$, see [1, 3, 5, 8, 15, 16]. These additional steps have complexities $O(K_n n^2)$ and $O(n^4)$ respectively, where $K_n$ is the maximal value of the inverse Christoffel function $\sum_{j=1}^n |L_j|^2$ which might grow more than linearly with
$n$ for certain choices of spaces $V_n$. Results similar to Lemma 3 have been obtained in the above references, with $(L_j)_{j=1,...,n}$ replaced by $(\tilde{L}_j)_{j=1,...,n}$.

One could stop at this point and compute the approximation $\tilde{u} = P_{V_n}^Y u$, which satisfies error bounds both in expectation when comparing to $e_n(u)$, see Lemma 3, or uniformly when comparing to $e_n(u)_{L^\infty}$, see Theorem 1 (iii) in [6]. Once the measurements of $u$ are performed, the computation of $\tilde{u}$ requires to solve a $n \times n$ linear system as in Remark 2, so the online stage takes a time $O(\tau n \log n + n^3)$, where $\tau$ is the cost of each measurement of $u$.

However, in applications where the evaluation cost $\tau$ becomes very high (for example when each evaluation $x \mapsto u(x)$ requires solving a PDE by some numerical code, or running a physical experiment), further reduction of the size of the sample may prove interesting, and justifies the interest for sparsification methods. The greedy selection method from [4], which is used in [23] and leads to (58), has a complexity in $O(mn^3) = O(n^4 \log n)$, but it can only be applied to the worst-case setting, with the uniform error bound $e_n(u)_{L^\infty}$.

On the other hand, the iterative splitting method that we have used in this paper following the ideas from [17, 19] is not easily implemented, and one obvious method consists in testing all partitions of $\{1, ..., m\}$ into sets $S_1$ and $S_2$ when applying Lemma 4. Note that this lemma is in practice used $L$ times, with $L = O(\log \log n)$ since $2^L = O(m/n) = O(\log n)$. The algorithm consisting in subdividing the sample $L$ times, each time checking that the Gram matrices corresponding to $S_1$ and $S_2$ are well conditioned, and keeping one such subset at random, thus has an exponential complexity $O(2^m n^3) = O(n^{cn})$. Having a different strategy that would produce the random sample in polynomial time is currently an open problem to us.

We summarise these computational observations in the following table, which illustrates the conflicts between reducing the sampling budget, ensuring optimal approximation results, and maintaining a reasonable cost for sample generation.

<table><thead><tr><th>sampling complexity</th><th>sample cardinality m</th><th>offline complexity</th><th>E(||u - ū||<sup>2</sup>) &le; Ce<sub>n</sub>(u)<sup>2</sup></th><th>||u - ū||<sup>2</sup> &le; Ce<sub>n</sub>(u)<sup>2</sup><sub>L<sup>\infty</sup></sub></th></tr></thead><tbody><tr><td>conditionned<br/>ρ<sup>⊗m</sup> | E</td><td>10n log(4n)</td><td>O(n<sup>3</sup> log n)</td><td>✓</td><td>✓</td></tr><tr><td>+ deterministic<br/>sparsification [4]</td><td>(1 + ε)n</td><td>O(n<sup>4</sup> log n)</td><td>✗</td><td>✓</td></tr><tr><td>+ random<br/>sparsification [17]</td><td>Cn</td><td>O(n<sup>cn</sup>) → O(n<sup>r</sup>) ?</td><td>✓</td><td>✓</td></tr></tbody></table>
As a final remark, let us to emphasize that although the results presented in our paper are mainly theoretical and not practically satisfactory, due both to the computational complexity of the sparsification, and to the high values of the numerical constants *C* and *K* in Theorem 1, they provide some intuitive justification to the boosted least-squares methods presented in [8], which consist in removing points from the initial sample as long as the corresponding Gram matrix *G*<sub>*X*</sub> remains well conditioned. For instance, Lemma 4 allows to keep splitting the sample even after *L* steps, if one still has a framing 1/2 *I* ≤ *G*<sub>*X*</sub> ≤ 3/2 *I* and a sufficiently large ratio |*X*|/n. Nevertheless, it would be of much interest to find a randomized version of [4] giving a bound of the form (50), since this would give algorithmic tractability, smaller values for *C* and *K*, and the possibility to balance these constants in Theorem 1.

References

[1] B. Adcock and J.M. Cardenas, *Near-optimal sampling strategies for multi-variate function approximation on general domains*, SIAM Journal on Mathematics of Data Science, 2 (2020), pp. 607-630.

[2] R. Ahlswede and A. Winter, *Strong converse for identification via quantum channels*, IEEE Transactions on Information Theory (2002), 48(3), pp. 569-579.

[3] B. Arras, M. Bachmayr and A. Cohen, *Sequential sampling for optimal weighted least squares approximations in hierarchical spaces*, SIAM Journal on Mathematics of Data Science, 1 (2019), pp. 189-207.

[4] J. Batson, D.A. Spielman, and N. Srivastava, *Twice-Ramanujan Sparsifiers*, SIAM Review, 56 (2014), 315-334.

[5] A. Cohen and M. Dolbeault, *Optimal sampling and Christoffel functions on general domains*, arXiv: 2010.11040 (2020).

[6] A. Cohen and G. Migliorati, *Optimal weighted least squares methods*, SMAI Journal of Computational Mathematics 3, 181–203, 2017.
[7] A. Doostan and J. Hampton, *Coherence motivated sampling and convergence analysis of least squares polynomial Chaos regression*, Computer Methods in Applied Mechanics and Engineering, 290 (2015), pp. 73-97.

[8] C. Haberstich, A. Nouy, and G. Perrin, *Boosted optimal weighted least-squares*, arXiv:1912.07075 (2019).

[9] J.D. Jakeman, A. Narayan, and T. Zhou, *A Christoffel function weighted least squares algorithm for collocation approximations*, Mathematics of Computation, 86 (2017), pp. 1913-1947.

[10] L. Kaemmerer, T. Ullrich, and T. Volkmer *Worst case recovery guarantees for least squares approximation using random samples*, arXiv: 1911.10111 (2019).

[11] D. Krieg, *Optimal Monte Carlo methods for L2 -approximation*, Constructive Approximation 49 (2019), pp. 385-403.

[12] D. Krieg and M. Ullrich, *Function values are enough for L2-approximation*, arXiv: 1905.02516 (2020), to appear in Foundations of Computational Mathematics.

[13] D. Krieg and M. Ullrich, *Function values are enough for L2-approximation: Part II*, Journal of Complexity (2021) https://doi.org/10.1016/j.jco.2021.101569.

[14] I. Limonova and V. Temlyakov, *On sampling discretization in L₂*, arXiv:2009.10789 (2020).

[15] G. Migliorati, *Adaptive approximation by optimal weighted least-squares methods*, SIAM Journal on Numerical Analysis, 7 (2019), pp. 2217-2245.

[16] G. Migliorati, *Multivariate approximation of functions on irregular domains by weighted least-squares methods*, IMA journal of numerical analysis (2020), https://doi.org/10.1093/imanum/draa023.

[17] A. Marcus, D. Spielman and N. Srivastava, *Interlacing families II: Mixed characteristic polynomials and the Kadison-Singer problem*, Annals of Mathematics (2015), http://dx.doi.org/10.4007/annals.2015.182.1.8.
[18] Shahaf Nitzan, Alexander Olevskii and Alexander Ulanovskii, *Exponential frames on unbounded sets*, Proceedings of the American Mathematical Society (2016), http://dx.doi.org/10.1090/proc/12868.

[19] N. Nagel, M. Schäfer, and T. Ullrich, *A new upper bound for sampling numbers*, arXiv: 2010.00327 (2020).

[20] E. Novak and H. Woźniakowski, *Tractability of multivariate problems. Volume I: Linear information*, EMS, Zürich, 2008.

[21] E. Novak and H. Woźniakowski, *Tractability of multivariate problems. Volume III: Standard information for operators*, EMS, Zürich, 2012.

[22] J. Tropp, *User-Friendly tail bounds for sums of random matrices*, Foundations of Computational mathematics, 12 (2012), pp. 389-434.

[23] V. N. Temlyakov, *On optimal recovery in $L^2$*, arXiv: 2010.03103 (2020).

[24] G.W. Wasilkowski and H. Woźniakowski, *The power of standard information for multivariate approximation in the randomized setting*, Mathematics of Computation, 76 (2007), pp. 965-988.