2021

Felix Bartel

# Stability and Error Guarantees for Least Squares Approximation with Noisy
Sampels

<EMAIL_ADDRESS>Chemnitz University of Technology,
Faculty of Mathematics, 09107 Chemnitz, Germany

###### Abstract

Given $n$ samples of a function $f\colon D\to\mathds{C}$ in random points
drawn with respect to a measure $\nu$ we develop theoretical analysis of the
$L_{2}(D,\mu)$-approximation error. We show that the weighted least squares
method from finite dimensional function spaces $V_{m}$, $\dim(V_{m})=m<\infty$
is stable and optimal up to a multiplicative constant when given exact samples
with logarithmic oversampling. Further, for noisy samples, our bounds describe
the bias-variance trade off depending on the dimension $m$ of the
approximation space $V_{m}$. All results hold with high probability.

For demonstration, we consider functions defined on the $d$-dimensional cube
given in unifom random samples. We analyze polynomials, the half-perid cosine,
and a bounded orthonormal basis of the non-periodic Sobolev space
$H_{\mathrm{mix}}^{2}$. Overcoming numerical issues of this
$H_{\mathrm{mix}}^{2}$ basis, this gives a novel stable approximation method
with quadratic error decay. Numerical experiments indicate the applicability
of our results.

###### keywords:

individual function approximation, least squares, sampling theory, unit cube,
polynomial approximation

###### pacs:

[

MSC Classification] 41A10, 41A25, 41A60, 41A63, 42C10, 65Txx, 65F22, 65D15,
94A20

##  1 Introduction

In this paper we study the reconstruction of complex-valued functions on a
$d$-dimensional domain $D\subset\mathds{R}^{d}$ from possibly noisy function
values

$\displaystyle\bm{y}=\bm{f}+\bm{\varepsilon}=(f(\bm{x}^{1})+\varepsilon_{1},\dots,f(\bm{x}^{n})+\varepsilon_{n})^{\mathsf{T}}\,,$

which are sampled in random points $\bm{x}^{1},\dots,\bm{x}^{n}\in D$. We
consider error bounds for the weighted least squares method for individual
functions, which is common in, e.g. partial differential equations HD15  or
uncertainty quantification HNTW20 . In this setting, the samples are drawn
after the function is fixed in contrast to worst-case or minmax-bounds, which
hold for a class of functions and usually do not include noise in the samples.
For individual function approximation the majority of error bounds are stated
in expectation, cf. (Baraud02, , Thm. 1.1) for penalized least-squares,
(CDL13, , Thm. 3) for plain least-squares, or (CM17, , Thm. 2.1), (HNP22, ,
Thm. 4.1), and (KUV21, , Thm. 6.1) for weighted least squares. Stronger
bounds, which hold with high probability, are only known in specific settings
and do not include noise: cf. (MNST14, , Thm. 3) for polynomial approximation
and (LPU21, , Thms. 3.20 & 3.21) for wavelet approximation. The contribution
and novelty of this paper is twofold:

  * •

We use concentration inequalities to show error bounds in the $L_{2}$\- and
$L_{\infty}$-norm which hold with high probability, include noisy samples, and
are applicable in the general setting of the known bounds in expectation. We
do not loose the optimality of the weighted least-squares method when samples
are noise-free, i.e., we only need logarithmic oversampling to achieve the
best possible error. (Note, there exists a set of points with linear
oversampling achieving the optimal error but this is not constructive nor
random, cf. DC22 .) Including noise, our bounds reflect the typical bias-
variance trade off which one wants to balance to prevent over- or
underfitting. The results enable to give performance guarantees for model
selection strategies like the balancing principle LMP20 ; PL13  or cross-
validation BaHiPo19 ; BH21 .

  * •

For an application we have a look at approximation on the $d$-dimensional unit
cube $[0,1]^{d}$ when samples are distributed uniform according to the
Lebesgue measure. Here, we improve on the results in MNST14  for polynomial
approximation in the noise-free setting by a logarithmic factor and, more
importantly, the approximation error there is related to the coarser best
possible error in $L_{\infty}$, similar to the general results in Temlyakov21
, instead of $L_{2}$. A drawback of polynomials is the need for quadratic
oversampling, which we show for the Legendre polynomials but holds in general,
cf. MNST14 . To circumvent this, we use the eigenfunctions of the embedding
$\mathrm{Id}\colon H^{s}\to L_{2}$ into the Sobolev space $H^{s}$ for $s=1,2$
which allow for logarithmic oversampling. The $H^{1}$ basis, also known as
half-period cosine, has become the standard in many applications and is
researched thouroughly, cf. WW08 ; IN08 ; Adcock10 ; AH11 ; DDP14 ; SNC16 ;
CKNS16 ; KMNN21 . But also for functions in Sobolev spaces $H^{s}$ of higher
smoothness their convergence is limited to be linear in theory (the rate $3/2$
can be observed in practice). This can be improved by using the $H^{2}$ basis,
examined theoretically in (AIN12, , Section 3) to have quadratic convergence.
But it is not used as it is prone to numerical errors and unusable for higher
degree approximation. Here, we propose an approximation and prove its accuracy
which leads to a numerically stable way for approximating non-periodic uniform
data with quadratic convergence.

For a more detailed formulation we need some notation. Given an
$m$-dimensional function space $V_{m}\subset L_{2}$, we define the best
possible approximation (projection) to $f\colon D\to\mathds{C}$ in $V_{m}$ and
its error:

$\displaystyle P(f,V_{m},L_{p})=\operatorname*{arg\,min}_{g\in
V_{m}}\|f-g\|_{L_{p}}\text{ and
}e(f,V_{m},L_{p})_{L_{q}}=\|f-P(f,V_{m},L_{p})\|_{L_{q}}$

for $p,q\in\\{2,\infty\\}$. Note, since $V_{m}$ is finite-dimensional the
minimum is actually attained. Following Baraud02 ; CDL13 ; MNST14 ; CM17 ;
KUV21 ; LPU21 , we use weighted least squares $S_{m}$, defined in (3.1), as
underlying approximation method. Because of its linearity, the approximation
error $\|f-S_{m}\bm{y}\|_{L_{2}}$ splits as follows:

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$
$\displaystyle=e(f,V_{m},L_{2})_{L_{2}}^{2}+\|P(f,V_{m},L_{2})-S_{m}\bm{y}\|_{L_{2}}^{2}$
$\displaystyle\leq\underbrace{e(f,V_{m},L_{2})_{L_{2}}^{2}}_{\text{truncation
error}}+2\underbrace{\|P(f,V_{m},L_{2})-S_{m}\bm{f}\|_{L_{2}}^{2}}_{\text{discretization
error}}+2\underbrace{\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}}_{\text{noise
error}}\,.$

For fixed number of points $n$, we have a look at the behaviour with respect
to $m$, the dimension of the approximation space $V_{m}$. The truncation error
is the best possible benchmark and usually has polynomial decay $m^{-s}$ for
some rate $s\geq 1$ depending on $f$ and the choice of $V_{m}$. We show, that
the discretization error obeys the same decay as the truncation error. Thus,
in the noiseless case we obtain the best possible error up to a multiplicative
constant, cf. Theorem 3.2.

Including noise, we show that we get an additional summand growing linear in
$m$, cf. Theorems 1.1. This resambles the well-known bias-variance trade off
modeling the over- and undersmoothing effects which one wants to balance, cf.
GKKW02 ; PL13 . This linear behaviour in $m$ is approved by (LMP20, , Thm.
4.9) (by using the regularization $g_{\lambda}(\sigma)=1/(\lambda+\sigma)$
with $\lambda=0$). An example of that behaviour for $D=[0,1]$ and
$\mu=\mathrm{d}x$ being the Lebesgue measure is depicted in Figure 4.1 where
the detailed example is found in Section 4.

$0$$1$$0$$1$$0$$1$$0$$1$$0$$1$$0$$1$$m=2$(underfitting)$m=4$$m=16$(overfitting)$m$$\|f-S_{m}\bm{y}\|_{\mathrm{L}_{2}}^{2}$
Figure 1.1: One-dimensional approximation on the unit-interval for three
different choices of $m$ and the schematic behaviour of the
$L_{2}$-approximation error $\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$ (solid) split into
$\|f-S_{m}\bm{f}\|_{L_{2}}^{2}+\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}$ (dashed)
with respect to $m$.

Our central theorem, complying this behaviour, looks as follows:

######  Theorem 1.1.

Let $f\colon D\to\mathds{C}$, $\bm{x}^{1},\dots,\bm{x}^{n}$, $n\in\mathds{N}$
be points drawn according to a probability measure
$\mathrm{d}\nu=\varrho\,\mathrm{d}\mu$ and
$\bm{y}=\bm{f}+\bm{\varepsilon}=(f(\bm{x}^{1})+\varepsilon_{1},\dots,f(\bm{x}^{n})+\varepsilon_{n})^{\mathsf{T}}$
noisy function values where
$\bm{\varepsilon}=(\varepsilon_{1},\dots,\varepsilon_{n})^{\mathsf{T}}$ is a
vector of independent complex-valued mean-zero random variables satisfying
$\mathds{E}(\lvert\varepsilon_{i}\rvert^{2})\leq\sigma^{2}$ and
$\lvert\varepsilon_{i}\rvert\leq B$ for $i=1,\dots,n$. Let further, $t\geq 0$,
$V_{m}$ be an $m$-dimensional function space with an
$L_{2}(D,\mu)$-orthonormal basis $\eta_{0},\dots,\eta_{m-1}$ satisfying

$\displaystyle
10\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}(\log(m)+t)\leq
n\quad\text{with}\quad
N(V_{m},\cdot)=\sum_{k=0}^{m-1}\lvert\eta_{k}(\cdot)\rvert^{2}\,.$

Then, for $S_{m}$ the weighted least squares method defined in (3.1) with
$\omega_{i}=1/\varrho(\bm{x}^{i})$, we have with joint probability exceeding
$1-3\exp(-t)$:

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$ $\displaystyle\leq
14\left(e(f,V_{m},L_{2})_{L_{2}}+\sqrt{\frac{t}{n}}e(f,V_{m},L_{2})_{L_{\infty}}\right)^{2}$
$\displaystyle\quad+\|4/\varrho\|_{\infty}\left(\frac{m}{n}\left(14B\sqrt{t\sigma^{2}}+\sigma^{2}\right)+\frac{128B^{2}t}{n}\right)\,,$

where $L_{2}=L_{2}(D,\mu)$ and $L_{\infty}=L_{\infty}(D,\mu)$.

The first line of the bound corresponds to the truncation error and
discretization error, decaying in $m$, and the second line the error due to
noise, increasing in $m$, cf. Figure 1.1. This theorem extends to the
$L_{\infty}$ case, cf. Theorem 3.4. The estimation of the noise error is using
a Hanson-Wright concentration inequality, which can be found using different
assumptions. Thus, we can replace the noise model by general Bernstein
conditions, cf. Lemma 2.3, or sub-Gaussian noise, cf. RV13

The sampling measure $\nu(E)\coloneqq\int_{E}\varrho\;d\mu$, induced by the
probability distribution $\varrho$, may differ from the error measure $\mu$,
which is known as the _change of measure_ and has applications in domain
adaptation, cf. PY10 . Note, the choice $\varrho$ affects the maximal size of
$V_{m}$ in the assumption and the amplification of the noise in bound. There
are two extremal cases:

  1. (i)

Having $\varrho(\bm{x})=N(V_{m},\bm{x})/m$, as it was done in HD15 ; NJZ17 ;
CM17 ; KUV21 , we obtain the assumption

$\displaystyle
10\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}(\log(m)+t)=10m(\log(m)+t)\leq
n\,,$

which allows for the biggest choice of $m$ possible. But this spoils
$\|1/\varrho\|_{\infty}$ in the error bound when the Christoffel function
attains small values.

  2. (ii)

For domains $D$ with bounded measure, we may choose
$\varrho(\bm{x})=1/\mu(D)$, as it was done in CDL13 ; CM17 ; LPU21 . As all
weights $\omega_{i}=\mu(D)$, $S_{m}$ becomes the plain least squares method.
In this case, $\|1/\varrho\|_{\infty}$ is minimal and noise is amplified the
least. But this choice spoils the assumption on the choice of $m$ when the
Christoffel function $N(V_{m},\bm{x})$ attains big values. This effect is
controllable, for instance, when working with a bounded orthonormal system
(BOS) ($\|\eta_{k}\|_{\infty}\leq B$ for some $B>0$ and all $k$). Then

$\displaystyle N(V_{m})\leq\sum_{k=0}^{m-1}\|\eta_{k}\|_{\infty}^{2}\leq
mB^{2}$

and the assumption on the size of $V_{m}$ can be replaced by

$\displaystyle
10\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}(\log(m)+t)\leq
10\mu(D)Bm(\log(m)+t)\leq n\,.$

An interesting example, where these effects occur, is the approximation of
functions on the unit interval $D=[0,1]$ from samples given in uniformly
random points.

  * •

When using algebraic polynomials and the Lebesgue error measure
$\mu=\mathrm{d}x$ we have to choose $\varrho\equiv 1$ to obtain uniform random
points. Orthogonalizing algebraic polynomials with respect to the Lebesgue
measure, we obtain $\eta_{k}=P_{k}/\|P_{k}\|_{L_{2}((0,1),\mathrm{d}x)}$
Legendre polynomials for our approximation space $V_{m}$. Since
$\|P_{k}\|_{L_{2}((0,1),\mathrm{d}x)}^{2}=2k+1$ and $P_{k}(0)=1$, we have

$\displaystyle N(V_{m},0)=\sum_{k=0}^{m-1}\frac{\lvert
P_{k}(0)\rvert^{2}}{\|P_{k}\|_{L_{2}}^{2}}=\sum_{k=0}^{m-1}(2k+1)=m^{2}\,.$
(1.1)

Thus, this case falls into category (i) from above and spoils our choice of
$m\leq\sqrt{n}$, i.e., quadratic oversampling as in MNST14 .

  * •

When using algebraic polynomials and the Chebyshev error measure
$\mu=(1-(2x-1)^{2})^{-1/2}\;\mathrm{d}x$ we have to choose
$\varrho(x)=\frac{4}{\pi}\sqrt{1-(2x-1)^{2}}$ to obtain uniform random points.
Orthogonalizing algebraic polynomials with respect to the Chebyshev measure,
we obtain Chebyshev polynomials $\eta_{k}(x)=T_{k}(x)=\cos(k\arccos(2x-1))$
for our approximation space $V_{m}$. These are a BOS, but the distribution
$\varrho$ spoils both the assumption on $m$ and the error bound, since
$1/\varrho(x)$ diverges at the border (this effect can be circumvented using a
padding technique at the border as we discuss in Section 4).

  * •

In Section 4.1 we construct orthogonal functions with respect to the Sobolev
space inner product

$\displaystyle\langle f,g\rangle_{H^{s}(0,1)}=\langle
f,g\rangle_{L_{2}((0,1),\mathrm{d}x)}+\langle
f^{(s)},g^{(s)}\rangle_{L_{2}((0,1),\mathrm{d}x)}\,,$

for $s=1$ and $s=2$, which are orthogonal with respect to the
$L_{2}((0,1),\mathrm{d}x)$ inner product as well. For $s=1$, these functions
were already considered in WW08 ; IN08 ; Adcock10 ; AH11 ; DDP14 ; SNC16 ;
CKNS16 ; KMNN21  and for $s=2$ in (AIN12, , Section 3), where also higher
orders can be found.

We show that they form an BOS and, by (ii) above, this basis is then suitable
for approximation in uniform random points on $D=[0,1]$ using plain least
squares and only logarithmic oversampling. The $H^{2}$ basis is prone to
numerical errors. To overcome this, we propose a numerically stable
approximation and proof its accuracy.

As for the Structure of this paper, we start with some tools from probability
theory in Section 2. In Sections 3 we show error bounds for the weighted least
squares method. The construction of the $H^{1}$ and $H^{2}$ basis mentioned
above are found in Section 4 along with a comparison to the Legendre and
Chebyshev polynomials. To indicate the applicability of our error bounds and
the proposed basis, we conduct numerical experiments in one and five
dimensions.

##  2 Tools from probabilty theory

For the later analysis we need concentration inequalities starting with
Bernstein’s inequality, which is found in the standard literature, cf. (FR13,
, Corollary 7.31) or (SC08, , Theorem 6.12).

######  Theorem 2.1 (Bernstein).

Let $\xi_{1},\dots,\xi_{n}$ be independent real-valued mean-zero random
variables satisfying $\mathds{E}(\xi_{i}^{2})\leq\sigma^{2}$ and
$\|\xi_{i}\|_{\infty}\leq B$ for $i=1,\dots,n$ and real numbers $\sigma^{2}$
and $B$. Then

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\xi_{i}\leq\frac{2Bt}{3n}+\sqrt{\frac{2\sigma^{2}t}{n}}$

with probability exceeding $1-\exp(-t)$.

Bernstein’s inequality gives a concentration bound for the sum of independent
random variables. We need similar bounds for quadratic forms in random
vectors, which are known as Hanson-Wright inequalities. To formulate them, we
need to introduce the spectral norm and the Frobenius norm of a matrix
$\bm{A}\in\mathds{C}^{m\times n}$

$\displaystyle\|\bm{A}\|_{2\to
2}=\sqrt{\lambda_{\max}(\bm{A}^{\ast}\bm{A})}=\sigma_{\max}(\bm{A})\quad\text{and}\quad\|\bm{A}\|_{F}=\sqrt{\sum_{k=1}^{m}\sum_{i=1}^{n}\lvert
a_{k,i}\rvert^{2}}\,,$

where $\lambda_{\max}$ and $\sigma_{\max}$ denote the largest eigenvalues and
singular values, respectively. The following result is such an inequality with
a Bernstein condition on the random variables and was shown in (Bellec19, ,
Theorem 3).

######  Theorem 2.2 (Hanson-Wright).

Let $\bm{\xi}=(\xi_{1},\dots,\xi_{n})^{\mathsf{T}}$ be a vector of independent
mean-zero random variables such that for all integers $p\geq 1$

$\displaystyle\mathds{E}(\lvert\xi_{i}\rvert^{2p})\leq
p!B^{2p-2}\sigma_{i}^{2}/2$ (2.1)

for real numbers $B\geq 0$, $\sigma_{i}\geq 0$, and all $i=1,\dots,n$. Let
further $\bm{A}\in\mathds{C}^{n\times n}$ and
$m=\mathds{E}(\bm{\xi}^{\ast}\bm{A}\bm{\xi})$. Then
$\bm{D}_{\sigma}=\operatorname*{diag}(\sigma_{1},\dots,\sigma_{n})$

$\displaystyle\bm{\xi}^{\ast}\bm{A}\bm{\xi}-m\leq
256B^{2}\|\bm{A}\bm{D}_{\sigma}\|_{2\to
2}t+8\sqrt{3}B\|\bm{A}\bm{D}_{\sigma}\|_{F}\sqrt{t}$

with probability exceeding $1-\exp(-t)$.

The following is a special case of the above Hanson-Wright inequality for
Hermitian positive semi-definite matrices and random variables with known
variance $\mathds{E}(\lvert\xi_{i}\rvert^{2})$ and a uniform bound
$\|\xi_{i}\|_{\infty}$.

######  Corollary 2.3.

Let $\bm{\xi}=(\xi_{1},\dots,\xi_{n})^{\mathsf{T}}$ be a vector of independent
complex-valued mean-zero random variables satisfying
$\mathds{E}(\lvert\xi_{i}\rvert^{2})\leq\sigma^{2}$ and
$\|\xi_{i}\|_{\infty}\leq B$ for $i=1,\dots,n$. Then for all
$\bm{A}\in\mathds{C}^{m\times n}$

$\displaystyle\|\bm{A}\bm{\xi}\|_{2}^{2}\leq 128B^{2}\|\bm{A}\|_{2\to
2}^{2}t+(8\sqrt{3}B\sqrt{t\sigma^{2}}+\sigma^{2})\|\bm{A}\|_{F}^{2}$

with probability exceeding $1-\exp(-t)$.

Proof: Since
$\|\bm{A}\bm{\xi}\|_{2}^{2}=\bm{\xi}^{\ast}\bm{A}^{\ast}\bm{A}\bm{\xi}$ is a
quadratic form we want to apply Theorem 2.2 on $\bm{A}^{\ast}\bm{A}$. For that
we check the moment condition (2.1) on $\xi_{1}^{2},\dots,\xi_{n}^{2}$. For
$p=1$ it is fulfilled for constants $B/\sqrt{2}$ and
$(\sqrt{2}\sigma_{i})^{2}$. For $p\geq 2$, we have $p!\geq 2^{p-1}$ and obtain

$\displaystyle\mathds{E}(\lvert a_{i}\xi_{i}\rvert^{2p})$
$\displaystyle\leq\|\xi_{i}\|_{\infty}^{2p-2}\mathds{E}(\lvert\xi_{i}\rvert^{2})$
$\displaystyle\leq(B)^{2p-2}\sigma^{2}$ $\displaystyle\leq
p!\left(\frac{B}{\sqrt{2}}\right)^{2p-2}\frac{(\sqrt{2}\sigma)^{2}}{2}\,.$

Therefore, Theorem 2.2 is applicable.

It is left to estimate the expected value. Since $\xi_{1},\dots,\xi_{n}$ are
independent and have bounded variance, we obtain

$\displaystyle\mathds{E}\left(\|\bm{A}\bm{\xi}\|_{2}^{2}\right)$
$\displaystyle=\sum_{k=1}^{m}\sum_{i=1}^{n}\sum_{j=1}^{n}a_{i,k}\overline{a_{j,k}}\,\mathds{E}(\xi_{i}\overline{\xi_{j}})$
$\displaystyle=\sum_{k=1}^{m}\left(\sum_{i=1}^{n}\sum_{j\neq
i}a_{i,k}\overline{a_{j,k}}\,\mathds{E}(\xi_{i}\overline{\xi_{j}})\right)+\sum_{i=1}^{n}\lvert
a_{i,k}\rvert^{2}\mathds{E}(\lvert\xi_{i}\rvert^{2})$
$\displaystyle\leq\sigma^{2}\|\bm{A}\|_{F}^{2}\,.$

The following tool is a concentration bound on the maximal singular values of
random matrices which was shown in (Tr12, , Theorem 1.1).

######  Lemma 2.4 (Matrix Chernoff).

For a finite sequence $\bm{A}_{1},\dots,\bm{A}_{n}\in\mathds{C}^{m\times m}$
of independent, Hermitian, positive semi-definite random matrices satisfying
$\lambda_{\max}(\bm{A}_{i})\leq R$ almost surely it holds

$\displaystyle\mathds{P}\left(\lambda_{\min}\left(\sum_{i=1}^{n}\bm{A}_{i}\right)\leq(1-t)\mu_{\min}\right)$
$\displaystyle\leq{m}\exp\left(-\frac{\mu_{\min}}{R}(t+(1-t)\log(1-t))\right)$
$\displaystyle\leq{m}\exp\left(-\frac{\mu_{\min}t^{2}}{2R}\right)$

and

$\displaystyle\mathds{P}\left(\lambda_{\max}\left(\sum_{i=1}^{n}\bm{A}_{i}\right)\geq(1+t)\mu_{\max}\right)$
$\displaystyle\leq{m}\exp\left(-\frac{\mu_{\max}}{R}(-t+(1+t)\log(1+t))\right)$
$\displaystyle\leq{m}\exp\left(-\frac{\mu_{\max}t^{2}}{3R}\right)$

for $t\in[0,1]$ where
$\mu_{\min}\coloneqq\lambda_{\min}(\sum_{i=1}^{n}\mathds{E}(\bm{A}_{i}))$ and
$\mu_{\max}\coloneqq\lambda_{\max}(\sum_{i=1}^{n}\mathds{E}(\bm{A}_{i}))$.

Proof: The first estimates are provided by (Tr12, , Theorem 1.1). Based on the
Taylor expansion

$\displaystyle(1+t)\log(1+t)=t+\sum_{k=2}^{\infty}\frac{(-1)^{k}}{k(k-1)}t^{k}\,,$

which holds true for $t\in[-1,1]$, we further derive the inequalities

$\displaystyle t+(1-t)\log(1-t)$
$\displaystyle=\sum_{k=2}^{\infty}\frac{1}{k(k-1)}t^{k}\geq\frac{t^{2}}{2}$
and $\displaystyle-t+(1+t)\log(1+t)$
$\displaystyle=\sum_{k=2}^{\infty}\frac{(-1)^{k}}{k(k-1)}t^{k}\geq\frac{t^{2}}{2}-\frac{t^{3}}{6}\geq\frac{t^{2}}{3}$

for the range $t\in[0,1]$.

##  3 Error bounds for least squares

In this section we develop $L_{2}$\- and $L_{\infty}$-error bounds for the
least squares method. To this end we introduce some notation and the method
itself. Let $\eta_{0},\dots,\eta_{m-1}\colon D\to\mathds{C}$ be an
$L_{2}$-orthonormal basis of $V_{m}$,

$\displaystyle
N(V_{m},\bm{x})=\sum_{k=0}^{m-1}\lvert\eta_{k}(\bm{x})\rvert^{2}\quad\text{and}\quad
N(V_{m})=\sup_{\bm{x}\in D}N(V_{m},\bm{x})$

be the _Christoffel function_ and its supremum. For our approximation method
$S_{m}$ we use the _weighted least squares approximation_ depending on
$\eta_{0},\dots,\eta_{m-1}$ and $\bm{x}^{1},\dots,\bm{x}^{n}$:

$\displaystyle(S_{m}\bm{y})(\bm{x})=\sum_{k=0}^{m-1}\hat{g}_{k}\eta_{k}(\bm{x})\quad\text{with}\quad$
$\displaystyle\bm{\hat{g}}=\operatorname*{arg\,min}_{\bm{\hat{a}}\in\mathds{C}^{m}}\|\bm{L}\bm{\hat{a}}-\bm{y}\|_{\bm{W}}^{2}\,,$
$\displaystyle\bm{L}=\begin{bmatrix}\eta_{0}(\bm{x}^{1})&\dots&\eta_{m-1}(\bm{x}^{1})\\\
\vdots&\ddots&\vdots\\\
\eta_{0}(\bm{x}^{n})&\dots&\eta_{m-1}(\bm{x}^{n})\end{bmatrix}\in\mathds{C}^{n\times
m}$ $\displaystyle\text{, and}\quad\bm{W}=\begin{bmatrix}\omega_{1}\\\
&\ddots\\\ &&\omega_{n}\end{bmatrix}\in\mathds{[}0,\infty)^{n\times n}$ (3.1)

where
$\|\bm{L}\bm{\hat{a}}-\bm{y}\|_{\bm{W}}^{2}=(\bm{L}\bm{\hat{a}}-\bm{y})^{\ast}\bm{W}(\bm{L}\bm{\hat{a}}-\bm{y})$.
If all weights are equal we speak of _plain least squares approximation_.

The coefficients $\bm{\hat{g}}$ of the approximation $S_{m}\bm{y}$ are found
by solving the normal equation

$\displaystyle\bm{\hat{g}}=(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}\bm{y}\,.$

Doing this by the means of an iterative solver, the stability and the
iteration count for a desired precision are determined by the spectral
properties of the above matrix, cf. (Greenbaum97, , Theorem 3.1.1). However,
these are fully determined by the spectral properties of $\bm{W}^{1/2}\bm{L}$,
since for a singular value decomposition
$\bm{W}^{1/2}\bm{L}=\bm{U}\bm{\Sigma}\bm{V}^{\ast}$, we obtain

$\displaystyle(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}^{1/2}=\bm{V}(\bm{\Sigma}^{\ast}\bm{\Sigma})^{-1}\bm{\Sigma}\bm{U}^{\ast}\,.$
(3.2)

For $f=\sum_{k=0}^{m-1}\hat{f}_{k}\eta_{k}\in V_{m}$, the singular values of
$\bm{W}^{1/2}\bm{L}$ relate the coefficients
$\bm{\hat{f}}=(\hat{f}_{0},\dots,\hat{f}_{m-1})^{\mathsf{T}}$ with the samples
$\bm{f}=(f(\bm{x}^{1}),\dots,f(\bm{x}^{n}))^{\mathsf{T}}=\bm{L}\bm{\hat{f}}$.
Such connections can be understood as $L_{2}$-Marcinkiewicz-Zygmund inequality
for $V_{m}$. In the following lemma this is done for random points and was
used before in (NSU21, , Theorem 2.3), (MU21, , Theorem 5.1), or (BSU22, ,
Theorem 2.1).

######  Lemma 3.1.

Let $t\geq 0$, $n\in\mathds{N}$, $\bm{x}^{1},\dots,\bm{x}^{n}$ be points drawn
according to a probability measure $\mathrm{d}\nu=\varrho\,\mathrm{d}\mu$. Let
further, $V_{m}$ be an $m$-dimensional function space with an orthonormal
basis $\eta_{0},\dots,\eta_{m-1}$ in $L_{2}$ with $m$ satisfying

$\displaystyle
10\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}(\log(m)+t)\leq
n$

and $\bm{L},\bm{W}$ be as in (3.1) with $\omega_{i}=1/\varrho(\bm{x}^{i})$.
Then

$\displaystyle\frac{n}{2}\|\bm{\hat{g}}\|_{2}^{2}\leq\|\bm{W}^{1/2}\bm{L}\bm{\hat{g}}\|_{2}^{2}\leq\frac{3n}{2}\|\bm{\hat{g}}\|_{2}^{2}\quad\text{for
all}\quad\bm{\hat{g}}\in\mathds{C}^{m},$

where each inequality holds with probability exceeding $1-\exp(-t)$,
repectively.

Proof: The result is a direct consequence of Tropp’s result in Lemma 2.4. For
a randomly chosen point $\bm{x}^{i}$ we define the random rank-one matrix
$\bm{A}_{i}=\frac{1}{n}\varrho^{-1}(\bm{x}^{i})(\bm{y}^{i}\otimes\bm{y}^{i})$
with
$\bm{y}^{i}=(\eta_{0}(\bm{x}^{i}),\dots,\eta_{m-1}(\bm{x}^{i}))^{\mathsf{T}}$.
By construction, it holds

$\displaystyle\sum_{i=1}^{n}\bm{A}_{i}=\bm{L}^{\ast}\bm{W}\bm{L}$

and by the orthogonality of $\eta_{k}$

$\displaystyle\left(\mathds{E}(\bm{A}_{i})\right)_{k,l}=\frac{1}{n}\int_{D}\eta_{k}(\bm{x})\overline{\eta_{l}(\bm{x})}\varrho^{-1}(\bm{x})\varrho(\bm{x})\;d\mu(\bm{x})=\frac{\delta_{k,l}}{n}\,,$

which gives
$\mathds{E}\left(\sum_{i=1}^{n}\bm{A}_{i}\right)=\mathrm{Id}_{m\times m}$ and,
therefore, $\mu_{\max}=\mu_{\min}=1$. Further, we have

$\displaystyle\lambda_{\max}\left(\frac{1}{n}\varrho^{-1}(\bm{x}^{i})(\bm{y}^{i}\otimes\bm{y}^{i})\right)=\frac{1}{n}\varrho(\bm{x}^{i})^{-1}\|\bm{y}^{i}\|_{2}^{2}\leq\frac{1}{n}\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}\,.$

Lemma 2.4 with $t=1/2$ then gives the lower bound

$\displaystyle\mathds{P}\left(\lambda_{\min}\left(\frac{1}{n}\bm{L}^{\ast}\bm{W}\bm{L}\right)\leq\frac{1}{2}\right)$
$\displaystyle\leq{m}\exp\left(-\frac{n}{10}\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}^{-1}\right)\,,$

which is smaller than $\exp(-t)$ by the assumption on $m$.

The bound for the largest eigenvalue works analogue.

We now formulate a bound on the $L_{2}$-error of the weighted least squares
method for exact function values. This result is heavily based on (LPU21, ,
Theorem 3.20) which extends to a more general setting.

######  Theorem 3.2 ($L_{2}$-error bound without noise).

Let $f\colon D\to\mathds{C}$, $\bm{x}^{1},\dots,\bm{x}^{n}$, $n\in\mathds{N}$
be points drawn according to a probability measure
$\mathrm{d}\nu=\varrho\,\mathrm{d}\mu$ and
$\bm{y}=(f(\bm{x}^{1}),\dots,f(\bm{x}^{n}))^{\mathsf{T}}$ exact function
values. Let further, $t\geq 0$, $V_{m}$ be an $m$-dimensional function space
with an orthonormal basis $\eta_{0},\dots,\eta_{m-1}$ satisfying

$\displaystyle
10\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}(\log(m)+t)\leq
n\,.$

Then, for $S_{m}$ the weighted least squares method defined in (3.1) with
$\omega_{i}=1/\varrho(\bm{x}^{i})$, we have with probability exceeding
$1-2\exp(-t)$:

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}\leq
8\left(e(f,V_{m},L_{2})_{L_{2}}+\sqrt{\frac{t}{n}}e(f,V_{m},L_{2})_{L_{\infty}}\right)^{2}\,.$

Proof: For abbreviation, we use $e_{2}=e(f,V_{m},L_{2})_{L_{2}}$ and
$e_{\infty}=e(f,V_{m},L_{2})_{L_{\infty}}$. Further, we define the event

$\displaystyle A\coloneqq\left\\{\bm{x}^{1},\dots,\bm{x}^{n}\in
D:\frac{n}{2}\leq\|\bm{W}^{1/2}\bm{L}\|_{2\to 2}^{2}\right\\}$ (3.3)

which has probability $\mathds{P}(A)\geq 1-\exp(-t)$ by Lemma 3.1 and the
assumption on $V_{m}$. We split the approximation error

$\displaystyle\|f-S_{m}\bm{f}\|_{L_{2}}^{2}$
$\displaystyle=e_{2}^{2}+\|P(f,V_{m},L_{2})-S_{m}\bm{f}\|_{L_{2}}^{2}\,.$

Due to the invariance of $S_{m}$ to functions in $V_{m}$, we pull it in front
and use compatibility of the operator norm to obtain

$\displaystyle\|f-S_{m}\bm{f}\|_{L_{2}}^{2}$ $\displaystyle\leq
e_{2}^{2}+\|S_{m}\|_{2\to
2}^{2}\sum_{i=1}^{n}\frac{1}{\varrho(\bm{x}^{i})}\lvert(f-P(f,V_{m},L_{2}))(\bm{x}^{i})\rvert^{2}\,.$

By (3.2) and the event (3.3), we have $\|S_{m}\|_{2\to
2}^{2}=\|\bm{W}^{1/2}\bm{L}\|_{2\to 2}^{-1}\leq 2/n$. Thus,

$\displaystyle\|f-S_{m}\bm{f}\|_{L_{2}}^{2}$ $\displaystyle\leq
3e_{2}^{2}+\frac{2}{n}\sum_{i=1}^{n}\left\lvert\omega_{i}\lvert(f-P(f,V_{m},L_{2}))(\bm{x}^{i})\rvert^{2}-e_{2}^{2}\right\rvert\,.$

It remains to estimate the latter summand. We define

$\displaystyle\xi_{i}=\frac{1}{\varrho(\bm{x}^{i})}\left\lvert(f-P_{m}f)(\bm{x}^{i})\right\rvert^{2}-e_{2}^{2},$

which is mean-zero since we sample with respect to the distribution $\varrho$.
Further, we have

$\displaystyle\mathds{E}(\xi_{i}^{2})=\mathds{E}\left(\frac{\lvert(f-P_{m}f)(\bm{x}^{i})\rvert^{4}}{\varrho^{2}(\bm{x}^{1})}\right)-e_{2}^{4}\leq\|f-P_{m}f\|_{L_{\infty}}^{2}e_{2}^{2}-e_{2}^{4}\leq
e_{2}^{2}(e_{2}+e_{\infty})^{2}\,,$

and

$\displaystyle\|\xi_{i}\|_{\infty}\leq\sup_{x\in
D}\left\lvert\frac{\lvert(f-P_{m}f)(\bm{x})\rvert^{2}}{\varrho(\bm{x})}-e_{2}^{2}\right\rvert\leq
e_{\infty}^{2}+e_{2}^{2}\,.$

Thus, the conditions in order to apply Bernstein are fulfilled:

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\xi_{i}$
$\displaystyle\leq\frac{2t}{3n}\left(e_{2}^{2}+e_{\infty}^{2}\right)+\sqrt{\frac{2t}{n}}\left(e_{\infty}e_{2}+e_{2}^{2}\right)$
$\displaystyle\leq\left(\frac{2}{3}+\sqrt{2}\right)e_{2}^{2}+\sqrt{\frac{2t}{n}}e_{\infty}e_{2}+\frac{2t}{3n}e_{\infty}^{2}$
(3.4)

with probability $1-\exp(-t)$, where $t\leq n$ was used in the last
inequality. Thus,

$\displaystyle\|f-S_{m}\bm{f}\|_{L_{2}}^{2}$
$\displaystyle\leq\left(\frac{13}{3}+2\sqrt{2}\right)e_{2}^{2}+\sqrt{\frac{8t}{n}}e_{\infty}e_{2}+\frac{4t}{3n}e_{\infty}^{2}$
$\displaystyle\leq\left(\frac{13}{3}+2\sqrt{2}\right)\left(e_{2}+\sqrt{\frac{t}{n}}e_{\infty}\right)^{2}\,.$

By union bound we obtain the overall probability exceeding the sum of the
probabilities of events given by (3.3) and (3.4).

######  Remark 3.3.

For any $g=\sum_{k=1}^{m}\langle g,\eta_{k}\rangle_{L_{2}}\eta_{k}\in V_{m}$
the Hölder-inequality gives an estimate on the $L_{\infty}$-norm in terms of
the $L_{2}$-norm:

$\displaystyle\|g\|_{L_{\infty}}=\left\|\sum_{k=1}^{m}\langle
g,\eta_{k}\rangle_{L_{2}}\eta_{k}\right\|_{L_{\infty}}\leq\left\|\sqrt{\sum_{k=1}^{m}\lvert\langle
g,\eta_{k}\rangle_{L_{2}}\rvert^{2}}\sqrt{\sum_{k=1}^{m}\lvert\eta_{k}\rvert^{2}}\right\|_{L_{\infty}}=\sqrt{N(V_{m})}\|g\|_{L_{2}}\,.$
(3.5)

Using the assumption on $V_{m}$ from Theorem 3.2, we have

$\displaystyle\sqrt{\frac{t}{n}}e(f,V_{m},L_{2})_{L_{\infty}}$
$\displaystyle\leq\sqrt{\frac{tN(V_{m})}{n}}e(f,V_{m},L_{2})_{L_{2}}$
$\displaystyle\leq\sqrt{\frac{t}{10(\log(m)+t)}\frac{N(V_{m})}{\|N(V_{m},\cdot)/\varrho(\cdot)\|_{\infty}}}e(f,V_{m},L_{2})_{L_{2}}\,,$

which converges to zero whenever
$N(V_{m})/\|N(V_{m},\cdot)/\varrho(\cdot)\|_{\infty}$ is bounded. Thus,
asymptotically, Theorem 3.2 provides the best possible error
$e(f,V_{m},L_{2})_{L_{2}}$ up to a multiplicative constant of eight.

Next, we prove the central theorem which includes the noisy case.

Proof: [Proof of Theorem 1.1] We split the approximation error

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$ $\displaystyle\leq
e(f,V_{m},L_{2})_{L_{2}}+2\|f-S_{m}\bm{f}\|_{L_{2}}^{2}+2\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}$

and bound the first two summands as in the proof of Theorem 3.2 with the
events given by (3.3) and (3.4). Note, the constant changes from
$13/3+2\sqrt{2}$ to $23/3+4\sqrt{2}\leq 14$. Now, we focus on the third
summand. Applying Corollary 2.3 gives

$\displaystyle\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}=\|(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}\bm{\varepsilon}\|_{2}^{2}$
$\displaystyle\leq
128\|\bm{\varepsilon}\|_{\infty}^{2}\|(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}\|_{2\to
2}^{2}t+(8\sqrt{3}\|\bm{\varepsilon}\|_{\infty}\sqrt{t\sigma^{2}}+\sigma^{2})\|(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}\|_{F}^{2}$
(3.6)

with probability $1-\exp(-t)$. Since
$\bm{L}^{\ast}\bm{W}\bm{L}\in\mathds{C}^{m\times m}$, the matrix
$(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}\bm{W}^{1/2}$ has rank at most $m$ and,
thus, we use
$\|\bm{A}\|_{F}^{2}\leq\operatorname*{rank}(\bm{A})\|\bm{A}\|_{2\to 2}^{2}$ to
obtain

$\displaystyle\|(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}\|_{F}^{2}\leq\|1/\varrho\|_{\infty}m\|(\bm{L}^{\ast}\bm{W}\bm{L})^{-1}\bm{L}^{\ast}\bm{W}^{1/2}\|_{2\to
2}^{2}\leq\|1/\varrho\|_{\infty}\frac{2m}{n}$

where the last inequality follows from (3.2) and event (3.3). Therefore,

$\displaystyle\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}\leq\|1/\varrho\|_{\infty}\left(128\|\bm{\varepsilon}\|_{\infty}^{2}\frac{2}{n}t+(8\sqrt{3}\|\bm{\varepsilon}\|_{\infty}\sqrt{t\sigma^{2}}+\sigma^{2})\frac{2m}{n}\right)\,.$

By union bound we obtain the overall probability exceeding the sum of the
probabilities of the events given by (3.3), (3.4), and (3.6).

The following theorem bounds the approximation error in the $L_{\infty}$-norm.

######  Theorem 3.4 ($L_{\infty}$-error bound with noise).

Let the assumptions of Theorem 1.1 hold. Then, for $S_{m}$ the weighted least
squares method defined in (3.1) with $\omega_{i}=1/\varrho(\bm{x}^{i})$, we
have with probability exceeding $1-3\exp(-t)$:

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{\infty}}$
$\displaystyle\leq\left(1+\sqrt{5N(V_{m})}\right)\left(e(f,V_{m},L_{\infty})_{L_{\infty}}+\sqrt{\frac{t}{n}}e(f,V_{m},L_{\infty})_{L_{2}}\right)$
$\displaystyle\quad+\sqrt{\|2/\varrho\|_{\infty}N(V_{m})}\sqrt{\frac{m}{n}\left(14B\sqrt{t\sigma^{2}}+\sigma^{2}\right)+\frac{128B^{2}t}{n}}\,.$

Proof: For abbreviation, we use $e_{2}=e(f,V_{m},L_{\infty})_{L_{2}}$ and
$e_{\infty}=e(f,V_{m},L_{\infty})_{L_{\infty}}$. Using (3.5) we reduce the
$L_{\infty}$-case to the $L_{2}$-case which we already covered. We split the
approximation error

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{\infty}}$
$\displaystyle\leq\|f-P(f,V_{m},L_{\infty})\|_{L_{\infty}}+\|P(f,V_{m},L_{\infty})-S_{m}\bm{f}\|_{L_{\infty}}+\|S_{m}\bm{\varepsilon}\|_{L_{\infty}}$
$\displaystyle\leq
e_{\infty}+\sqrt{N(V_{m})}\|P(f,V_{m},L_{\infty})-S_{m}\bm{f}\|_{L_{2}}+\sqrt{N(V_{m})}\|S_{m}\bm{\varepsilon}\|_{L_{2}}\,.$

Analogously to (3.4) we obtain

$\displaystyle\|P(f,V_{m},L_{\infty})-S_{m}\bm{f}\|_{L_{2}}^{2}$
$\displaystyle\leq\left(\frac{2}{3}+\sqrt{2}\right)e_{\infty}^{2}+\sqrt{\frac{2t}{n}}e_{\infty}e_{2}+\frac{2t}{3n}e_{2}^{2}$
$\displaystyle\leq\left(\frac{2}{3}+\sqrt{2}\right)\left(e_{\infty}+\sqrt{\frac{t}{n}}e_{2}\right)^{2}\,,$

where the last inequality follows from $t\leq n$. Thus,

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{\infty}}$
$\displaystyle\leq\left(1+\sqrt{\frac{4+6\sqrt{2}}{3}N(V_{m})}\right)\left(e_{\infty}+\sqrt{\frac{t}{n}}e_{2}\right)+\sqrt{N(V_{m})}\|S_{m}\bm{\varepsilon}\|_{L_{2}}\,.$

Using the same bound as in Theorem 1.1 for $\|S_{m}\bm{\varepsilon}\|_{L_{2}}$
we obtain the assertion.

##  4 Application on the unit cube

In this section we are interested in approximating functions on the
$d$-dimensional unit cube $D=[0,1]^{d}$ when sample points are drawn
uniformly, i.e., with respect to the Lebesgue measure $\mathrm{d}\bm{x}$. The
deterministic equivalent to uniform sampling are equispaced points. When using
these for polynomial interpolation, Runge already knew in 1901, that higher
degree polynomials lead to oscillatory behaviour towards the border which
spoil the approximation error. Even though, we do not interpolate, we will
observe similar behaviour using Legendre and Chebyshev polynomials. We propose
alternative bases, which are stable for large $m=\dim(V_{m})$ as well.

Throughout this section we have $L_{2}=L_{2}((0,1)^{d},\mathrm{d}\bm{x})$
unless stated otherwise.

###  4.1 Sobolev spaces on the unit interval

Let $d=1$, $D=[0,1]$ be the unit interval equipped with the Lebesgue measure
$\mathrm{d}x$. In order to get hold on appropriate finite-dimensional function
spaces $V_{m}$ for approximation, we have a look at Sobolev spaces
$H^{s}=H^{s}(0,1)$ with integer smoothness $s\geq 0$. The inner product of
these Hilbert spaces is given by

$\displaystyle\langle f,g\rangle_{H^{s}}=\langle f,g\rangle_{L_{2}}+\langle
f^{(s)},g^{(s)}\rangle_{L_{2}}\,.$

Since $\|f\|_{L_{2}}^{2}\leq\|f\|_{H^{s}}^{2}=\langle f,f\rangle_{H^{s}}$, the
embedding operator $\mathrm{Id}\colon H^{s}\hookrightarrow L_{2}$ is compact.
Thus, $W=\mathrm{Id}^{\ast}\circ\mathrm{Id}\colon H^{s}\to H^{s}$ is compact
and self-adjoint. Applying the spectral theorem gives for $f\in H^{s}$

$\displaystyle W(f)=\sum_{k=0}^{\infty}\sigma_{k}\langle
f,e_{k}\rangle_{H^{s}}e_{k}$

where $(\sigma_{k})_{k=0}^{\infty}$ is the non-increasing rearrangement of the
singular values of $W$ and $(e_{k})_{k=0}^{\infty}\subset H^{s}$ the
corresponding system of eigenfunctions forming an orthonormal basis in
$H^{s}$. Since

$\displaystyle\langle
e_{k},e_{l}\rangle_{L_{2}}=\langle\mathrm{Id}(e_{k}),\mathrm{Id}(e_{l})\rangle_{L_{2}}=\langle
W(e_{k}),e_{l}\rangle_{H^{s}}=\sigma_{k}^{2}\langle
e_{k},e_{l}\rangle_{H^{s}}=\sigma_{k}^{2}\delta_{k,l}\,,$

the functions $\eta_{k}=\sigma_{k}^{-1}e_{k}$ form an orthonormal system in
$L_{2}$. Setting $V_{m}=\operatorname*{span}\\{\eta_{k}\\}_{k=0}^{m-1}$, we
obtain for $H^{s}\ni f=\sum_{k=0}^{\infty}\langle f,e_{k}\rangle_{H^{s}}e_{k}$

$\displaystyle e(f,V_{m},L_{2})_{L_{2}}^{2}=\left\|\sum_{k=m}^{\infty}\langle
f,e_{k}\rangle_{H^{s}}e_{k}\right\|_{L_{2}}^{2}=\sum_{k=m}^{\infty}\left\lvert\langle
f,e_{k}\rangle_{H^{s}}\sigma_{k}\right\rvert^{2}\leq\|f\|_{H^{s}}^{2}\sigma_{m}^{2}\,.$
(4.1)

Thus, the eigenfunctions corresponding to the largest singular values are a
canonical candidate for the approximation space $V_{m}$. To put this into
concrete terms, in the next two theorems, we give the singular values and
eigenfunctions for $H^{1}$ and $H^{2}$.

######  Theorem 4.1.

The operator $W=\mathrm{Id}^{\ast}\circ\mathrm{Id}\colon H^{1}\to H^{1}$ has
singular values $\sigma_{k}^{2}=\frac{1}{1+\pi^{2}k^{2}}$ with corresponding
$L_{2}$-normalized eigenfunctions

$\displaystyle\eta_{k}(x)=\begin{cases}1&\text{for }k=0\\\ \sqrt{2}\cos(\pi
kx)&\text{for }k\geq 1\,.\end{cases}$

Proof: For $\sigma$ a singular value of $W$ with corresponding eigenfunction
$\eta\in H^{1}$ and $\varphi\in H^{1}$ a test function, we have

$\displaystyle\langle\eta,\varphi\rangle_{L_{2}}=\langle\mathrm{Id}(\eta),\mathrm{Id}(\varphi)\rangle_{L_{2}}=\langle
W(\eta),\varphi\rangle_{H^{1}}=\sigma^{2}\langle\eta,\varphi\rangle_{H^{1}}=\sigma^{2}\left(\langle\eta,\varphi\rangle_{L_{2}}+\langle\eta^{\prime},\varphi^{\prime}\rangle_{L_{2}}\right)\,.$

Partial differentiation yields
$\langle\eta^{\prime},\varphi^{\prime}\rangle_{L_{2}}=\eta^{\prime}(1)\varphi(1)-\eta^{\prime}(0)\varphi(0)-\langle\eta^{\prime\prime},\varphi\rangle_{L_{2}}$.
Thus

$\displaystyle\left\langle\frac{1-\sigma^{2}}{\sigma^{2}}\eta+\eta^{\prime\prime},\varphi\right\rangle_{L_{2}}=\eta^{\prime}(1)\varphi(1)-\eta^{\prime}(0)\varphi(0)\,.$

Since this has to hold for all test functions $\varphi\in H^{1}$, we obtain
the differential equation

$\displaystyle\frac{1-\sigma^{2}}{\sigma^{2}}\eta=-\eta^{\prime\prime}\quad\text{with}\quad\eta(0)^{\prime}=\eta(1)^{\prime}=0\,.$

The proposed functions are exactly the ones fulfilling this differential
equation.

The $H^{1}$ basis above was already considered in (WW08, , Lemma 4.1) with the
same proof technique, in IN08  as a modified Fourier expansion, and in SNC16
in the context of samples along tent-transformed rank-1 lattices, and further
in Adcock10 ; AH11 ; DDP14 ; CKNS16 ; KMNN21 . The following $H^{2}$ basis was
already posed in (AIN12, , Section 3), where higher-order Sobolev-spaces are
found as well. The proof of the following theorem is found in Appendix A.

######  Theorem 4.2.

The operator $W=\mathrm{Id}^{\ast}\circ\mathrm{Id}\colon H^{2}\to H^{2}$ has
singular values $\sigma_{0}^{2}=1$ with corresponding $L_{2}$-normalized
eigenfunctions

$\displaystyle\eta_{0}(x)=1\quad\text{and}\quad\eta_{1}(x)=2\sqrt{3}x-\sqrt{3}$

and for $k\geq 2$, $\sigma_{k}^{2}=\frac{1}{1+t_{k}^{4}}$ with $t_{k}>0$ the
solutions of $\cosh(t_{k})\cos(t_{k})=1$ ($t_{k}\approx\frac{2k-1}{2}\pi$, cf.
Lemma A.1) and

$\displaystyle\eta_{k}(x)=\cosh(t_{k}x)+\cos(t_{k}x)-\displaystyle{\frac{\cosh(t_{k})-\cos(t_{k})}{\sinh(t_{k})-\sin(t_{k})}(\sinh(t_{k}x)+\sin(t_{k}x))}\,.$

Further, it holds

$\displaystyle\|\eta_{k}\|_{\infty}\leq\begin{cases}1&\text{for }k=0\\\
\sqrt{3}&\text{for }k=1\\\ \sqrt{6}&\text{for }k\geq 2\,.\end{cases}$

The singular values $\sigma_{k}$ for $H^{2}$ decay quadraticly in contrast to
linearly for $H^{1}$. Thus, approximating a twice differentiable function,
$m=\dim(V_{m})$ can be chosen smaller when using the $H^{2}$ basis whilst
achieving the same truncation error $e(f,V_{m},L_{2})_{L_{2}}$. Furthermore,
as noise enters with the factor $m/n$, cf. Theorem 1.1, this helps prevent
overfitting as well and leads to a smaller approximation error.

However, as $\cosh$ and $\sinh$ both grow exponentially, the representation of
the $H^{2}$ basis in Theorem 4.2 is prone to cancellations and, therefore,
numerical unstable. In the next theorem we pose an approximation which is
numerically stable.

######  Theorem 4.3.

For $0<t_{2}<t_{3}<\dots$ fulfilling $\cosh(t_{k})\cos(t_{k})=1$, let
$\eta_{k}$ be as in Theorem 4.2. Further, for $n\geq 2$, let
$\tilde{t}_{k}=\pi(2k-1)/2$ and

$\displaystyle\tilde{\eta}_{k}(x)=\sqrt{2}\cos\left(\tilde{t}_{k}x+\pi/4\right)$
$\displaystyle\quad+\mathds{1}_{[0,1/2]}(x)\exp\left(-\tilde{t}_{k}x\right)+\mathds{1}_{[1/2,1]}(x)(-1)^{k}\exp\left(-\tilde{t}_{k}(1-x)\right)\,.$

Then $\lvert\eta_{k}(x)-\tilde{\eta}_{k}(x)\rvert\leq\varepsilon$ for
$k\geq\frac{2}{\pi}\log(16/\varepsilon)+1$. In particular, the approximation
$\tilde{\eta}_{k}$ is exact up to machine precision $\varepsilon=10^{-16}$ for
$k\geq 27$.

For the proof see Appendix A. For the numerical experiments we use the exact
representation from Theorem 4.2 for $m<10$ and the approximation from Theorem
4.3 for $m\geq 10$.

###  4.2 Polynomial approximation on the unit interval

Next, we examine how the $H^{1}$ and $H^{2}$ bases compares to polynomial
approximation when points are distributed uniformly, i.e.,
$V_{m}=\Pi_{m}=\operatorname*{span}\\{1,x,\dots,x^{m-1}\\}$ and
$\nu(x)=\varrho(x)d\mu(x)=\mathrm{d}x$. Polynomial approximation results often
assume $f\in X^{s}$ with

$\displaystyle
X^{s}\coloneqq\\{f\colon[0,1]\to\mathds{C}:f,\dots,f^{(s-1)}\text{ absolute
continuous},\,f^{(s)}\in BV([0,1])\\}\,,$

where $BV([0,1])$ are all functions with bounded variation. This assumtion is
stronger than assuming $f\in H^{s}$ as the following remark shows.

######  Remark 4.4 ($X^{s}\hookrightarrow H^{s+1/2-\varepsilon}$).

For a rigorous investigation of the relation of $X^{s}$ and $H^{s}$, we need
to define the Besov space $B_{p,q}^{s}$ for $p=1$, $q=\infty$, and integer
smoothness $s$

$\displaystyle B_{1,\infty}^{s}\coloneqq\left\\{f\in L_{1}:\sup_{h\neq
0}\frac{\|\Delta_{h}^{2}f^{(s-1)}\|_{L_{1}}}{\lvert h\rvert}<\infty\right\\}$

with the finite difference $(\Delta_{h}f)(x)\coloneqq f(x+h)-f(x)$ and
$\Delta_{h}^{2}=\Delta_{h}\circ\Delta_{h}$, cf. (Triebel92, , Section 1.2.5).

For $f\in X^{s}$ the derivative $f^{(s)}$ is of bounded variation. Thus, also
the finite difference $\Delta_{h}^{2}f$ is of bounded variation. In
particular, $f^{(s)}\in L_{1}$ and, therefore, $f\in B_{1,\infty}^{s+1}$. By
(Triebel92, , (2.3.2/23)), we further have $B_{1,\infty}^{s+1}\hookrightarrow
B_{1,1}^{s+1-\varepsilon}$ for any $\varepsilon>0$. Thus,

$\displaystyle X^{s}\hookrightarrow B_{1,\infty}^{s+1}\hookrightarrow
B_{1,1}^{s+1-\varepsilon}\hookrightarrow
W_{2}^{s+1/2-\varepsilon}=H^{s+1/2-\varepsilon}\,,$

where the third embedding follows from the Sobolev inequality, cf. (Triebel83,
, (2.7.1/1)), and the Sobolev space for non-integer smoothness $s$ consists of
functiosn $f$ such that $\langle f,\eta_{k}\rangle\leq Ck^{-s}$ for some
constant $C<\infty$.

Assuming $f\in X^{s}$, we have a look into approximating with Legendre- and
Chebyshev polynomials:

  * •

The canonical choice for the target measure is $\mu=\mathrm{d}x$ and
$\varrho\equiv 1$ such that $\nu(x)=\varrho(x)d\mu(x)=\mathrm{d}x$.
Orthogonalizing the first $m-1$ monomials with respect to
$\mu(x)=\mathrm{d}x$, we obtain the Legendre polynomials $P_{k}$.

For the error of the projection, assuming $f\in H^{s}$, the following was
shown in (Wang21, , Theorem 3.5):

$\displaystyle
e(f,V_{m},L_{2})_{L_{2}}^{2}\leq\frac{2V}{\pi(s+1/2)(m-s)^{2s+1}}\,,$ (4.2)

where $V$ is the total variation of $f^{(s)}$. With this stronger assumption
$f\in X^{s}$ half an order is gained by polynomial approximation in contrast
to the $H^{s}$ bases. This is expected as we require half an order of
smoothness in $L_{2}$ more, cf. Remark 4.4. (In the later numerical
experiments, we observe the gain of half an order for $H^{1}$ and $H^{2}$ as
well.)

A drawback comes with the Christoffel function $N(V_{m},\cdot)$. Since
$N(V_{m},0)=m^{2}$, cf. (1.1), this spoils the choice of $m$ to quadratic
oversampling:

$\displaystyle 10m^{2}(\log(m)+t)\leq n\,,$

which is usual for polynomial approximation in uniform points, cf. MNST14 .

  * •

When using the Chebyshev measure $\mu(x)=(1-(2x-1)^{2})^{-1/2}\;\mathrm{d}x$
we have to compensate with $\varrho(x)=\sqrt{1-(2x-1)^{2}}$ to obtain uniform
random samples as well. Orthogonalizing the first $m-1$ monomials with respect
to the Chebyshev measure, we obtain the Chebyshev polynomials $T_{k}$, which
are a BOS.

As for the error, assuming $f\in X^{s}$, we use (Trefethen13, , Theorem 7.1)
or (PPST18, , Theorem 6.16) to obtain the same bound as for Legendre
polynomials (4.2) but with respect to the
$L_{2}((0,1),(1-(2x-1)^{2})^{-1/2}\mathrm{d}x)$ norm.

As $\|1/\varrho\|_{\infty}$ diverges at the border, this spoils the choice of
the polynomial degree $m$ and our bound. Note, when we exclude some area
around the border for sampling, it does not diverge and the resulting error
can be controlled. This is called padding and was done in (PS22, , Section
4.1.2)

Thus, with polynomial approximation we assume half an order of smoothness
more, cf. Remark 4.4, which we also see in the approximation rate
$\mathcal{O}(m^{s+1/2})$.

######  Remark 4.5.

Note, when using the Chebyshev polynomials and samples with respect to the
Chebyshev measure, we have $\varrho\equiv 1$. Since the Chebyshev polynomials
are a BOS, this does not spoil our bounds.

Furthermore, using the Legendre polynomials ($\mu=\mathrm{d}x$) and samples
with respect to the Chebyshev measure ($\varrho=1/\pi(1-(2x-1)^{2})^{-1/2}$)
works as well. To see this we use (RW12, , Lemma 5.1):

$\displaystyle\sqrt{1-(2x-1)^{2}}\lvert
P_{k}(x)\rvert^{2}\leq\frac{2}{\pi}\left(2+\frac{1}{k}\right)$

for $k\geq 1$. Thus, $\|N(V_{m},\cdot)/\varrho(\cdot)\|_{\infty}$ and
$\|1/\varrho(\cdot)\|_{\infty}$ are bounded and do not spoil the choice of
polynomial degree $m$ nor the error bound.

###  4.3 Numerics on the unit interval

To support our findings, we give a numerical example. As a test function we
use

$\displaystyle f(x)=B_{2}^{\text{cut}}(x)\quad\text{with}\quad
B_{2}^{\text{cut}}(x)=\begin{cases}-x^{2}+3/4&\text{for }x\in[0,1/2]\\\
x^{2}/2-3/2x+9/8&\text{for }x\in[1/2,1]\\\ \end{cases}$ (4.3)

which was already considered in PV15 ; NP21 . The function
$B_{2}^{\text{cut}}$ is shown in Figure 1.1 and is a cutout of the B-spline of
order two. It and its first derivative are absolute continuous and the second
derivative is of bounded variation. Therefore $f\in X^{3}$ and the polynomial
approximation bounds from above are applicable. According to Remark 4.4 we
further have $f\in H^{5/2-\varepsilon}$ for any $\varepsilon>0$, i.e., there
exists $C\geq 0$ such that for $k\geq 0$ it holds $\langle
f,\eta_{k}\rangle_{L_{2}}\leq Ck^{-5/2+\varepsilon}$. In particular, $f\in
H^{2}$ and (4.1) is applicable for approximating with the $H^{1}$ and $H^{2}$
bases.

We sample $f$ in $10\,000$ uniformly random nodes and add $0.1\%M$ Gaussian
noise to obtain $\bm{y}=\bm{f}+\bm{\varepsilon}$, where
$M=\max_{x\in[0,1]}f(x)-\min_{x\in[0,1]}f(x)=5/8$. For $V_{m}$ we consider the
four choices from above: The Chebyshev polynomials
$V_{m}=\operatorname*{span}\\{T_{k}\\}_{k=0}^{m-1}$
($\mu=(1-(2x-1)^{2})^{-1/2}\;\mathrm{d}x$ and $\varrho=2/\pi$); the Legendre
polynomials $V_{m}=\operatorname*{span}\\{P_{k}\\}_{k=0}^{m-1}$
($\mu=\mathrm{d}x$ and $\varrho\equiv 1$); the first $m$ $H^{1}$ functions
from Theorem 4.1, and the the first $m$ $H^{2}$ functions from Theorems 4.2
and 4.3 ($\mu=\mathrm{d}x$ and $\varrho\equiv 1$ as well);. For
$m=\dim(V_{m})$ up to $1\,000$ we do the following:

  1. (i)

Compute the minimal and maximal singular values of
$1/\sqrt{n}\bm{W}^{1/2}\bm{L}$, with $\bm{W}$ and $\bm{L}$ given in (3.1).

  2. (ii)

We use least squares with $20$ iterations to obtain the approximation
$S_{m}\bm{y}=\sum_{k=0}^{m-1}\hat{g}_{k}\eta_{k}$, defined in (3.1).

  3. (iii)

We compute the $L_{2}$-error by using Parseval’s equality:

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}=\|f\|_{L_{2}}^{2}-\sum_{k=0}^{m-1}\lvert\hat{f}_{k}\rvert^{2}+\sum_{k=0}^{m-1}\lvert\hat{f}_{k}-\hat{g}_{k}\rvert^{2}\,,$

where the coefficients $\hat{f}_{k}=\langle f,\eta_{k}\rangle_{L_{2}}$ can be
computed analyticaly.

  4. (iv)

We compute the split approximation error:

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}\leq
2\|f-S_{m}\bm{f}\|_{L_{2}}^{2}+2\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}\,,$

where we compute both quantities separately, again, using Parseval’s equality.

The results are depicted in Figure 4.1.

$10^{-2}$$10^{0}$$10^{1}$$10^{2}$$10^{3}$$10^{-2}$$10^{0}$$10^{1}$$10^{2}$$10^{3}$$10^{-2}$$10^{0}$$10^{1}$$10^{2}$$10^{3}$$10^{-2}$$10^{0}$$10^{1}$$10^{2}$$10^{3}$$10^{-10}$$10^{-6}$$10^{-2}$$10^{1}$$10^{2}$$10^{3}$$10^{-10}$$10^{-6}$$10^{-2}$$10^{1}$$10^{2}$$10^{3}$$10^{-10}$$10^{-6}$$10^{-2}$$10^{1}$$10^{2}$$10^{3}$$10^{-10}$$10^{-6}$$10^{-2}$$10^{1}$$10^{2}$$10^{3}$$m$Chebyshev$m$Legendre$m$$\mathrm{H}^{1}$basis$m$$\mathrm{H}^{2}$basis$m$$m$$m$$m$
Figure 4.1: One-dimensional experiment for different choices of $V_{m}$. Top
row: minimal and maximal singular value of $1/\sqrt{n}\bm{W}^{1/2}\bm{L}$.
Bottom row: Solid: $\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$. Dashed:
$e(f,V_{m},L_{2})^{2}+\|P(f,V_{m},L_{2})-S_{m}\bm{f}\|_{L_{2}}^{2}$ and
$\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}$.

  * •

The smallest singular values for the Chebyshev polynomials and the Legendre
polynomials decay rapidly for bigger $m$. This coincides with the violation of
the assumtion in Lemma 3.1 for small $m$:

$\displaystyle
10\left\|\frac{N(V_{m},\cdot)}{\varrho(\cdot)}\right\|_{\infty}(\log(m)+t)\leq
n,$

where $\|N(V_{m},\cdot)/\varrho(\cdot)\|_{\infty}$ is unbounded in the
Chebyshev case and quadratic in the Legendre case, cf. (1.1). In this
experiment, for $m=1\,000$ the condition number
$\sigma_{\max}(\bm{W}^{1/2}\bm{L})/\sigma_{\min}(\bm{W}^{1/2}\bm{L})$ exceeded
$10^{29}$ for the algebraic polynomials and was below 14 for the $H^{s}$
basis.

  * •

The error for exact function values $\|f-S_{m}\bm{f}\|_{L_{2}}^{2}$ has decay
$3/2$ for $H^{1}$ and $5/2$ for the other bases. This conforms with the theory
for the polynomial bases. For the $H^{1}$ and $H^{2}$ bases the theory
predicted only decay rate $1$ and $2$, cf. Theorems 4.1, 4.2, and (4.1).

  * •

For the noise error $\|S_{m}\bm{\varepsilon}\|_{L_{2}}^{2}$ we observe linear
growth in $m=\dim(V_{m})$ as predicted in Theorem 1.1. Furthermore, this error
is bigger by a factor of around $40$ in the Chebyshev case compared to the
others. The maximal weight $\|\bm{W}\|_{\infty}$ in this case is around $40$
as well. The error due to noise in our bound has the factor
$\|1/\varrho\|_{\infty}$ which can be replaced by $\|\bm{W}\|_{\infty}$ to
sharpen the bound and explain this effect.

This numerical experiment and the earlier theoretical discussion shows, that
the $H^{1}$ and the $H^{2}$ bases are suitable for approximating functions on
the unit interval given in uniform random samples. They are numerically stable
in contrast to polynomial approximation with Chebyshev or Legendre. In
particular, the least squares matrix is well-conditioned and we can limit the
iterations when using an iterative solver, cf. (Greenbaum97, , Theorem 3.1.1).

###  4.4 Sobolev spaces with dominating mixed smoothness on the unit cube

The ideas from Subsection 4.1 can be extended to higher dimensions using the
concept of dominating mixed smoothness. We focus on the case of $H^{1}$ and
$H^{2}$, but the same can be done for polynomials as well, cf. (STW11, ,
Section 8.5.1).

Let $D=[0,1]^{d}$ be the $d$-dimensional unit cube equipped with the Lebesgue
measure $\mathrm{d}\bm{x}$. The Sobolev space with dominating mixed smoothness
of integer degree $s\geq 0$ is given by
$H_{\mathrm{mix}}^{s}=H_{\mathrm{mix}}^{s}(0,1)^{d}=H^{s}(0,1)\otimes\dots\otimes
H^{s}(0,1)$. The inner product of these Hilbert spaces is given by

$\displaystyle\langle
f,g\rangle_{H_{\mathrm{mix}}^{s}}=\sum_{\bm{j}\in\\{0,s\\}^{d}}\langle
D^{(\bm{j})}f,D^{(\bm{j})}g\rangle_{L_{2}}\,.$

With $\sigma_{k}$ and $\eta_{k}$ the singular values and eigenfunctions of
$H^{s}$, the singular values and eigenfunctions of
$W=\mathrm{Id}^{\ast}\circ\mathrm{Id}\colon H_{\mathrm{mix}}^{s}\to
H_{\mathrm{mix}}^{s}$ extend as follows:

$\displaystyle\sigma_{\bm{k}}^{2}=\prod_{j=1}^{d}\sigma_{k_{j}}^{2}\quad\text{and}\quad\eta_{\bm{k}}(\bm{x})=\prod_{j=1}^{d}\eta_{k_{j}}(x_{j})\,.$

To obtain the eigenfunctions corresponding to the smallest singular values, we
now work with multi-indices $\bm{k}$. The indices corresponding to the largest
singular values lie on a, so called, hyperbolic cross

$\displaystyle
I_{R}(H_{\mathrm{mix}}^{s})\coloneqq\left\\{\bm{k}\in\mathds{N}^{d}:\prod_{j=1}^{d}\sigma_{k_{j}}^{2}\geq
R\right\\}\,.$

For $V_{m}=\operatorname*{span}\\{\eta_{\bm{k}}:\bm{k}\in
I_{R}(H_{\mathrm{mix}}^{s})\\}$ and $f\in H_{\mathrm{mix}}^{s}$, we obtain by
(4.1)

$\displaystyle e(f,V_{m},L_{2})_{L_{2}}\leq
R\|f\|_{H_{\mathrm{mix}}^{s}}^{2}\,.$

In Figure 4.2 we have equally sized index sets for $H_{\mathrm{mix}}^{1}$ and
$H_{\mathrm{mix}}^{2}$. Note, that $R$ is smaller for $H_{\mathrm{mix}}^{2}$
as its singular values decay faster, cf. Theorems 4.1 and 4.2.

$I_{R}(H_{\mathrm{mix}}^{1})$ with $254$ frequencies and $R=5.3\cdot 10^{-5}$
$I_{R}(H_{\mathrm{mix}}^{2})$ with $254$ frequencies and $R=8.3\cdot 10^{-8}$
Figure 4.2: Hyperbolic cross in three dimensions.

###  4.5 Numerics on the unit cube

For a numerical experiment we do the same as in the one-dimensional case but
only consider the $H_{\mathrm{mix}}^{2}$ case. For our test function we
tensorize the B-Spline cutout

$\displaystyle f(\bm{x})=\prod_{j=1}^{d}B_{2}^{\text{cut}}(x_{j})$

where $B_{2}^{\text{cut}}$ was defined in (4.3). We increase the dimension to
$d=5$ and the number of samples to $1\,000\,000$ and use Gaussian noise with
variance $\sigma^{2}\in\\{0.00,\,0.01M,\,0.03M\\}$ where
$M=\max_{\bm{x}\in[0,1]^{d}}f(\bm{x})-\min_{\bm{x}\in[0,1]^{d}}f(\bm{x})=5/8$.

Let $V_{m}=\operatorname*{span}\\{\eta_{\bm{k}}:\bm{k}\in
I_{R}(H_{\mathrm{mix}}^{2})\\}$ of size $m$ with $\eta_{\bm{k}}$ the
tensorized $H_{\mathrm{mix}}^{2}$ basis, cf. Theorems 4.2 and 4.3. Since the
$H_{\mathrm{mix}}^{2}$ basis is a BOS, we obtain

$\displaystyle\frac{N(V_{m})}{m}\leq 6\,.$

With $t=6$, we satisfy the assumptions of Theorem 1.1 for $m\leq 12\,250$ and
obtain a probability exceeding $0.99$ for the error bound in Theorem 1.1. For
$m=\dim(V_{m})$ up to $10\,000$ we do the following:

  1. (i)

We use plain least squares with $20$ iterations to obtain the approximation
$S_{m}\bm{y}=\sum_{k=0}^{m-1}\hat{g}_{k}\eta_{k}$, defined in (3.1).

  2. (ii)

We compute the $L_{2}$-error by using Parseval’s equality analog to the one-
dimensional case.

  3. (iii)

We compute our bound: Applying Theorem 1.1 using (3.5), $t=6$, and
$n=1\,000\,000$, we obtain

$\displaystyle\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$ $\displaystyle\leq
14\left(1+\sqrt{\frac{6N(V_{m})}{n}}\right)^{2}e(f,V_{m},L_{2})_{L_{2}}^{2}$
$\displaystyle\quad+\frac{m}{n}\left(138B\sqrt{\sigma^{2}}+4\sigma^{2}\right)+0.0031B^{2}$

with probability exceeding $0.99$ and all the remaining quantities known in
our experiment.

$10^{-6}$$10^{-4}$$10^{-2}$$10^{0}$$10^{2}$$10^{3}$$10^{4}$$10^{-4}$$10^{-2}$$10^{0}$$10^{2}$$10^{3}$$10^{4}$$10^{-4}$$10^{-2}$$10^{0}$$10^{2}$$10^{3}$$10^{4}$$m$$\sigma^{2}=0.00$$m$$\sigma^{2}=0.01M$$m$$\sigma^{2}=0.03M$
Figure 4.3: Five-dimensional experiment for $H_{\mathrm{mix}}^{2}$. The solid
lines represent the $L_{2}$-error $\|f-S_{m}\bm{y}\|_{L_{2}}^{2}$ and the
dashed lines the bound from Theorems 3.2 and1.1.

The results are depicted in Figure 4.3.

  * •

The bounds capture the error behaviour well. But it seems that there is room
for improvement in the constants, especially in the experiments with noise.
Here, improving constants in the Hanson-Wright inequality in Theorem 2.3 could
be a starting point.

  * •

Furthermor, this experiment shows, that the $H^{2}$ basis is easily suitable
for high-dimensional approximation as well.

## Acknowledgement

I would like to thank Tino Ullrich for many helpful discussions and his
expertise with function spaces. Further, I would like to thank Michael
Schmischke for several insights on approximation on the unit cube, Sergei
Pereverzyev for his knowledge on regularization theory as well as Ralf
Hielscher and Daniel Potts for further discussions. Next to that I thank the
Deutscher Akademischer Austauschdienst (DAAD) for funding my research
scholarship.

## References

  *   * (1) Hampton, J., Doostan, A.: Coherence motivated sampling and convergence analysis of least squares polynomial Chaos regression. Comput. Methods Appl. Mech. Engrg. 290, 73–97 (2015). https://doi.org/10.1016/j.cma.2015.02.006
  * (2) Haji-Ali, A.-L., Nobile, F., Tempone, R., Wolfers, S.: Multilevel weighted least squares polynomial approximation. ESAIM Math. Model. Numer. Anal. 54(2), 649–677 (2020). https://doi.org/10.1051/m2an/2019045
  * (3) Baraud, Y.: Model selection for regression on a random design. ESAIM Probab. Stat. 6, 127–146 (2002). https://doi.org/10.1051/ps:2002007
  * (4) Cohen, A., Davenport, M.A., Leviatan, D.: On the stability and accuracy of least squares approximations. Found. Comput. Math. 13(5), 819–834 (2013). https://doi.org/10.1007/s10208-013-9142-3
  * (5) Cohen, A., Migliorati, G.: Optimal weighted least-squares methods. SMAI J. Comput. Math. 3, 181–203 (2017). https://doi.org/10.5802/smai-jcm.24
  * (6) Haberstich, C., Nouy, A., Perrin, G.: Boosted optimal weighted least-squares. Math. Comp. (2022). https://doi.org/10.1090/mcom/3710
  * (7) Kämmerer, L., Ullrich, T., Volkmer, T.: Worst-case recovery guarantees for least squares approximation using random samples. Constr. Approx. 54(2), 295–352 (2021). https://doi.org/10.1007/s00365-021-09555-0
  * (8) Migliorati, G., Nobile, F., von Schwerin, E., Tempone, R.: Analysis of discrete ${L}^{2}$ projection on polynomial spaces with random evaluations. Found. Comput. Math. (2014). https://doi.org/10.1007/s10208-013-9186-4
  * (9) Lippert, L., Potts, D., Ullrich, T.: Fast hyperbolic wavelet regression meets ANOVA. arXiv:2108.13197 (2021). https://doi.org/10.48550/arXiv.2108.13197
  * (10) Dolbeault, M., Cohen, A.: Optimal pointwise sampling for ${L}^{2}$ approximation. Journal of Complexity 68, 101602 (2022). https://doi.org/10.1016/j.jco.2021.101602
  * (11) Lu, S., Mathé, P., Pereverzev, S.V.: Balancing principle in supervised learning for a general regularization scheme. Appl. Comput. Harmon. Anal. 48(1), 123–148 (2020). https://doi.org/10.1016/j.acha.2018.03.001
  * (12) Pereverzyev, S.V., Lu, S.: Regularization Theory for Ill-Posed Problems. Selected Topics, p. 287. DeGruyter, Berlin, Boston (2013). https://doi.org/10.1515/9783110286496
  * (13) Bartel, F., Hielscher, R., Potts, D.: Fast cross-validation in harmonic approximation. Appl. Comput. Harmon. Anal. 49(2), 415–437 (2020). https://doi.org/10.1016/j.acha.2020.05.002
  * (14) Bartel, F., Hielscher, R.: Concentration inequalities for cross-validation in scattered data approximation. J. Approx. Theory 277, 105715 (2022). https://doi.org/10.1016/j.jat.2022.105715
  * (15) Temlyakov, V.: On optimal recovery in ${L}_{2}$. Journal of Complexity 65, 101545 (2021). https://doi.org/10.1016/j.jco.2020.101545
  * (16) Werschulz, A.G., Woźniakowski, H.: Tractability of multivariate approximation over a weighted unanchored sobolev space. Constr. Approx. 30(3), 395–421 (2009). https://doi.org/10.1007/s00365-009-9066-y
  * (17) Iserles, A., Nørsett, S.P.: From high oscillation to rapid approximation. I. Modified Fourier expansions. IMA J. Numer. Anal. 28(4), 862–887 (2008). https://doi.org/10.1093/imanum/drn006
  * (18) Adcock, B.: Multivariate modified Fourier series and application to boundary value problems. Numer. Math. 115(4), 511–552 (2010). https://doi.org/10.1007/s00211-010-0287-6
  * (19) Adcock, B., Huybrechs, D.: Multivariate modified Fourier expansions. In: Spectral and High Order Methods for Partial Differential Equations. Lect. Notes Comput. Sci. Eng., vol. 76, pp. 85–92. Springer, Heidelberg (2011). https://doi.org/10.1007/978-3-642-15337-2_5. https://doi.org/10.1007/978-3-642-15337-2_5
  * (20) Dick, J., Nuyens, D., Pillichshammer, F.: Lattice rules for nonperiodic smooth integrands. Numer. Math. 126(2), 259–291 (2014). https://doi.org/10.1007/s00211-013-0566-0
  * (21) Suryanarayana, G., Nuyens, D., Cools, R.: Reconstruction and collocation of a class of non-periodic functions by sampling along tent-transformed rank-1 lattices. J. Fourier Anal. Appl. 22(1), 187–214 (2016). https://doi.org/10.1007/s00041-015-9412-3
  * (22) Cools, R., Kuo, F.Y., Nuyens, D., Suryanarayana, G.: Tent-transformed lattice rules for integration and approximation of multivariate non-periodic functions. J. Complexity 36, 166–181 (2016). https://doi.org/10.1016/j.jco.2016.05.004
  * (23) Kuo, F.Y., Migliorati, G., Nobile, F., Nuyens, D.: Function integration, reconstruction and approximation using rank-1 lattices. Math. Comp. 90(330), 1861–1897 (2021). https://doi.org/10.1090/mcom/3595
  * (24) Adcock, B., Iserles, A., Nørsett, S.P.: From high oscillation to rapid approximation II: expansions in Birkhoff series. IMA J. Numer. Anal. 32(1), 105–140 (2012). https://doi.org/10.1093/imanum/drq038
  * (25) Györfi, L., Kohler, M., Krzyżak, A., Walk, H.: A Distribution-free Theory of Nonparametric Regression. Springer Series in Statistics, p. 647. Springer, New York (2002). https://doi.org/10.1007/b97848. https://doi.org/10.1007/b97848
  * (26) Rudelson, M., Vershynin, R.: Hanson-wright inequality and sub-gaussian concentration. Electronic Communications in Probability 18(none) (2013). https://doi.org/10.1214/ecp.v18-2865
  * (27) Pan, S.J., Yang, Q.: A survey on transfer learning. IEEE Transactions on Knowledge and Data Engineering 22(10), 1345–1359 (2010). https://doi.org/10.1109/TKDE.2009.191
  * (28) Narayan, A., Jakeman, J.D., Zhou, T.: A Christoffel function weighted least squares algorithm for collocation approximations. Math. Comp. 86(306), 1913–1947 (2017). https://doi.org/10.1090/mcom/3192
  * (29) Foucart, S., Rauhut, H.: A Mathematical Introduction to Compressive Sensing. Applied and Numerical Harmonic Analysis, p. 625. Birkhäuser/Springer, New York (2013). https://doi.org/10.1007/978-0-8176-4948-7
  * (30) Steinwart, I., Christmann, A.: Support Vector Machines. Information Science and Statistics, p. 601. Springer, New York (2008)
  * (31) Bellec, P.C.: Concentration of quadratic forms under a Bernstein moment assumption. arXiv:1901.08736 (2019) 1901.08736
  * (32) Tropp, J.A.: User-friendly tail bounds for sums of random matrices. Found. Comput. Math. 12(4), 389–434 (2012). https://doi.org/10.1007/s10208-011-9099-z
  * (33) Greenbaum, A.: Iterative Methods for Solving Linear Systems. Frontiers in Applied Mathematics, vol. 17, p. 220. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA (1997). https://doi.org/10.1137/1.9781611970937
  * (34) Nagel, N., Schäfer, M., Ullrich, T.: A new upper bound for sampling numbers. Found. Comp. Math. (2021). https://doi.org/10.1007/s10208-021-09504-0
  * (35) Moeller, M., Ullrich, T.: $L_{2}$-norm sampling discretization and recovery of functions from RKHS with finite trace. Sampl. Theory Signal Process. Data Anal. 19(2), 13 (2021). https://doi.org/10.1007/s43670-021-00013-3
  * (36) Bartel, F., Schäfer, M., Ullrich, T.: Constructive subsampling of finite frames with applications in optimal function recovery. arXiv:2202.12625 (2022) 2202.12625
  * (37) Triebel, H.: Theory of Function Spaces. II. Monographs in Mathematics, vol. 84, p. 370. Birkhäuser Verlag, Basel (1992). https://doi.org/10.1007/978-3-0346-0419-2
  * (38) Triebel, H.: Theory of Function Spaces. Modern Birkhäuser Classics, p. 285. Birkhäuser Verlag, Basel (2010)
  * (39) Wang, H.: New error bounds for Legendre approximations of differentiable functions. arXiv:2111.03833 (2021). https://doi.org/10.48550/arXiv.2111.03833
  * (40) Trefethen, L.N.: Approximation Theory and Approximation Practice, p. 305. Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA (2013). https://doi.org/10.1137/1.9781611975949
  * (41) Plonka, G., Potts, D., Steidl, G., Tasche, M.: Numerical Fourier Analysis. Applied and Numerical Harmonic Analysis, p. 168. Birkhäuser/Springer, Cham (2018). https://doi.org/10.1007/978-3-030-04306-3
  * (42) Potts, D., Schmischke, M.: Learning multivariate functions with low-dimensional structures using polynomial bases. J. Comput. Appl. Math. 403, 113821 (2022). https://doi.org/10.1016/j.cam.2021.113821
  * (43) Rauhut, H., Ward, R.: Sparse Legendre expansions via $\ell_{1}$-minimization. J. Approx. Theory 164(5), 517–533 (2012). https://doi.org/10.1016/j.jat.2012.01.008
  * (44) Potts, D., Volkmer, T.: Fast and exact reconstruction of arbitrary multivariate algebraic polynomials in Chebyshev form. In: 2015 International Conference on Sampling Theory and Applications (SampTA), pp. 392–396 (2015). https://doi.org/10.1109/SAMPTA.2015.7148919
  * (45) Nasdala, R., Potts, D.: A note on transformed Fourier systems for the approximation of non-periodic signals. arXiv:2102.00717 (2021). https://doi.org/10.48550/arXiv.2102.00717
  * (46) Shen, J., Tang, T., Wang, L.-L.: Spectral Methods. Springer Series in Computational Mathematics, vol. 41, p. 470. Springer, Heidelberg (2011). https://doi.org/10.1007/978-3-540-71041-7. Algorithms, analysis and applications

##  Appendix A Calculations for the $H^{2}(0,1)$ basis

Proof: [Proof of the first part of Theorem 4.2] Analogously to Theorem 4.1,
for $\sigma$ a singular value of $W$ with corresponding eigenfunction $\eta\in
H_{\mathrm{mix}}^{2}$, we obtain the following differential equation

$\displaystyle\frac{1-\sigma^{2}}{\sigma^{2}}\eta=\eta^{(4)}\quad\text{with}\quad\eta^{(2)}(0)=\eta^{(2)}(1)=\eta^{(3)}(0)=\eta^{(3)}(1)=0.$

Now we distinguish three cases for the value of $\sigma^{2}$:

First case. Let us assume $\sigma^{2}=1$. The ansatz function becomes

$\displaystyle\eta(x)=A+Bx+Cx^{2}+Dx^{3}.$

From the conditions $\eta^{(2)}(0)=\eta^{(3)}(0)=0$ we obtain $D=C=0$. The two
remaining degrees of freedom are restricted by demanding
$L_{2}(0,1)$-orthonormality. By simple calculus we obtain the proposed
eigenfunctions $\eta_{0}$ and $\eta_{1}$.

Second case. Lets assume
$(1-\sigma^{2})/\sigma^{2}>0\Leftrightarrow\sigma^{2}<1$. Introducing
$t\coloneqq\sqrt[4]{(1-\sigma^{2})/\sigma^{2}}$, we use the ansatz

$\displaystyle\eta(x)=A\cos(tx)+B\sin(tx)+C\cosh(tx)+D\sinh(tx).$

The conditions $\eta^{(2)}(0)=\eta^{(3)}(0)=0$ transform to $A=C$ and $B=D$,
respectively. The conditions $\eta^{(2)}(1)=\eta^{(3)}(1)=0$ can be put into a
system of equations:

$\displaystyle\begin{bmatrix}\cosh(t)-\cos(t)&\sinh(t)-\sin(t)\\\
\sinh(t)+\sin(t)&\cosh(t)-\cos(t)\end{bmatrix}\begin{bmatrix}A\\\
B\end{bmatrix}=\bm{0}$

or, by using $\cosh^{2}(t)-\sinh^{2}(t)=\cos^{2}(t)+\sin^{2}(t)=1$,
equivalently

$\displaystyle\begin{bmatrix}\cosh(t)-\cos(t)&\sinh(t)-\sin(t)\\\
0&1-\cosh(t)\cos(t)\end{bmatrix}\begin{bmatrix}A\\\ B\end{bmatrix}=\bm{0}.$

For non-trivial solutions we need non-regularity of that matrix which
transforms to the condition $\cosh(t)\cos(t)=1$. With the leftover degree of
freedom we choose

$\displaystyle A=C=1\quad\text{and}\quad
B=D=-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}$

and obtain $\eta_{k}$ for $k\geq 2$ as proposed in the theorem.

For the $L_{2}$-norm we obtain

$\displaystyle\int_{0}^{1}\lvert\eta_{n}\rvert^{2}\;dx\quad=\int_{0}^{1}(\cosh(tx)+\cos(tx))^{2}\;dx+B^{2}\int_{0}^{1}(\sinh(tx)+\sin(tx))^{2}\;dx$
$\displaystyle\qquad+2B\int_{0}^{1}(\cosh(tx)+\cos(tx))(\sinh(tx)+\sin(tx))\;dx$
$\displaystyle\quad=1+\frac{\sin(2t)+\sinh(2t)+4\cos(t)\sinh(t)+4\sin(t)\cosh(t)}{4t}$
$\displaystyle\qquad+B^{2}\frac{-\sin(2t)+\sinh(2t)-4\cos(t)\sinh(t)+4\sin(t)\cosh(t)}{4t}+2B\frac{(\sin(t)+\sinh(t))^{2}}{2t}$
$\displaystyle\quad=1+\frac{1+B^{2}}{4t}(\sinh(2t)+4\sin(t)\cosh(t))+\frac{1-B^{2}}{4t}(\sin(2t)+4\cos(t)\sinh(t))$
$\displaystyle\qquad+B\frac{(\sin(t)+\sinh(t))^{2}}{t}$
$\displaystyle\quad=1+\frac{1+B^{2}}{2t}\cosh(t)(\sinh(t)+2\sin(t))+\frac{1-B^{2}}{2t}\cos(t)(\sin(t)+2\sinh(t))$
$\displaystyle\qquad+B\frac{(\sin(t)+\sinh(t))^{2}}{t}\,.$

Using $\cos(t)\cosh(t)=1$, we obtain

$\displaystyle 1+B^{2}=2\frac{\sinh(t)}{\sinh(t)-\sin(t)}\quad\text{and}\quad
1-B^{2}=-2\frac{\sin(t)}{\sinh(t)-\sin(t)}\,.$ (A.1)

Thus,

$\displaystyle\int_{0}^{1}\lvert\eta_{n}\rvert^{2}\;dx$
$\displaystyle=1+\frac{1}{t(\sinh(t)-\sin(t))}\Bigg{[}\sinh(t)\cosh(t)(\sinh(t)+2\sin(t))$
$\displaystyle\quad-\sin(t)\cos(t)(\sin(t)+2\sinh(t))-(\cosh(t)-\cos(t))(\sin(t)+\sinh(t))^{2}\Bigg{]}$
$\displaystyle=1+\frac{\cos(t)\sinh^{2}(t)-\cosh(t)\sin^{2}(t)}{t(\sinh(t)-\sin(t))}$
$\displaystyle=1+\frac{\cos(t)\cosh^{2}(t)-\cos(t)-\cosh(t)+\cosh(t)\cos^{2}(t)}{t(\sinh(t)-\sin(t))}$

where $\cos^{2}(t)+\sin^{2}(t)=\cosh^{2}(t)-\sinh^{2}(t)=1$ was used in the
last equality. Using $\cosh(t)\cos(t)=1$, the latter summand evaluates to zero
and we have proven the $L_{2}$-normality.

Third case. Assume $\sigma^{2}>1$. Set
$t\coloneqq\sqrt[4]{(\sigma^{2}-1)/(\sigma^{2})}$. The ansatz becomes

$\displaystyle\eta(x)$ $\displaystyle=A\cosh(tx)\cos(tx)+B\cosh(tx)\sin(tx)$
$\displaystyle\phantom{=}+C\sinh(tx)\cos(tx)+D\sinh(tx)\sin(tx).$

The conditions $\eta^{(2)}(0)=\eta^{(3)}(0)=0$ transform to $D=0$ and $B=C$.
The two remaining degrees of freedom are fixed by the conditions
$\eta^{(2)}(1)=\eta^{(3)}(1)=0$ which, in matrix form, look as follows

$\displaystyle\begin{bmatrix}-\sinh(t)\sin(t)&\sinh(t)\cos(t)-\cosh(t)\sin(t)\\\
-\sinh(t)\cos(t)-\cosh(t)\sin(t)&-2\sinh(t)\sin(t)\end{bmatrix}\begin{bmatrix}A\\\
B\end{bmatrix}=\begin{bmatrix}0\\\ 0\end{bmatrix}.$

For a non-trivial solution we need that matrix to be non-regular. To achieve
that we have a look at the roots of its determinant:

$\displaystyle
2\sinh^{2}(t)\sin^{2}(t)+\sinh^{2}(t)\cos^{2}(t)-\cosh^{2}(t)\sin^{2}(t)\overset{!}{=}0.$

Using $\sin^{2}(t)+\cos^{2}(t)=\cosh^{2}(t)-\sinh^{2}(t)=1$ we have

$\displaystyle\sinh^{2}(t)-\sin^{2}(t)=\frac{1}{2}\cosh(2t)+\frac{1}{2}\cos(2t)-1\overset{!}{=}0$

which is only fulfilled for $t=0$, or equivalently, $\sigma^{2}=1$. Hence,
there are no eigenvalues bigger than 1.

######  Lemma A.1.

For $0<t_{2}<t_{3}<\dots$ fulfilling $\cosh(t_{k})\cos(t_{k})=1$ and
$\tilde{t}_{k}=\frac{2k-1}{2}\pi$, we have

$\displaystyle\frac{3}{2}\pi<t_{2}\quad\text{and}\quad\left\lvert\tilde{t}_{k}-t_{k}\right\rvert\leq\varepsilon$

for $k\geq\frac{1}{\pi}\log(\pi/\varepsilon)$. In particular
$\lvert\tilde{t}_{k}-t_{k}\rvert\leq\pi\exp(-2\pi)$ for all $k\geq 2$.

Proof:

$x$$y$$\pi/2$$3\pi/2$$5\pi/2$ Figure A.1: $\cos(t)$ and $1/\cosh(t)$

Since $0<1/\cosh(t)<1$ for $t>0$ and the oscillating behaviour of $\cos(t)$,
as depicted in Figure A.1, we obtain

$\displaystyle
t_{k}\in\begin{cases}\left(\frac{2k-1}{2}\pi,\frac{2k}{2}\pi\right)&\text{for
}k\text{ even}\\\ \left(\frac{2k-2}{2}\pi,\frac{2k-1}{2}\pi\right)&\text{for
}k\text{ odd}\,.\end{cases}$

In particular, $\frac{3}{2}\pi<t_{2}$. Furthermore, for even $k$ and
$t\in\left(\frac{2k-1}{2}\pi,\frac{2k}{2}\pi\right)$ we have

$\displaystyle\frac{1}{\cosh(t)}\leq 2\exp(-t)\leq
2\exp\left(-\frac{2k-1}{2}\pi\right)\quad\text{and}\quad\cos(t)\geq\frac{t-\frac{2k-1}{2}\pi}{\pi/2}.$

The function bounds intersect for a value larger than $t_{k}$, which we use to
refine the interval:

$\displaystyle
t_{k}\in\left(\tilde{t}_{k},\tilde{t}_{k}+\pi\exp\left(-\frac{2k-1}{2}\pi\right)\right).$

Similarly, for odd $k$ and
$t\in\left(\frac{2k-2}{2}\pi,\frac{2k-1}{2}\pi\right)$ we obtain

$\displaystyle
t_{k}\in\left(\tilde{t}_{k}-\pi\exp\left(-\frac{2k-2}{2}\pi\right),\tilde{t}_{k}\right).$

Thus, for $k\geq 2$ we have
$\left\lvert\tilde{t}_{k}-t_{k}\right\rvert\leq\pi\exp(-(k-1)\pi)$, which is
smaller than $\varepsilon$ for

$\displaystyle
k\geq\frac{\log(\pi/\varepsilon)}{\pi}+1\geq\frac{\log(\pi/\varepsilon)}{\pi}\,.$

######  Lemma A.2.

For $0<t_{2}<t_{3}<\dots$ fulfilling $\cosh(t_{k})\cos(t_{k})=1$, we have that
$\eta_{t_{k}}$ defined by

$\displaystyle\eta_{t_{k}}^{I}(x)=\cosh(t_{k}x)-\frac{\cosh(t_{k})-\cos(t_{k})}{\sinh(t_{k})-\sin(t_{k})}\sinh(t_{k}x)$
(A.2)

is convex and non-negative for all even $k$ and monotone for all odd $k$.

Proof: Step 1. We distinguish for different values of
$B=B(t)\coloneqq(\cosh(t)-\cos(t))/(\sinh(t)-\sin(t))$. For $B<1$ we have

$\displaystyle\eta_{t}^{I}(x)=\cosh(tx)-B\sinh(tx)\geq\cosh(tx)-\sinh(tx)\geq
0$

and by the same argument $(\eta_{t}^{I}(x))^{(2)}=t^{2}\eta_{t}^{I}(x)\geq 0$
for all $x\geq 0$. Thus, $\eta_{t}^{I}(x)$ is convex and non-negative.

For $B>1$ we obtain

$\displaystyle(\eta_{t}^{{I}})^{\prime}=t(\sinh(tx)-B\cosh(tx))\leq
t(\sinh(tx)-\cosh(tx))\leq 0$

for all $x\geq 0$. Thus, $\eta_{t}^{{I}}$ is monotone.

Step 2. It is left to show for which $k$’s $B(t_{k})$ attains a value smaller
or bigger than one:

$\displaystyle B(t_{k})\lessgtr 1$
$\displaystyle\quad\Leftrightarrow\quad\cosh(t_{k})-\cos(t_{k})\lessgtr\sinh(t_{k})-\sin(t_{k})$
$\displaystyle\quad\Leftrightarrow\quad\exp(-t_{k})-\sqrt{2}\cos(t_{k}+\pi/4)\lessgtr
0.$

We will show that $\exp(-t_{k})-\sqrt{2}\cos(t_{k}+\pi/4)$ has the same sign
as $(-1)^{k+1}$ and, thus, are finished. We do this by estimating their
difference by a quantity smaller than one. With
$\tilde{t}_{k}=\frac{2k-1}{2}\pi$ we obtain

$\displaystyle\lvert\exp(-t_{k})-\sqrt{2}\cos(t_{k}+\pi/4)-(-1)^{k+1}\rvert$
$\displaystyle\quad=\lvert\exp(-t_{k})-\sqrt{2}\cos(t_{k}+\pi/4)+\sqrt{2}\cos(\tilde{t}_{k}+\pi/4)\rvert.$

Using that $\cos$ is Lipschitz-continuous with constant $1$ and Lemma A.1 we
estimate the above by

$\displaystyle\lvert\exp(-t_{k})-\sqrt{2}\cos(t_{k}+\pi/4)-(-1)^{k+1}\rvert$
$\displaystyle\leq\lvert\exp(-t_{k})\rvert+\sqrt{2}\lvert
t_{k}-\tilde{t}_{k}\rvert$
$\displaystyle\leq\exp(-3/2\pi)+\sqrt{2}\pi\exp\left(-2\pi\right),$

which is certainly smaller than one.

######  Lemma A.3.

For $0<t_{2}<t_{3}<\dots$ fulfilling $\cosh(t_{k})\cos(t_{k})=1$, we have that
$\eta_{t_{k}}^{I}$ defined in (A.2) is even with respect to the axis $x=1/2$
for all even $k$ and vice versa.

Proof: Step 1. We will show that $\eta_{t_{k}}^{I}$ has any symmetry around
$x=1/2$. We shift the function and split it into an odd and an even part. For
$B=B(t)=(\cosh(t)-\cos(t))/(\sinh(t)-\sin(h))$, we obtain

$\displaystyle\eta_{t}^{I}(x+1/2)=\cosh(tx+t/2)-B\sinh(tx+t/2)$
$\displaystyle\quad=\underbrace{(\cosh(t/2)-B\sinh(t/2))}_{\eqqcolon\alpha}\cosh(tx)+\underbrace{(\sinh(t/2)-B\cosh(t/2))}_{\eqqcolon\beta}\sinh(tx).$

Multiplying the two factors $\alpha$ and $\beta$ in front of $\cosh(tx)$ and
$\sinh(tx)$, we obtain

$\displaystyle\alpha\cdot\beta$
$\displaystyle=-B\cosh^{2}(t/2)-B\sinh^{2}(t/2)+(1+B^{2})\cosh(t/2)\sinh(t/2)$
$\displaystyle=-B\frac{\cosh(t)-1}{2}-B\frac{\cosh(t)+1}{2}+(1+B^{2})\frac{\sinh(t)}{2}$
$\displaystyle=-B\cosh(t)+(1+B^{2})\frac{\sinh(t)}{2}\,.$

Using (A.1), $\cosh(t)\cos(t)=1$, and $1=\cosh^{2}(t)-\sinh^{2}(t)$ this
evaluates to

$\displaystyle\alpha\cdot\beta=-\frac{\cosh^{2}(t)-1}{\sinh(t)-\sin(t)}+\frac{\sinh^{2}(t)}{\sinh(t)-\sin(t)}=0\,.$

And since we are not dealing with the zero function either $\alpha$ or $\beta$
is zero. Thus, $x\mapsto\eta_{t}^{I}(x+1/2)$ obeys a symmetry.

Step 2. It remains to specify the kind of symmetry. By Lemma A.2 we have that
$\eta_{t}^{I}$ is convex for even $k$. Since a convex non-constant function
cannot be odd it has to be even. Also by Lemma A.2 we have that $\eta_{t}^{I}$
is monotone for odd $k$. Since a monotone non-zero function cannot be even it
has to be odd.

Proof: [Proof of the second part of Theorem 4.2] The cases $k\in\\{0,1\\}$ are
clear. For $k\geq 2$ we split the function into $\eta_{t}^{I}$ defined in
(A.2) and

$\displaystyle\eta_{t}^{{II}}(x)=\cos(tx)-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\sin(tx).$

We will show that each of these is bounded by $1.01\sqrt{2}$ and, thus, obtain
the assertion.

Step 1. In order to bound $\eta_{t_{k}}^{{I}}$ we firstly have a look at the
boundary points $x\in\\{0,1\\}$. With Lemma A.3 we obtain

$\eta_{t_{k}}^{{I}}(0)=\left\lvert\eta_{t_{k}}^{I}(1)\right\rvert=1.$ (A.3)

By Lemma A.2 $\eta_{t_{k}}^{I}$ is either non-negative and convex or monotone
and, thus, cannot exceed its values on the boundary.

Step 2. In order to bound $\eta_{t}^{{II}}$ we define

$\displaystyle
B\coloneqq\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\quad\text{and}\quad\vartheta=\arg(1+Bi).$

Next, we use the exponential definition of sine and cosine and the polar
representation of complex numbers to obtain

$\displaystyle\eta_{t}^{{II}}(x)$ $\displaystyle=\cos(tx)-B\sin(tx)$
$\displaystyle=\frac{\exp(itx)+\exp(-itx)}{2}+Bi\frac{\exp(itx)-\exp(-itx)}{2}$
$\displaystyle=\frac{(1+Bi)\exp(itx)+(1-Bi)\exp(-itx)}{2}$
$\displaystyle=\sqrt{1+B^{2}}\frac{\exp(i(tx+\vartheta))+\exp(-i(tx+\vartheta))}{2}$
$\displaystyle=\sqrt{1+B^{2}}\cos(tx+\vartheta)$

Thus, by (A.1)

$\displaystyle\lvert\eta_{t}^{{II}}(x)\rvert\leq\sqrt{1+B^{2}}=\sqrt{\frac{2}{1-\sin(t)/\sinh(t)}}\leq\sqrt{\frac{2}{1-1/\sinh(t)}}$
(A.4)

From Lemma A.1 we use $t\geq 3/2\pi$ in combination with the monotonicity in
(A.4) we have $\lvert\eta_{t}^{{II}}(x)\rvert\leq 1.01\sqrt{2}$.

######  Lemma A.4.

For $t\geq\max\\{2\log(4/\varepsilon),3/2\pi\\}$ we have for $x\in[0,1/2]$

$\displaystyle\left\lvert\left(1-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\right)\sinh(tx)\right\rvert\leq\varepsilon\,.$

Proof: We use $\cosh(t)-\sinh(t)=\exp(-t)$ and
$\cos(t)-\sin(t)=\sqrt{2}\cos(t+\pi/4)$ to estimate

$\displaystyle\left\lvert\left(1-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\right)\sinh(tx)\right\rvert$
$\displaystyle=\lvert\sqrt{2}\cos(t+\pi/4)-\exp(-t)\rvert\left\lvert\frac{\sinh(tx)}{\sinh(t)-\sin(t)}\right\rvert.$

Since $x\leq 1/2$, $\sinh$ strictly monotone growing, and $t\geq 3/2\pi$ by
Lemma A.1, we further estimate

$\displaystyle\left\lvert\left(1-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\right)\sinh(tx)\right\rvert$
$\displaystyle\leq
2\left\lvert\frac{\sinh(t/2)}{\sinh(t)-\sin(t)}\right\rvert$
$\displaystyle=2\left\lvert\frac{1}{2\cosh(t/2)}\frac{1}{1-\sin(t)/\sinh(t)}\right\rvert$

Using $1-\sin(t)/\sinh(t)>1/2$ for $t>3/2\pi$, we obtain

$\displaystyle\left\lvert\left(1-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\right)\sinh(tx)\right\rvert\leq\frac{2}{\cosh(t/2)}\leq\frac{4}{\exp(t/2)}\,,$

which is smaller than $\varepsilon$ for $t\geq 2\log(4/\varepsilon)$.

######  Lemma A.5.

For $0<t_{2}<t_{3}<\dots$ fulfilling $\cosh(t_{k})\cos(t_{k})=1$, we have

$\displaystyle\left\lvert\eta_{t_{k}}^{{II}}(x)-\sqrt{2}\cos(t_{k}x+\pi/4)\right\rvert$
$\displaystyle\leq\varepsilon\quad\text{for}\quad x\in[0,1]$ and
$\displaystyle\left\lvert\eta_{t_{k}}^{{I}}(x)-\exp(-tx)\right\rvert$
$\displaystyle\leq\varepsilon\quad\text{for}\quad x\in[0,1/2]$

for $k\geq\frac{2}{\pi}\log(4/\varepsilon)+1$.

Proof: Step 1. For the first inequality we use

$\displaystyle\sqrt{2}\cos(tx+\pi/4)=\cos(tx)-\sin(tx)$

to obtain

$\displaystyle\left\lvert\eta_{t}^{{II}}-\sqrt{2}\cos(tx+\pi/4)\right\rvert=\left\lvert\left(1-\frac{\cosh(t)-\cos(t)}{\sinh(t)-\sin(t)}\right)\sin(tx)\right\rvert$

which is smaller than $\varepsilon$ for
$t>\max\\{2\log(4/\varepsilon),3/2\pi\\}$ by Lemma A.4.

The second inequality follows analogously from $\exp(-tx)=\cosh(tx)-\sinh(tx)$
and Lemma A.4.

Step 2. It is left to show the condition
$t\geq\max\\{2\log(4/\varepsilon),3/2\pi\\}$ from Step 1. By Lemma A.1 we have
$t_{k}\geq 3/2\pi$. Further, by assumption, we have

$\displaystyle
k\geq\frac{2}{\pi}\log\left(\frac{4}{\varepsilon}\right)+1\geq\frac{2}{\pi}\log\left(\frac{4}{\varepsilon}\right)+\exp(-2\pi)+\frac{1}{2}$

Thus,

$\displaystyle
2\log\left(\frac{4}{\varepsilon}\right)\leq\frac{2k-1}{2}\pi-\pi\exp(-2\pi)\leq
t_{k}$

where the last inequality follows from Lemma A.1.

Proof: [Proof of Theorem 4.3] Because of the symmetry shown in Lemma A.3 we
assume without loss of generality $x\in[0,1/2]$. Then

$\displaystyle\lvert\eta_{n}(x)-\tilde{\eta}_{n}(x)\rvert\leq\left\lvert\eta_{n}^{I}(x)-\exp(-t_{k}x)\right\rvert+\left\lvert\eta_{n}^{{II}}(x)-\sqrt{2}\cos(t_{k}x+\pi/4)\right\rvert$
$\displaystyle\quad+\left\lvert\exp(-t_{k}x)-\exp(-\tilde{t}_{k}x)\right\rvert+\left\lvert\sqrt{2}\cos(t_{k}x+\pi/4)-\sqrt{2}\cos(\tilde{t}_{k}x+\pi/4)\right\rvert.$

By Lemma A.5, the first two summands are each smaller than $\varepsilon/4$
each for $n>\frac{2}{\pi}\log(16/\varepsilon)+1$. We estimate the two latter
summands as follows.

Since $\cos$ is Lipschitz continuous with constant one we have

$\displaystyle\left\lvert\sqrt{2}\cos(t_{k}x+\pi/4)-\sqrt{2}\cos\left(\tilde{t}_{k}x+\frac{\pi}{4}\right)\right\rvert\leq\left\lvert\sqrt{2}\left(t_{k}-\tilde{t}_{k}\right)\right\rvert$

which, by Lemma A.1 is smaller than $\varepsilon/4$ for
$k>\frac{1}{\pi}\log(4\pi/(\sqrt{2}\varepsilon))$.

Since $\exp$ is Lipschitz continuous with constant $1$ on $(-\infty,0)$, we
have

$\displaystyle\left\lvert\exp(-t_{k}x)-\exp(-\tilde{t}_{k}x)\right\rvert\leq\left\lvert
t_{k}-\tilde{t}_{k}\right\rvert$

which, by Lemma A.1 is smaller than $\varepsilon/4$ for
$k>\frac{1}{\pi}\log(16/\varepsilon)$.

Overall, we obtain
$\lvert\eta_{n}(x)-\tilde{\eta}_{n}(x)\rvert<4\frac{\varepsilon}{4}=\varepsilon$
for

$\displaystyle
k>\max\left\\{\frac{2}{\pi}\log(16/\varepsilon)+1,\,\frac{1}{\pi}\log(4\pi/(\sqrt{2}\varepsilon)),\,\frac{1}{\pi}\log(16/\varepsilon)\right\\}=\frac{2}{\pi}\log(16/\varepsilon)+1\,.$