Nonparametric regression using deep neural networks with ReLU activation function

Introduction

Multilayer (or deep) neural networks have been successfully trained recently to achieve impressive results for complicated tasks such as object detection on images and speech recognition. Deep learning is now considered to be the state-of-the art for these tasks. But there is a lack of mathematical understanding. One problem is that fitting a neural network to data is highly non-linear in the parameters. Moreover, the function class is non-convex and different regularization methods are combined in practice.

This article is inspired by the idea to build a statistical theory that provides some understanding of these procedures. As the full method is too complex to be theoretically tractable, we need to make some selection of important characteristics that we believe are crucial for the success of the procedure.

The statistical analysis for the ReLU activation function is quite different from earlier approaches and we discuss this in more detail in the overview on related literature in Section 6. Viewed as a nonparametric method, ReLU networks have some surprising properties. To explain this, notice that deep networks with ReLU activation produce functions that are piecewise linear in the input. Nonparametric methods which are based on piecewise linear approximations are typically not able to capture higher-order smoothness in the signal and are rate-optimal only up to smoothness index two. Interestingly, ReLU activation combined with a deep network architecture achieves near minimax rates for arbitrary smoothness of the regression function.

The number of hidden layers of state-of-the-art network architectures has been growing over the past years, cf. . There are versions of the recently developed deep network ResNet which are based on $152$ layers, cf. . Our analysis indicates that for the ReLU activation function the network depth should be scaled with the sample size. This suggests, that for larger samples, additional hidden layers should be added.

Recent deep architectures include more network parameters than training samples. The well-known AlexNet for instance is based on $60$ million network parameters using only 1.2 million samples. We account for high-dimensional parameter spaces in our analysis by assuming that the number of potential network parameters is much larger than the sample size. For noisy data generated from the nonparametric regression model, overfitting does not lead to good generalization errors and incorporating regularization or sparsity in the estimator becomes essential. In the deep networks literature, one option is to make the network thinner assuming that only few parameters are non-zero (or active), cf. , Section 7.10. Our analysis shows that the number of non-zero parameters plays the role of the effective model dimension and - as is common in non-parametric regression - needs to be chosen carefully.

Existing statistical theory often requires that the size of the network parameters tends to infinity as the sample size increases. In practice, estimated network weights are, however, rather small. We can incorporate small parameters in our theory, proving that it is sufficient to consider neural networks with all network parameters bounded in absolute value by one.

Multilayer neural networks are typically applied to high-dimensional input. Without additional structure in the signal besides smoothness, nonparametric estimation rates are then slow because of the well-known curse of dimensionality. This means that no statistical procedure can do well regarding pointwise reconstruction of the signal. Multilayer neural networks are believed to be able to adapt to many different structures in the signal, therefore avoiding the curse of dimensionality and achieving faster rates in many situations. In this work, we stick to the regression setup and show that deep ReLU networks can indeed attain faster rates under a hierarchical composition assumption on the regression function, which includes (generalized) additive models and the composition models considered in .

Parts of the success of multilayer neural networks can be explained by the fast algorithms that are available to estimate the network weights from data. These iterative algorithms are based on minimization of some empirical loss function using stochastic gradient descent. Because of the non-convex function space, gradient descent methods might get stuck in a saddle point or converge to one of the potentially many local minima. derives a heuristic argument showing that the risk of most of the local minima is not much larger than the risk of the global minimum. Despite the huge number of variations of the stochastic gradient descent, the common objective of all approaches is to reduce the empirical loss. In our framework we associate to any network reconstruction method a parameter quantifying the expected discrepancy between the achieved empirical risk and the global minimum of the energy landscape. The main theorem then states that a network estimator is minimax rate optimal (up to log factors) if and only if the method almost minimizes the empirical risk.

We also show that wavelet series estimators are unable to adapt to the underlying structure under the composition assumption on the regression function. By deriving lower bounds, it is shown that the rates are suboptimal by a polynomial factor in the sample size $n.$ This provides an example of a function class for which fitting a neural network outperforms wavelet series estimators.

Our setting deviates in two aspects from the computer science literature on deep learning. Firstly, we consider regression and not classification. Secondly, we restrict ourselves in this article to multilayer feedforward artificial neural networks, while most of the many recent deep learning applications have been obtained using specific types of networks such as convolutional or recurrent neural networks.

The article is structured as follows. Section 2 introduces multilayer feedforward artificial neural networks and discusses mathematical modeling. This section also contains the definition of the network classes. The considered function classes for the regression function and the main result can be found in Section 3. Specific structural constraints such as additive models are discussed in Section 4. In Section 5 it is shown that wavelet estimators can only achieve suboptimal rates under the composition assumption. We give an overview of relevant related literature in Section 6. The proof of the main result together with additional discussion can be found in Section 7.

Notation: Vectors are denoted by bold letters, e.g. $\mathbf{x}:=(x_{1},\ldots,x_{d})^{\top}.$ As usual, we define $|\mathbf{x}|_{p}:=(\sum_{i=1}^{d}|x_{i}|^{p})^{1/p},$ $|\mathbf{x}|_{\infty}:=\max_{i}|x_{i}|,$ $|\mathbf{x}|_{0}:=\sum_{i}\mathbf{1}(x_{i}\neq 0),$ and write $\|f\|_{p}:=\|f\|_{L^{p}(D)}$ for the $L^{p}$-norm on $D,$ whenever there is no ambiguity of the domain $D.$ For two sequences $(a_{n})_{n}$ and $(b_{n})_{n},$ we write $a_{n}\lesssim b_{n}$ if there exists a constant $C$ such that $a_{n}\leq Cb_{n}$ for all $n.$ Moreover, $a_{n}\asymp b_{n}$ means that $a_{n}\lesssim b_{n}$ and $b_{n}\lesssim a_{n}.$ We denote by $\log_{2}$ the logarithm with respect to the basis two and write $\lceil x\rceil$ for the smallest integer $\geq x.$

Mathematical definition of multilayer neural networks

In computer science, neural networks are more commonly introduced via their representation as directed acyclic graphs, cf. Figure 1. Using this equivalent definition, the nodes in the graph (also called units) are arranged in layers. The input layer is the first layer and the output layer the last layer. The layers that lie in between are called hidden layers. The number of hidden layers corresponds to $L$ and the number of units in each layer generates the width vector $\mathbf{p}.$ Each node/unit in the graph representation stands for a scalar product of the incoming signal with a weight vector which is then shifted and applied to the activation function.

Mathematical modeling of deep network characteristics: Given a network function $f(\mathbf{x})=W_{L}\sigma_{\mathbf{v}_{L}}W_{L-1}\sigma_{\mathbf{v}_{L-1}}\cdots W_{1}\sigma_{\mathbf{v}_{1}}W_{0}\mathbf{x},$ the network parameters are the entries of the matrices $(W_{j})_{j=0,\ldots,L}$ and vectors $(\mathbf{v}_{j})_{j=1,\ldots,L}.$ These parameters need to be estimated/learned from the data.

The aim of this article is to consider a framework that incorporates essential features of modern deep network architectures. In particular, we allow for large depth $L$ and large number of potential network parameters. For the main result, no upper bound on the number of network parameters is needed. Thus we consider high-dimensional settings with more parameters than training data.

Another characteristic of trained networks is that the size of the learned network parameters is typically not very large. Common network initialization methods initialize the weight matrices $W_{j}$ by a (nearly) orthogonal random matrix if two successive layers have the same width, cf. , Section 8.4. In practice, the trained network weights are typically not far from the initialized weights. Since in an orthogonal matrix all entries are bounded in absolute value by one, also the trained network weights will not be large.

Existing theoretical results, however, often require that the size of the network parameters tends to infinity. If large parameters are allowed, one can for instance easily approximate step functions by ReLU networks. To be more in line with what is observed in practice, we consider networks with all parameters bounded by one. This constraint can be easily build into the deep learning algorithm by projecting the network parameters in each iteration onto the interval $.$

If $\|W_{j}\|_{\infty}$ denotes the maximum-entry norm of $W_{j},$ the space of network functions with given network architecture and network parameters bounded by one is

with the convention that $\mathbf{v}_{0}$ is a vector with all components equal to zero.

In deep learning, sparsity of the neural network is enforced through regularization or specific forms of networks. Dropout for instance sets randomly units to zero and has the effect that each unit will be active only for a small fraction of the data, cf. , Section 7.2. In our notation this means that each entry of the vectors $\sigma_{\mathbf{v}_{k}}W_{k-1}\cdots W_{1}\sigma_{\mathbf{v}_{1}}W_{0}\mathbf{x},$ $k=1,\ldots,L$ is zero over a large range of the input space $\mathbf{x}\in^{d}.$ Convolutional neural networks filter the input over local neighborhoods. Rewritten in the form (2) this essentially means that the $W_{i}$ are banded Toeplitz matrices. All network parameters corresponding to higher off-diagonal entries are thus set to zero.

In this work we model the network sparsity assuming that there are only few non-zero/active network parameters. If $\|W_{j}\|_{0}$ denotes the number of non-zero entries of $W_{j}$ and $\||f|_{\infty}\|_{\infty}$ stands for the sup-norm of the function $\mathbf{x}\mapsto|f(\mathbf{x})|_{\infty},$ then the $s$-sparse networks are given by

The upper bound on the uniform norm of $f$ is most of the time dispensable and therefore omitted in the notation. We consider cases where the number of network parameters $s$ is small compared to the total number of parameters in the network.

In deep learning, it is common to apply variations of stochastic gradient descent combined with other techniques such as dropout to the loss induced by the log-likelihood (see Section 6.2.1.1 in ). For nonparametric regression with normal errors, this coincides with the least-squares loss (in machine learning terms this is the cross-entropy for this model, cf. , p.129). The common objective of all reconstruction methods is to find networks $f$ with small empirical risk $\tfrac{1}{n}\sum_{i=1}^{n}(Y_{i}-f(\mathbf{X}_{i}))^{2}.$ For any estimator $\widehat{f}_{n}$ that returns a network in the class $\mathcal{F}(L,\mathbf{p},s,F)$ we define the corresponding quantity

To evaluate the statistical performance of an estimator $\widehat{f}_{n},$ we derive bounds for the prediction error

with $\mathbf{X}\stackrel{{\scriptstyle\mathcal{D}}}{{=}}\mathbf{X}_{1}$ being independent of the sample $(\mathbf{X}_{i},Y_{i})_{i}.$

The term $\Delta_{n}(\widehat{f}_{n},f_{0})$ can be related via empirical process theory to constant $\times(R(\widehat{f}_{n},f_{0})-R(\widehat{f}_{n}^{\operatorname{ERM}},f_{0}))+$ remainder, with $\widehat{f}_{n}^{\operatorname{ERM}}$ an empirical risk minimizer. Therefore, $\Delta_{n}(\widehat{f}_{n},f_{0})$ is the key quantity that together with the minimax estimation rate sharply determines the convergence rate of $\widehat{f}_{n}$ (up to $\log n$-factors). Determining the decay of $\Delta_{n}(\widehat{f}_{n},f_{0})$ in $n$ for commonly employed methods such as stochastic gradient descent is an interesting problem in its own. We only sketch a possible proof strategy here. Because of the potentially many local minima and saddle points of the loss surface or energy landscape, gradient descent based methods have only a small chance to reach the global minimum without getting stuck in a local minimum first. By making a link to spherical spin glasses, provide a heuristic suggesting that the loss of any local minima lies in a band that is lower bounded by the loss of the global minimum. The width of the band depends on the width of the network. If the heuristic argument can be made rigorous, then the width of the band provides an upper bound for $\Delta_{n}(\widehat{f}_{n},f_{0})$ for all methods that converge to a local minimum. This would allow us then to study deep learning without an explicit analysis of the algorithm. For more on the energy landscape, see .

Main results

The theoretical performance of neural networks depends on the underlying function class. The classical approach in nonparametric statistics is to assume that the regression function is $\beta$-smooth. The minimax estimation rate for the prediction error is then $n^{-2\beta/(2\beta+d)}.$ Since the input dimension $d$ in neural network applications is very large, these rates are extremely slow. The huge sample sizes often encountered in deep learning applications are by far not sufficient to compensate the slow rates.

We therefore consider a function class that is natural for neural networks and exhibits some low-dimensional structure that leads to input dimension free exponents in the estimation rates. We assume that the regression function $f_{0}$ is a composition of several functions, that is,

with $g_{i}:[a_{i},b_{i}]^{d_{i}}\rightarrow[a_{i+1},b_{i+1}]^{d_{i+1}}.$ Denote by $g_{i}=(g_{ij})_{j=1,\ldots,d_{i+1}}^{\top}$ the components of $g_{i}$ and let $t_{i}$ be the maximal number of variables on which each of the $g_{ij}$ depends on. Thus, each $g_{ij}$ is a $t_{i}$-variate function. As an example consider the function $f_{0}(x_{1},x_{2},x_{3})=g_{11}(g_{01}(x_{1},x_{3}),g_{02}(x_{1},x_{2}))$ for which $d_{0}=3,$ $t_{0}=2,$ $d_{1}=t_{1}=2,$ $d_{2}=1.$ We always must have $t_{i}\leq d_{i}$ and for specific constraints such as additive models, $t_{i}$ might be much smaller than $d_{i}.$ The single components $g_{0},\ldots,g_{q}$ and the pairs $(\beta_{i},t_{i})$ are obviously not identifiable. As we are only interested in estimation of $f_{0}$ this causes no problems. Among all possible representations, one should always pick one that leads to the fastest estimation rate in Theorem 1 below.

It is conceivable that for many of the problems for which neural networks perform well a hidden hierarchical input-output relationship of the form (6) is present with small values $t_{i},$ cf. . Slightly more specific function spaces, which alternate between summations and composition of functions, have been considered in . We provide below an example of a function class that can be decomposed in the form (6) but is not contained in these spaces.

Recall that a function has Hölder smoothness index $\beta$ if all partial derivatives up to order $\lfloor\beta\rfloor$ exist and are bounded and the partial derivatives of order $\lfloor\beta\rfloor$ are $\beta-\lfloor\beta\rfloor$ Hölder, where $\lfloor\beta\rfloor$ denotes the largest integer strictly smaller than $\beta.$ The ball of $\beta$-Hölder functions with radius $K$ is then defined as

We assume that each of the functions $g_{ij}$ has Hölder smoothness $\beta_{i}.$ Since $g_{ij}$ is also $t_{i}$-variate, $g_{ij}\in\mathcal{C}_{t_{i}}^{\beta_{i}}([a_{i},b_{i}]^{t_{i}},K_{i})$ and the underlying function space becomes

with $\mathbf{d}:=(d_{0},\ldots,d_{q+1}),$ $\mathbf{t}:=(t_{0},\ldots,t_{q}),$ $\bm{\beta}:=(\beta_{0},\ldots,\beta_{q}).$

For estimation rates in the nonparametric regression model, the crucial quantity is the smoothness of $f.$ Imposing smoothness on the functions $g_{i},$ we must then find the induced smoothness on $f.$ If, for instance, $q=1,$ $\beta_{0},\beta_{1}\leq 1,$ $d_{0}=d_{1}=t_{0}=t_{1}=1,$ then $f=g_{1}\circ g_{0}$ and $f$ has smoothness $\beta_{0}\beta_{1},$ cf. . We should then be able to achieve at least the convergence rate $n^{-2\beta_{0}\beta_{1}/(2\beta_{0}\beta_{1}+1)}.$ For $\beta_{1}>1,$ the rate changes. Below we see that the convergence of the network estimator is described by the effective smoothness indices

Recall the definition of $\Delta_{n}(\widehat{f}_{n},f_{0})$ in (5). We can now state the main result.

Consider the $d$-variate nonparametric regression model (1) for composite regression function (6) in the class $\mathcal{G}(q,\mathbf{d},\mathbf{t},\bm{\beta},K).$ Let $\widehat{f}_{n}$ be an estimator taking values in the network class $\mathcal{F}(L,(p_{i})_{i=0,\ldots,L+1},s,F)$ satisfying

$\sum_{i=0}^{q}\log_{2}(4t_{i}\vee 4\beta_{i})\log_{2}n\leq L\lesssim n\phi_{n},$

$n\phi_{n}\lesssim\min_{i=1,\ldots,L}p_{i},$

There exist constants $C,C^{\prime}$ only depending on $q,\mathbf{d},\mathbf{t},\bm{\beta},F,$ such that if $\Delta_{n}(\widehat{f}_{n},f_{0})\leq C\phi_{n}L\,\log^{2}n,$ then

and if $\Delta_{n}(\widehat{f}_{n},f_{0})\geq C\phi_{n}L\,\log^{2}n,$ then

In order to minimize the rate $\phi_{n}L\,\log^{2}n,$ the best choice is to choose $L$ of the order of $\log_{2}n.$ The rate in the regime $\Delta_{n}(\widehat{f}_{n},f_{0})\leq C\phi_{n}\,\log^{3}n$ becomes then

The convergence rate in Theorem 1 depends on $\phi_{n}$ and $\Delta_{n}(\widehat{f}_{n},f_{0}).$ Below we show that $\phi_{n}$ is a lower bound for the minimax estimation risk over this class. Recall that the term $\Delta_{n}(\widehat{f}_{n},f_{0})$ is large if $\widehat{f}_{n}$ has a large empirical risk compared to an empirical risk minimizer. Having this term in the convergence rate is unavoidable as it also appears in the lower bound in (9). Since for any empirical risk minimizer the $\Delta_{n}$-term is zero by definition, we have the following direct consequence of the main theorem.

Let $\widetilde{f}_{n}\in\mathop{\rm arg\min}_{f\in\mathcal{F}(L,\mathbf{p},s,F)}\sum_{i=1}^{n}(Y_{i}-f(\mathbf{X}_{i}))^{2}$ be an empirical risk minimizer. Under the same conditions as for Theorem 1, there exists a constant $C^{\prime},$ only depending on $q,\mathbf{d},\mathbf{t},\bm{\beta},F,$ such that

The number of network parameters in a fully connected network is of the order $\sum_{i=0}^{L}p_{i}p_{i+1}.$ This shows that Theorem 1 requires sparse networks. More specifically, the network has at least $\sum_{i=1}^{L}p_{i}-s$ completely inactive nodes, meaning that all incoming signal is zero. The choice $s\asymp n\phi_{n}\log n$ in condition $(iv)$ balances the squared bias and the variance. From the proof of the theorem convergence rates can also be derived if $s$ is chosen of a different order.

For convenience, Theorem 1 is stated without explicit constants. The proofs, however, are non-asymptotic although we did not make an attempt to minimize the constants. It is well-known that deep learning outperforms other methods only for large sample sizes. This indicates that the method might be able to adapt to underlying structure in the signal and therefore achieving fast convergence rates but with large constants or remainder terms which spoil the results for small samples.

The proof of the risk bounds in Theorem 1 is based on the following oracle-type inequality.

Consider the $d$-variate nonparametric regression model (1) with unknown regression function $f_{0},$ satisfying $\|f_{0}\|_{\infty}\leq F$ for some $F\geq 1.$ Let $\widehat{f}_{n}$ be any estimator taking values in the class $\mathcal{F}(L,\mathbf{p},s,F)$ and let $\Delta_{n}(\widehat{f}_{n},f_{0})$ be the quantity defined in (5). For any $\varepsilon\in(0,1],$ there exists a constant $C_{\varepsilon},$ only depending on $\varepsilon,$ such that with

One consequence of the oracle inequality is that the upper bounds on the risk become worse if the number of layers increases. In practice it also has been observed that too many layers lead to a degradation of the performance, cf. , , Section 4.4 and , Section 4. Residual networks can overcome this problem. But they are not of the form (2) and will need to be analyzed separately.

One may wonder whether there is anything special about ReLU networks compared to other activation functions. A close inspection of the proof shows that two specific properties of the ReLU function are used.

One of the advantages of deep ReLU networks is the projection property

that we can use to pass a signal without change through several layers in the network. This is important since the approximation theory is based on the construction of smaller networks for simpler tasks that might not all have the same network depth. To combine these subnetworks we need to synchronize the network depths by adding hidden layers that do not change the output. This can be easily realized by choosing the weight matrices in the network to be the identity (assuming equal network width in successive layers) and using (11), see also (18). This property is not only a theoretical tool. To pass an outcome without change to a deeper layer is also often helpful in practice and realized by so called skip connections in which case they do not need to be learned from the data. A specific instance are residual networks with ReLU activation function that are successfully applied in practice. The difference to standard feedforward networks is that if all networks parameters are set to zero in a residual network, the network becomes essentially the identity map. For other activation functions it is much harder to approximate the identity.

Another advantage of the ReLU activation is that all network parameters can be taken to be bounded in absolute value by one. If all network parameters are initialized by a value in $,$ this means that each network parameter only need to be varied by at most two during training. It is unclear whether other results in the literature for non-ReLU activation functions hold for bounded network parameters. An important step is the approximation of the square function $x\mapsto x^{2}.$ For any twice differentiable and non-linear activation function, the classical approach to approximate the square function by a network is to use rescaled second order differences $(\sigma(t+2xh)-2\sigma(t+xh)+\sigma(xh))/(h^{2}\sigma^{\prime\prime}(t))\rightarrow x^{2}$ for $h\rightarrow 0$ and a $t$ with $\sigma^{\prime\prime}(t)\neq 0$. To achieve a sufficiently good approximation, we need to let $h$ tend to zero with the sample size, making some of the network parameters necessarily very large.

The $\log^{2}n$-factor in the convergence rate $\phi_{n}L\log^{2}n$ is likely an artifact of the proof. Next we show that $\phi_{n}$ is a lower bound for the minimax estimation risk over the class $\mathcal{G}(q,\mathbf{d},\mathbf{t},\bm{\beta},K)$ in the interesting regime $t_{i}\leq\min(d_{0},\ldots,d_{i-1})$ for all $i.$ This means that no dimensions are added on deeper abstraction levels in the composition of functions. In particular, it avoids that $t_{i}$ is larger than the input dimension $d_{0}.$ Outside of this regime, it is hard to determine the minimax rate and in some cases it is even possible to find another representation of $f$ as a composition of functions which yields a faster convergence rate.

Consider the nonparametric regression model (1) with $\mathbf{X}_{i}$ drawn from a distribution with Lebesgue density on $^{d}$ which is lower and upper bounded by positive constants. For any non-negative integer $q,$ any dimension vectors $\mathbf{d}$ and $\mathbf{t}$ satisfying $t_{i}\leq\min(d_{0},\ldots,d_{i-1})$ for all $i,$ any smoothness vector $\bm{\beta}$ and all sufficiently large constants $K>0,$ there exists a positive constant $c,$ such that

where the $\inf$ is taken over all estimators $\widehat{f}_{n}.$

The $L^{2}$-minimax rate coincides in most regimes with the sup-norm rate obtained in Section 4.1 of for composition of two functions. But unlike the classical nonparametric regression model, the minimax estimation rates for $L^{2}$-loss and sup-norm loss differ for some setups by a polynomial power in the sample size $n.$

There are several recent results in approximation theory that provide lower bounds on the number of required network weights $s$ such that all functions in a function class can be approximated by a $s$-sparse network up to some prescribed error, cf. for instance . Results of this flavor can also be quite easily derived by combining the minimax lower bound with the oracle inequality. The argument is that if the same approximation rates would hold for networks with less parameters then we would obtain rates that are faster than the minimax rates, which clearly is a contradiction. This provides a new statistical route to establish approximation theoretic properties.

for some $\varepsilon\leq c_{2},$ then for any width vector $\mathbf{p}$ with $p_{0}=d$ and $p_{L+1}=1,$

A more refined argument using Lemma 4 instead of Theorem 2 yields also lower bounds for $L^{2}.$

Examples of specific structural constraints

In this section we discuss several well-studied special cases of compositional constraints on the regression function.

Additive models: In an additive model the regression function has the form

This can be written as a composition of functions

Suppose that $f_{i}\in\mathcal{C}_{1}^{\beta}(,K)$ for $i=1,\ldots,d.$ Then, $f:\,^{d}\stackrel{{\scriptstyle g_{0}}}{{\longrightarrow}}[-K,K]^{d}\stackrel{{\scriptstyle g_{1}}}{{\longrightarrow}}[-Kd,Kd].$ Since for any $\gamma>1,$ $g_{1}\in\mathcal{C}_{d}^{\gamma}([-K,K]^{d},(K+1)d),$

For network architectures $\mathcal{F}(L,\mathbf{p},s,F)$ satisfying $F\geq(K+1)d,$ $2\log_{2}(4(\beta\vee 2)d)\log n\leq L\lesssim\log n,$ $n^{1/(2\beta+1)}\lesssim\min_{i}p_{i}$ and $s\asymp n^{1/(2\beta+1)}\log n,$ we thus obtain by Theorem 1,

This coincides up to the $\log^{3}n$-factor with the minimax estimation rate.

Generalized additive models: Suppose the regression function is of the form

For network architectures satisfying the assumptions of Theorem 1, the bound on the estimation rate becomes

Theorem 3 shows that $n^{-2\beta(\gamma\wedge 1)/(2\beta(\gamma\wedge 1)+1)}+n^{-2\gamma/(2\gamma+1)}$ is also a lower bound. Let us also remark that for the special case $\beta=\gamma\geq 2$ and $\beta,\gamma$ integers, Theorem 2.1 of establishes the estimation rate $n^{-2\beta/(2\beta+1)}.$

Sparse tensor decomposition: Assume that the regression function $f_{0}$ has the form

For networks with architectures satisfying $3\log_{2}(4(\beta+1)(d+1)N)\log_{2}n\leq L\lesssim\log n,$ $n^{1/(2\beta+1)}\lesssim\min_{i}p_{i}$ and $s\asymp n^{1/(2\beta+1)}\log n,$ Theorem 1 yields the rate

and the exponent in the rate does not depend on the input dimension $d.$

Suboptimality of wavelet series estimators

As argued before the composition assumption in (6) is very natural and generalizes many structural constraints such as additive models. In this section, we show that wavelet series estimators are unable to take advantage from the underlying composition structure in the regression function and achieve in some setups much slower convergence rates.

More specifically, we consider general additive models of the form $f_{0}(\mathbf{x})=h(x_{1}+\ldots+x_{d})$ with $h\in\mathcal{C}_{1}^{\alpha}([0,d],K).$ This can also be viewed as a special instance of the single index model, where the aim is not to estimate $h$ but $f_{0}.$ Using (13), the prediction error of neural network reconstructions with small empirical risk and depth $L\asymp\log n$ is then bounded by $n^{-2\alpha/(2\alpha+1)}\log^{3}n.$ The lower bound below shows that wavelet series estimators cannot converge faster than with the rate $n^{-2\alpha/(2\alpha+d)}.$ This rate can be much slower if $d$ is large. Wavelet series estimators thus suffer in this case from the curse of dimensionality while neural networks achieve fast rates.

To construct a counterexample, it is enough to consider the nonparametric regression model $Y_{i}=f_{0}(\mathbf{X}_{i})+\varepsilon_{i},$ $i=1,\ldots,n$ with uniform design $\mathbf{X}_{i}:=(U_{i,1},\ldots,U_{i,d})\sim\operatorname{Unif}^{d}.$ The empirical wavelet coefficients are

Because of $E[\widehat{d}_{\bm{\lambda}_{1}\ldots\bm{\lambda}_{d}}(f_{0})]=d_{\bm{\lambda}_{1}\ldots\bm{\lambda}_{d}}(f_{0}),$ this gives unbiased estimators for the wavelet coefficients. By the law of total variance,

For the lower bounds we may assume that the smoothness indices are known. For estimation, we can truncate the series expansion on a resolution level that balances squared bias and variance of the total estimator. More generally, we study estimators of the form

for an arbitrary subset $I\subset\Lambda\times\ldots\times\Lambda$. Using that the design is uniform,

Let $q$ be as above and set $\nu:=\lceil\log_{2}d\rceil+1.$ For any $0<\alpha\leq 1$ and any $K>0,$ there exists a non-zero constant $c(\psi,d)$ only depending on $d$ and properties of the wavelet function $\psi$ such that for any $j,$ we can find a function $f_{j,\alpha}(\mathbf{x})=h_{j,\alpha}(x_{1}+\ldots+x_{d})$ with $h_{j,\alpha}\in\mathcal{C}_{1}^{\alpha}([0,d],K)$ satisfying

for all $p_{1},\ldots,p_{d}\in\{0,1,\ldots,2^{j-q-\nu}-1\}.$

If $\widehat{f}_{n}$ denotes the wavelet estimator (15) for a compactly supported wavelet $\psi$ and an arbitrary index set $I,$ then, for any $0<\alpha\leq 1$ and any Hölder radius $K>0,$

A close inspection of the proof shows that the theorem even holds for $0<\alpha\leq r$ with $r$ the smallest positive integer for which $\int x^{r}\psi(x)dx\neq 0.$

A brief summary of related statistical theory for neural networks

This section is intended as a condensed overview on related literature summarizing main proving strategies for bounds on the statistical risk. An extended summary of the work until the late 90’s is given in . To control the stochastic error of neural networks, bounds on the covering entropy and VC dimension can be found in the monograph . A challenging part in the analysis of neural networks is the approximation theory for multivariate functions. We first recall results for shallow neural networks, that is, neural networks with one hidden layer.

Shallow neural networks: A shallow network with one output unit and width vector $(d,m,1)$ can be written as

The universal approximation theorem states that a neural network with one hidden layer can approximate any continuous function $f$ arbitrarily well with respect to the uniform norm provided there are enough hidden units, cf. . If $f$ has a derivative $f^{\prime}$, then the derivative of the neural network also approximates $f^{\prime}.$ The number of required hidden units might be, however, extremely large, cf. and . There are several proofs for the universal approximation theorem based on the Fourier transform, the Radon transform and the Hahn-Banach theorem .

The proofs can be sharpened in order to obtain rates of convergence. In the convergence rate $n^{-2\beta/(2\beta+d+5)}$ is derived. Compared with the minimax estimation rate this is suboptimal by a polynomial factor. The reason for the loss of performance with this approach is that rewriting the function as a network requires too many parameters.

In a similar strategy is used to derive the rate $C_{f}(d\tfrac{\log n}{n})^{1/2}$ for the squared $L^{2}$-risk, where $C_{f}:=\int|\omega|_{1}|\mathcal{F}f(\omega)|d\omega$ and $\mathcal{F}f$ denotes the Fourier transform of $f.$ If $C_{f}<\infty$ and $d$ is fixed the rate is always $n^{-1/2}$ up to logarithmic factors. Since $\sum_{i}\|\partial_{i}f\|_{\infty}\leq C_{f},$ this means that $C_{f}<\infty$ can only hold if $f$ has Hölder smoothness at least one. This rate is difficult to compare with the standard nonparametric rates except for the special case $d=1,$ where the rate is suboptimal compared with the minimax rate $n^{-2/(2+d)}$ for $d$-variate functions with smoothness one. More interestingly, the rate $C_{f}(d\tfrac{\log n}{n})^{1/2}$ shows that neural networks can converge fast if the underlying function satisfies some additional structural constraint. The same rate can also be obtained by a Fourier series estimator, see , Section 1.7. In a similar fashion, studies abstract function spaces on which shallow networks achieve fast convergence rates.

Results for multilayer neural networks: In it is shown how to approximate a polytope by a neural network with two hidden layers. Based on this result, uses two-layer neural networks with sigmoidal activation function and achieves the nonparametric rate $n^{-2\beta/(2\beta+d)}$ up to $\log n$-factors for $\beta\leq 1.$ This is extended in to a composition assumption and further generalized to $\beta>1$ in the recent article . Unfortunately, the result requires that the activation function is at least as smooth as the signal, cf. Theorem 1 in and therefore rules out the ReLU activation function.

The activation function $\sigma(x)=x^{2}$ is not of practical relevance but has some interesting theory. Indeed, with one hidden layer, we can generate quadratic polynomials and with $L$ hidden layers polynomials of degree $2^{L}.$ For this activation function, the role of the network depth is the polynomial degree and we can use standard results to approximate functions in common function classes. A natural generalization is the class of activation functions satisfying $\lim_{x\rightarrow-\infty}x^{-k}\sigma(x)=0$ and $\lim_{x\rightarrow+\infty}x^{-k}\sigma(x)=1.$

If the growth is at least quadratic ($k\geq 2$), the approximation theory has been derived in for deep networks with number of layers scaling with $\log d.$ The same class has also been considered recently in . For the approximations to work, the assumption $k\geq 2$ is crucial and the same approach does not generalize to the ReLU activation function, which satisfies the growth condition with $k=1$, and always produces functions that are piecewise linear in the input.

Approximation theory for the ReLU activation function has been only recently developed in . The key observation is that there are specific deep networks with few units which approximate the square function well. In particular, the function approximation presented in is essential for our approach and we use a similar strategy to construct networks that are close to a given function. We are, however, interested in a somehow different question. Instead of deriving existence of a network architecture with good approximation properties, we show that for any network architecture satisfying the conditions of Theorem 1 good approximation rates are obtainable. An additional difficulty in our approach is the boundedness of the network parameters.

Proofs

For the approximation of a function by a network, we first construct smaller networks computing simpler objects. Let $\mathbf{p}=(p_{0},\ldots,p_{L+1})$ and $\mathbf{p}^{\prime}=(p_{0}^{\prime},\ldots,p_{L+1}^{\prime}).$ To combine networks, we make frequently use of the following rules.

Enlarging: $\mathcal{F}(L,\mathbf{p},s)\subseteq\mathcal{F}(L,\mathbf{q},s^{\prime})$ whenever $\mathbf{p}\leq\mathbf{q}$ componentwise and $s\leq s^{\prime}.$

Additional layers/depth synchronization: To synchronize the number of hidden layers for two networks, we can add additional layers with identity weight matrix, such that

Parallelization: Suppose that $f,g$ are two networks with the same number of hidden layers and the same input dimension, that is, $f\in\mathcal{F}(L,\mathbf{p})$ and $g\in\mathcal{F}(L,\mathbf{p}^{\prime})$ with $p_{0}=p_{0}^{\prime}.$ The parallelized network $(f,g)$ computes $f$ and $g$ simultaneously in a joint network in the class $\mathcal{F}(L,(p_{0},p_{1}+p_{1}^{\prime},\ldots,p_{L+1}+p_{L+1}^{\prime})).$

To see this, let $f(\mathbf{x})=W_{L}\sigma_{\mathbf{v}_{L}}W_{L-1}\ldots\sigma_{\mathbf{v}_{1}}W_{0}\mathbf{x}\in\mathcal{F}(L,\mathbf{p},s).$ If all entries of the $j$-th column of $W_{i}$ are zero, then we can remove this column together with the $j$-th row in $W_{i-1}$ and the $j$-th entry of $\mathbf{v}_{i}$ without changing the function. This shows then that $f\in\mathcal{F}(L,(p_{0},\ldots,p_{i-1},p_{i}-1,p_{i+1},\ldots,p_{L+1}),s).$ Because there are $s$ active parameters, we can iterate this procedure at least $p_{i}-s$ times for any $i=1,\ldots,L.$ This proves f\in\mathcal{F}\big{(}L,(p_{0},p_{1}\wedge s,p_{2}\wedge s,\ldots,p_{L}\wedge s,p_{L+1}),s\big{)}.

For any function $f\in\mathcal{C}_{r}^{\beta}(^{r},K)$ and any integers $m\geq 1$ and $N\geq(\beta+1)^{r}\vee(K+1)e^{r},$ there exists a network

The proof of the theorem is given in the supplement. The idea is to first build networks that for given input $(x,y)$ approximately compute the product $xy.$ We then split the input space into small hyper-cubes and construct a network that approximates a local Taylor expansion on each of these hyper-cubes.

Based on Theorem 5, we can now construct a network that approximates $f=g_{q}\circ\ldots\circ g_{0}.$ In a first step, we show that $f$ can always be written as composition of functions defined on hypercubes $^{t_{i}}.$ As in the previous theorem, let $g_{ij}\in\mathcal{C}_{t_{i}}^{\beta_{i}}([a_{i},b_{i}]^{t_{i}},K_{i})$ and assume that $K_{i}\geq 1.$ For $i=1,\ldots,q-1,$ define

Here, $2K_{i-1}\mathbf{x}-K_{i-1}$ means that we transform the entries by $2K_{i-1}x_{j}-K_{i-1}$ for all $j.$ Clearly,

Using the definition of the Hölder balls $C_{r}^{\beta}(D,K),$ it follows that $h_{0j}$ takes values in $,$ $h_{0j}\in\mathcal{C}_{t_{0}}^{\beta_{0}}(^{t_{0}},1),$ $h_{ij}\in\mathcal{C}_{t_{i}}^{\beta_{i}}(^{t_{i}},(2K_{i-1})^{\beta_{i}})$ for $i=1,\ldots,q-1,$ and $h_{qj}\in\mathcal{C}_{t_{q}}^{\beta_{q}}(^{t_{q}},K_{q}(2K_{q-1})^{\beta_{q}}).$ Without loss of generality, we can always assume that the radii of the Hölder balls are at least one, that is, $K_{i}\geq 1.$

Let $h_{ij}$ be as above with $K_{i}\geq 1$. Then, for any functions $\widetilde{h}_{i}=(\widetilde{h}_{ij})_{j}^{\top}$ with $\widetilde{h}_{ij}:^{t_{i}}\rightarrow,$

Define $H_{i}=h_{i}\circ\ldots\circ h_{0}$ and $\widetilde{H}_{i}=\widetilde{h}_{i}\circ\ldots\circ\widetilde{h}_{0}.$ If $Q_{i}$ is an upper bound for the Hölder semi-norm of $h_{ij},$ $j=1,\ldots,d_{i+1},$ we find using triangle inequality,

Together with the inequality $(y+z)^{\alpha}\leq y^{\alpha}+z^{\alpha}$ which holds for all $y,z\geq 0$ and all $\alpha\in,$ the result follows. ∎

It is enough to prove the result for all sufficiently large $n.$ Throughout the proof $C^{\prime}$ is a constant only depending on $(q,\mathbf{d},\mathbf{t},\bm{\beta},F)$ that may change from line to line. Combining Theorem 2 with the assumed bounds on the depth $L$ and the network sparsity $s,$ it follows for $n\geq 3,$

where we used $\varepsilon=1/2$ for the lower bound and $\varepsilon=1$ for the upper bound. Taking $C=8C^{\prime},$ we find that $\tfrac{1}{8}\Delta_{n}(\widehat{f}_{n},f_{0})\leq R(\widehat{f},f_{0})$ whenever $\Delta_{n}(\widehat{f}_{n},f_{0})\geq C\phi_{n}L\log^{2}n.$ This proves the lower bound in (9).

To derive the upper bounds in (8) and (9) we need to bound the approximation error. To do this, we rewrite the regression function $f_{0}$ as in (21), that is, $f_{0}=h_{q}\circ\ldots h_{0}$ with $h_{i}=(h_{ij})_{j}^{\top},$ $h_{ij}$ defined on $^{t_{i}},$ and for any $i<q,$ $h_{ij}$ mapping to $.$

We apply Theorem 5 to each function $h_{ij}$ separately. Take $m=\lceil\log_{2}n\rceil$ and let $L^{\prime}_{i}:=8+(\lceil\log_{2}n\rceil+5)(1+\lceil\log_{2}(t_{i}\vee\beta_{i})\rceil).$ This means there exists a network $\widetilde{h}_{ij}\in\mathcal{F}(L^{\prime}_{i},(t_{i},6(t_{i}+\lceil\beta_{i}\rceil)N,\ldots,6(t_{i}+\lceil\beta_{i}\rceil)N,1),s_{i})$ with $s_{i}\leq 141(t_{i}+\beta_{i}+1)^{3+t_{i}}N(\lceil\log_{2}n\rceil+6),$ such that

where $Q_{i}$ is any upper bound of the Hölder norms of $h_{ij}.$ If $i<q,$ then we apply to the output the two additional layers $1-(1-x)_{+}.$ This requires four additional parameters. Call the resulting network $h_{ij}^{*}\in\mathcal{F}(L^{\prime}_{i}+2,(t_{i},6(t_{i}+\lceil\beta_{i}\rceil)N,\ldots,6(t_{i}+\lceil\beta_{i}\rceil)N,1),s_{i}+4)$ and observe that $\sigma(h_{ij}^{*})=(\widetilde{h}_{ij}(x)\vee 0)\wedge 1.$ Since $h_{ij}(\mathbf{x})\in,$ we must have

If the networks $h_{ij}^{*}$ are computed in parallel, $h_{i}^{*}=(h_{ij}^{*})_{j=1,\ldots,d_{i+1}}$ lies in the class

with $r_{i}:=\max_{i}d_{i+1}(t_{i}+\lceil\beta_{i}\rceil).$ Finally, we construct the composite network $f^{*}=\widetilde{h}_{q1}\circ\sigma(h_{q-1}^{*})\circ\ldots\circ\sigma(h_{0}^{*})$ which by the composition rule in Section 7.1 can be realized in the class

with $E:=3(q-1)+\sum_{i=0}^{q}L^{\prime}_{i}.$ Observe that there is an $A_{n}$ that is bounded in $n$ such that $E=A_{n}+\log_{2}n(\sum_{i=0}^{q}\lceil\log_{2}(t_{i}\vee\beta_{i})\rceil+1).$ Using that $\lceil x\rceil<x+1,$ $E\leq\sum_{i=0}^{q}(\log_{2}(4)+\log_{2}(t_{i}\vee\beta_{i}))\log_{2}n\leq L$ for all sufficiently large $n.$ By (18) and for sufficiently large $n,$ the space (25) can be embedded into $\mathcal{F}(L,\mathbf{p},s)$ with $L,\mathbf{p},s$ satisfying the assumptions of the theorem by choosing $N=\lceil c\max_{i=0,\ldots,q}n^{\frac{t_{i}}{2\beta_{i}^{*}+t_{i}}}\rceil$ for a sufficiently small constant $c>0$ only depending on $q,\mathbf{d},\mathbf{t},\bm{\beta}.$ Combining Lemma 3 with (23) and (24)

For the approximation error in (22) we need to find a network function that is bounded in sup-norm by $F.$ By the previous inequality there exists a sequence $(\widetilde{f}_{n})_{n}$ such that for all sufficiently large $n,$ $\widetilde{f}_{n}\in\mathcal{F}(L,\mathbf{p},s)$ and $\|\widetilde{f}_{n}-f_{0}\|_{\infty}^{2}\leq 2C\max_{i=0,\ldots,q}c^{-2\beta_{i}^{*}/t_{i}}n^{-(2\beta_{i}^{*})/(2\beta_{i}^{*}+t_{i})}.$ Define $f^{*}_{n}=\widetilde{f}_{n}(\|f_{0}\|_{\infty}/\|\widetilde{f}_{n}\|_{\infty}\newline \wedge 1).$ Then, $\|f^{*}_{n}\|_{\infty}\leq\|f_{0}\|_{\infty}=\|g_{q}\|_{\infty}\leq K\leq F,$ where the last inequality follows from assumption $(i).$ Moreover, $f^{*}_{n}\in\mathcal{F}(L,\mathbf{p},s,F).$ Writing $f_{n}^{*}-f_{0}=(f_{n}^{*}-\widetilde{f}_{n})+(\widetilde{f}_{n}-f_{0}),$ we obtain $\|f_{n}^{*}-f_{0}\|_{\infty}\leq 2\|\widetilde{f}_{n}-f_{0}\|_{\infty}.$ This shows that (26) also holds (with constants multiplied by $8$) if the infimum is taken over the smaller space $\mathcal{F}(L,\mathbf{p},s,F).$ Together with (22) the upper bounds in (8) and (9) follow for any constant $C$. This completes the proof. ∎

2 Proof of Theorem 2

Several oracle inequalities for the least-squares estimator are know, cf. . The common assumption of bounded response variables is, however, violated in the nonparametric regression model with Gaussian measurement noise. Additionally we provide also a lower bound of the risk and give a proof that can be easily generalized to arbitrary noise distributions. Let $\mathcal{N}(\delta,\mathcal{F},\|\cdot\|_{\infty})$ be the covering number, that is, the minimal number of $\|\cdot\|_{\infty}$-balls with radius $\delta$ that covers $\mathcal{F}$ (the centers do not need to be in $\mathcal{F}$).

Consider the $d$-variate nonparametric regression model (1) with unknown regression function $f_{0}.$ Let $\widehat{f}$ be any estimator taking values in $\mathcal{F}.$ Define

and assume $\{f_{0}\}\cup\mathcal{F}\subset\{f:^{d}\rightarrow[-F,F]\}$ for some $F\geq 1.$ If $\mathcal{N}_{n}:=\mathcal{N}(\delta,\mathcal{F},\|\cdot\|_{\infty})\geq 3,$ then,

The proof of the lemma can be found in the supplement. Next, we prove a covering entropy bound. Recall the definition of the network function class $\mathcal{F}(L,\mathbf{p},s,F)$ in (4).

For a proof see the supplement. A related result is Theorem 14.5 in .

The assertion follows from Lemma 5 with $\delta=1/n,$ Lemma 4 and Remark 1 since $F\geq 1.$ ∎

3 Proof of Theorem 3

Throughout this proof, \|\cdot\|_{2}=\|\cdot\|_{L^{2}^{d}}. By assumption there exist positive $\gamma\leq\Gamma$ such that the Lebesgue density of $\mathbf{X}$ is lower bounded by $\gamma$ and upper bounded by $\Gamma$ on $^{d}.$ For such design, $R(\widehat{f}_{n},f_{0})\geq\gamma\|\widehat{f}_{n}-f_{0}\|_{2}^{2}.$ Denote by $P_{f}$ the law of the data in the nonparametric regression model (1). For the Kullback-Leibler divergence we have $\operatorname{KL}(P_{f},P_{g})=nE[(f(\mathbf{X}_{1})-g(\mathbf{X}_{1}))^{2}]\leq\Gamma n\|f-g\|_{2}^{2}.$ Theorem 2.7 in states that if for some $M\geq 1$ and $\kappa>0,$ $f_{(0)},\ldots,f_{(M)}\in\mathcal{G}(q,\mathbf{d},\mathbf{t},\bm{\beta},K)$ are such that

$\|f_{(j)}-f_{(k)}\|_{2}\geq\kappa\sqrt{\phi_{n}}$ for all $0\leq j<k\leq M$

$n\sum_{j=1}^{M}\|f_{(j)}-f_{(0)}\|_{2}^{2}\leq M\log(M)/(9\Gamma),$

then there exists a positive constant $c=c(\kappa,\gamma),$ such that

Hence $\psi_{\mathbf{u}}\in\mathcal{C}_{t^{*}}^{\beta^{*}}(^{t^{*}},(\beta^{*})^{t^{*}}t^{*}).$ For a vector $\mathbf{w}=(w_{\mathbf{u}})_{\mathbf{u}\in\mathcal{U}_{n}}\in\{0,1\}^{|\mathcal{U}_{n}|},$ define

By construction, $\psi_{\mathbf{u}}$ and $\psi_{\mathbf{u}^{\prime}},$ $\mathbf{u},\mathbf{u}^{\prime}\in\mathcal{U}_{n},$ $\mathbf{u}\neq\mathbf{u}^{\prime}$ have disjoint support. As a consequence $\phi_{\mathbf{w}}\in\mathcal{C}_{t^{*}}^{\beta^{*}}(^{t^{*}},2(\beta^{*})^{t^{*}}t^{*}).$

and $f_{\mathbf{w}}\in\mathcal{G}(q,\mathbf{d},\mathbf{t},\bm{\beta},K)$ provided $K$ is taken sufficiently large.

For all $\mathbf{u},$ $\|\psi_{\mathbf{u}}\|_{2}^{2}=h_{n}^{2\beta^{**}+t^{*}}\|K^{B}\|_{2}^{2t^{*}}.$ If $\operatorname{Ham}(\mathbf{w},\mathbf{w}^{\prime})=\sum_{\mathbf{u}\in\mathcal{U}_{n}}\mathbf{1}(w_{\mathbf{u}}\neq w_{\mathbf{u}^{\prime}})$ denotes the Hamming distance, we find

By the Varshamov - Gilbert bound (, Lemma 2.9) and because of $m_{n}^{t^{*}}=|\mathcal{U}_{n}|,$ we conclude that there exists a subset $\mathcal{W}\subset\{0,1\}^{m_{n}^{t^{*}}}$ of cardinality $|\mathcal{W}|\geq 2^{m_{n}^{t^{*}}/8},$ such that $\operatorname{Ham}(\mathbf{w},\mathbf{w}^{\prime})\geq m_{n}^{t^{*}}/8$ for all $\mathbf{w},\mathbf{w}^{\prime}\in\mathcal{W},$ $\mathbf{w}\neq\mathbf{w}^{\prime}.$ This implies that for $\kappa=\|K^{B}\|_{2}^{t^{*}}/(\sqrt{8}\rho^{\beta^{**}}),$

Using the definition of $h_{n}$ and $\rho,$

This shows that the functions $f_{\mathbf{w}}$ with $\mathbf{w}\in\mathcal{W}$ satisfy $(i)$ and $(ii).$ The assertion follows. ∎

4 Proof of Lemma 1

We will choose $c_{2}\leq 1.$ Since $\|f_{0}\|_{\infty}\leq K,$ it is therefore enough to consider the infimum over $\mathcal{F}(L,\mathbf{p},s,F)$ with $F=K+1.$ Let $\widetilde{f}_{n}$ be an empirical risk minimizer. Recall that $\Delta_{n}(\widetilde{f}_{n},f_{0})=0.$ Because of the minimax lower bound in Theorem 3, there exists a constant $c_{3}$ such that $c_{3}n^{-2\beta/(2\beta+d)}\leq\sup_{f_{0}\in\mathcal{C}_{1}^{\beta}(,K)}R(\widetilde{f}_{n},f_{0})$ for all sufficiently large $n.$ Because of $p_{0}=d$ and $p_{L+1}=1,$ Theorem 2 yields

for some constant $C.$ Given $\varepsilon,$ set $n_{\varepsilon}:=\lfloor(\sqrt{8}\varepsilon/\sqrt{c_{3}})^{-(2\beta+d)/\beta}\rfloor.$ Observe that for $\varepsilon\leq\sqrt{c_{3}/8},$ $n_{\varepsilon}^{-1}\leq 2(\sqrt{8}\varepsilon/\sqrt{c_{3}})^{(2\beta+d)/\beta}$ and $8\varepsilon^{2}/c_{3}\leq n_{\varepsilon}^{-2\beta/(2\beta+d)}.$ For sufficiently small $c_{2}>0$ and all $\varepsilon\leq c_{2},$ we can insert $n_{\varepsilon}$ in the previous inequality and find

for constants $C_{1},C_{2}$ depending on $K,\beta,$ and $d.$ The result follows using the condition $s\leq c_{1}\varepsilon^{-d/\beta}/(L\log(1/\varepsilon))$ and choosing $c_{1}$ small enough. ∎

5 Proofs for Section 5

For a real number $u,$ denote by $\{u\}$ the fractional part of $u.$

We need to study the cases $\mu_{0}\neq 0$ and $\mu_{0}=0$ separately. If $\mu_{0}\neq 0,$ define $g(u)=r^{-1}u^{r}\mathbf{1}_{[0,1/2]}(u)+r^{-1}(1-u)^{r}\mathbf{1}_{(1/2,1]}(u).$ Notice that $g$ is Lipschitz with Lipschitz constant one. Let $h_{j,\alpha}(u)=K2^{-j\alpha-1}g(\{2^{j-q-\nu}u\})$ with $q$ and $\nu$ as defined in the statement of the lemma. For a $T$-periodic function $u\mapsto s(u)$ the $\alpha$-Hölder semi-norm for $\alpha\leq 1$ can be shown to be $\sup_{u\neq v,|u-v|\leq T}|s(u)-s(v)|/|u-v|^{\alpha}.$ Since $g$ is $1$-Lipschitz, we have for $u,v$ with $|u-v|\leq 2^{q+\nu-j},$