Random Vector Functional Link Networks for Function Approximation on Manifolds

The learning speed of feed-forward neural networks is notoriously slow and has presented a bottleneck in deep learning applications for several decades. For instance, gradient-based learning algorithms, which are used extensively to train neural networks, tend to work slowly when all of the network parameters must be iteratively tuned. To counter this, both researchers and practitioners have tried introducing randomness to reduce the learning requirement. Based on the original construction of Igelnik and Pao, single layer neural-networks with random input-to-hidden layer weights and biases have seen success in practice, but the necessary theoretical justification is lacking. In this paper, we begin to fill this theoretical gap. We provide a (corrected) rigorous proof that the Igelnik and Pao construction is a universal approximator for continuous functions on compact domains, with approximation error decaying asymptotically like $O(1/\sqrt{n})$ for the number $n$ of network nodes. We then extend this result to the non-asymptotic setting, proving that one can achieve any desired approximation error with high probability provided $n$ is sufficiently large. We further adapt this randomized neural network architecture to approximate functions on smooth, compact submanifolds of Euclidean space, providing theoretical guarantees in both the asymptotic and non-asymptotic forms. Finally, we illustrate our results on manifolds with numerical experiments.


INTRODUCTION
In recent years, neural networks have once again triggered an increased interest among researchers in the machine learning community.So-called deep neural networks model functions using a composition of multiple hidden layers, each transforming (possibly non-linearly) the previous layer before building a final output representation, see [18,36,9,12,43].In machine learning parlance, these layers are determined by sets of weights and biases that can be tuned so that the network mimics the action of a complex function.In particular, a single layer feed-forward neural network (SLFN) with n nodes may be regarded as a parametric function f n : R N → R of the form Here, the function ρ : R → R is called an activation function and is potentially nonlinear.Some typical examples include the sigmoid function ρ(z) = 1 1+exp(−z) , ReLU ρ(z) = max{0, z}, and sign functions, among many others.The parameters of the SLFN are the number of nodes n ∈ N in the the hidden layer, the input-to-hidden layer weights and biases {w k } n k=1 ⊂ R N and {b k } n k=1 ⊂ R (resp.), and the hidden-to-output layer weights {v k } n k=1 ⊂ R. In this way, neural networks are fundamentally parametric families of functions whose parameters may be chosen to approximate a given function.
It has been shown that any compactly supported continuous function can be approximated with any given precision by a single layer neural network with a suitably chosen number of nodes [3], and harmonic analysis techniques have been used to study stability of such approximations [5].Other recent results that take a different approach directly analyze the capacity of neural networks from a combinatorial point of view [41,2].
While these results ensure existence of a neural network approximating a function, practical applications require construction of such an approximation.The parameters of the neural network can be chosen using optimization techniques to minimize the difference between the network and the function f : R N → R it is intended to model.In practice, the function f is usually not known, and we only have access to a set {(x k , f (x k ))} m k=1 of values of the function at finitely many points sampled from its domain, called a training set.The approximation error can be measured by comparing the training data to the corresponding network outputs when evaluated on the same set of points, and the parameters of the neural network f n can be learned by minimizing a given loss function L(x 1 , . . ., x m ); a typical loss function is the sum-of-squares error The SLFN which approximates f is then determined using an optimization algorithm, such as backpropagation, to find the network parameters which minimize L(x 1 , . . ., x m ).It is known that there exist weights and biases which make the loss function vanish when the number of nodes n is at least m, provided the activation function is bounded, nonlinear, and has at least one finite limit at either ±∞ [13].
Unfortunately, optimizing the parameters in SLFNs can be difficult.For instance, any non-linearity in the activation function can cause back-propagation to be very time consuming or get caught in local minima of the loss function [35].Moreover, deep neural networks can require massive amounts of training data, and so are typically unreliable for applications with very limited data availability, such as agriculture, healthcare, and ecology [26].
To address some of the difficulties associated with training deep neural networks, both researchers and practitioners have attempted to incorporate randomness in some way.Indeed, randomization-based neural networks that yield closed form solutions typically require less time to train and avoid some of the pitfalls of traditional neural networks trained using back-propagation [35,32,39].One of the popular randomizationbased neural network architectures is the Random Vector Functional Link (RVFL) network [29,14], which is a single layer feed-forward neural network in which the input-to-hidden layer weights and biases are selected randomly and independently from a suitable domain and the remaining hidden-to-output layer weights are learned using training data.
By eliminating the need to optimize the input-to-hidden layer weights and biases, RVFL networks turn supervised learning into a purely linear problem.To see this, define ρ(X) ∈ R n×m to be the matrix whose jth column is {ρ(⟨w k , x j ⟩ + b k )} n k=1 and f (X) ∈ R m the vector whose jth entry is f (x j ).Then the vector v ∈ R n of hidden-to-output layer weights is the solution to the matrix-vector equation f (X) = ρ(X) T v, which can be solved by computing the Moore-Penrose pseudoinverse of ρ(X) T .
Even though RVFL networks are proving their usefulness in practice, the supporting theoretical framework is currently lacking [ see 46].Most theoretical research into the approximation capabilities of deep neural networks centers around two main concepts: universal approximation of functions on compact domains and point-wise approximation on finite training sets.For instance, in the early 1990s it was shown that multi-layer feed-forward neural networks having activation functions that are continuous, bounded, and non-constant are universal approximators (in the L p sense for 1 ≤ p < ∞) of continuous functions on compact domains [11,20].The most notable result in the existing literature regarding the universal approximation capability of RVFL networks is due to B. Igelnik and Y.H. Pao in the mid-1990s, who showed that such neural networks can universally approximate continuous functions on compact sets [14]; the noticeable lack of results since has left a sizable gap between theory and practice.In this paper, we begin to bridge this gap by further improving the Igelnik and Pao result, and bringing the mathematical theory behind RFVL networks into the modern spotlight.Below, we introduce the notation that will be used throughout this paper, and describe our main contributions.

Notation
For a function f : R N → R, the set supp(f ) ⊂ R N denotes the support of f .We denote by C c (R N ) and C 0 (R N ) the classes of continuous functions mapping R N to R whose support sets are compact and vanish at infinity, respectively.Given a set S ⊂ R N , we define its radius to be rad(S) := sup x∈S ∥x∥ 2 ; moreover, if dµ denotes the uniform volume measure on S, then we write vol(S) := S dµ to represent the volume of S. For any probability distribution P : R N → [0, 1], a random variable X distributed according to P is denoted by X ∼ P , and we write its expectation as EX := R N XdP .The open ℓ p ball of radius r > 0 centered at x ∈ R N is denoted by B N p (x, r) for all 1 ≤ p ≤ ∞; the ℓ p unit-ball centered at the origin is abbreviated B N p .Given a fixed δ > 0 and a set S ⊂ R N , a minimal δ-net for S, which we denote C(δ, S), is the smallest subset of S satisfying S ⊂ ∪ x∈C(δ,S) B N 2 (x, δ); the δ-covering number of S is the cardinality of a minimal δ-net for S and is denoted N (δ, S) := |C(δ, S)|.

Main results
In this paper, we study the uniform approximation capabilities of RVFL networks.More specifically, we consider the problem of using RVFL networks to estimate a continuous, compactly supported function on N -dimensional Euclidean space.
The first theoretical result on approximating properties of RVFL networks, due to Igelnik and Pao, guarantees that continuous functions may be universally approximated on compact sets using RVFL networks, provided the number of nodes n ∈ N in the network goes to infinity [14].Moreover, it shows that the mean square error of the approximation vanishes at a rate proportional to 1/n.At the time, this result was state-of-the-art and justified how RVFL networks were used in practice.However, the original theorem is not technically correct.In fact, several aspects of the proof technique are flawed.Some of the minor flaws are mentioned in [21], but the subsequent revisions do not address the more significant issues which would make the statement of the result technically correct.We address these issues in this paper, see Remark 1.Thus, our first contribution to the theory of RVFL networks is a corrected version of the original Igelnik and Pao theorem: For any ε > 0, there exist distributions from which input weights {w k } n k=1 and biases {b k } n k=1 are drawn, and there exist hidden-to-output layer weights {v k } n k=1 ⊂ R that depend on the realization of weights and biases, such that the sequence of RVFL networks {f n } ∞ n=1 defined by as n → ∞.
For a more precise formulation of Theorem 1 and its proof, we refer the reader to Theorem 5 and Section 3.1.REMARK 1.
1.Even though in Theorem 1 we only claim existence of the distribution for input weights {w k } n k=1 and biases {b k } n k=1 , such a distribution is actually constructed in the proof.Namely, for any ε > 0, there exist constants α, Ω > 0 such that the random variables are independently drawn from their associated distributions, and b k := −⟨w k , y k ⟩ − αu k .
2.We note that, unlike the original theorem statement in [14], Theorem 1 does not show exact convergence of the sequence of constructed RVFL networks f n to the original function f .Indeed, it only ensures that the limit f n is ε-close to f .This should still be sufficient for practical applications since, given a desired accuracy level ε > 0, one can find values of α, Ω, n such that this accuracy level is achieved on average.Exact convergence can be proved if one replaces α and Ω in the distribution described above by sequences {α n } ∞ n=1 and {Ω n } ∞ n=1 of positive numbers, both tending to infinity with n.In this setting, however, there is no guaranteed rate of convergence; moreover, as n increases, the ranges of the random variables {w k } n k=1 and {u k } n k=1 become increasingly larger, which may cause problems in practical applications.
3.The approach we take to construct the RVFL network approximating a function f allows one to compute the output weights {v k } n k=1 exactly (once the realization of random parameters is fixed), in the case where the function f is known.For the details, we refer the reader to equations (6) and (8) in the proof of Theorem 1.If we only have access to a training set that is sufficiently large and uniformly distributed over the support of f , these formulas can be used to compute the output weights approximately, instead of solving the least squares problem.
4.Note that the normalization R ρ(z)dz = 1 of the activation function can be replaced by the condition R ρ(z)dz ̸ = 0. Indeed, in the case when ρ ∈ L 1 (R) ∩ L ∞ (R) and R ρ(z)dz / ∈ {0, 1}, one can simply use Theorem 1 to approximate 1 R ρ(z)dz f by a sequence of RVFL network with the activation function 1 R ρ(z)dz ρ.Mutatis mutandis in the case when R ρ ′ (z)dz ′ / ∈ {0, 1}.More generally, this trick allows any of our theorems to be applied in the case R ρ(z)dz ̸ = 0.
One of the drawbacks of Theorem 1 is that the mean square error guarantee is asymptotic in the number of nodes used in the neural network.This is clearly impractical for applications, and so it is desirable to have a more explicit error bound for each fixed number n of nodes used.To this end, we provide a new, non-asymptotic version of Theorem 1, which provides an error guarantee with high probability whenever the number of network nodes is large enough, albeit at the price of an additional Lipschitz requirement on the activation function: Suppose further that ρ is κ-Lipschitz on R for some κ > 0. For any ε > 0 and η ∈ (0, 1), suppose that n ≥ C(N, f, ρ)ε −1 log(η −1 /ε), where C(N, f, ρ) is independent of ε and η and depends on f , ρ, and superexponentially on N .Then there exist distributions from which input weights {w k } n k=1 and biases {b k } n k=1 are drawn, and there exist hidden-to-output layer weights {v k } n k=1 ⊂ R that depend on the realization of weights and biases, such that the RVFL network defined by For simplicity, the bound on the number n of the nodes on the hidden layer here is rough.For a more precise formulation of this result that contains a bound with explicit constant, we refer the reader to Theorem 6 in Section 3.2.We also note that the distribution of the input weight and bias here can be selected as described in Remark 1.
The constructions of RVFL networks presented in Theorems 1 and 2 depend heavily on the dimension of the ambient space R N .If N is small, this dependence does not present much of a problem.However, many modern applications require the ambient dimension to be large.Fortunately, a common assumption in practice is that support of the signals of interest lie on a lower-dimensional manifold embedded in R N .For instance, the landscape of cancer cell states can be modeled using nonlinear, locally continuous "cellular manifolds;" indeed, while the ambient dimension of this state space is typically high (e.g., single-cell RNA sequencing must account for approximately 20,000 gene dimensions), cellular data actually occupies an intrinsically lower dimensional space [4].Likewise, while the pattern space of neural population activity in the brain is described by an exponential number of parameters, the spatiotemporal dynamics of brain activity lie on a lower-dimensional subspace or "neural manifold" [25].In this paper, we propose a new RVFL network architecture for approximating continuous functions defined on smooth compact manifolds that allows to replace the dependence on the ambient dimension N with dependence on the manifold intrinsic dimension.We show that RVFL approximation results can be extended to this setting.More precisely, we prove the following analog of Theorem 2. THEOREM 3. Let M ⊂ R N be a smooth, compact d-dimensional manifold with finite atlas {(U j , ϕ j )} j∈J and f ∈ C(M).Fix any activation function ρ ∈ L 1 (R) ∩ L ∞ (R) with R ρ(z)dz = 1 such that ρ is κ-Lipschitz on R for some κ > 0. For any ε > 0 and η ∈ (0, 1), suppose n ≥ C(d, f, ρ)ε −1 log(η −1 /ε), where C(d, f, ρ) is independent of ε and η and depends on f , ρ, and superexponentially on d.Then there exists an RVFL-like approximation f n of the function f with a parameter selection similar to the Theorem 1 construction that satisfies For a the construction of the RVFL-like approximation f n and a more precise formulation of this result and an analog of Theorem 1 applied to manifolds, we refer the reader to Section 3.3.1 and Theorems 7 and 8.We note that the approximation f n here is not obtained as a single RVFL network construction, but rather as a combination of several RVFL networks in local manifold coordinates.

Organization
The remaining part of the paper is organized as follows.In Section 2, we discuss some theoretical preliminaries on concentration bounds for Monte-Carlo integration and on smooth compact manifolds.Monte-Carlo integration is an essential ingredient in our construction of RVFL networks approximating a given function, and we use the results listed in this section to establish approximation error bounds.Theorem 1 is proven in Section 3.1, where we break down the proof into four main steps, constructing a limit-integral representation of the function to be approximated in Lemmas 3 and 4, then using Monte-Carlo approximation of the obtained integral to construct an RVFL network in Lemma 5, and, finally, establishing approximation guarantees for the constructed RVFL network.The proofs of Lemmas 3, 4, and 5 can be found in Sections 3.5.1,3.5.2,and 3.5.3,respectively.We further study properties of the constructed RVFL networks and prove the non-asymptotic approximation result of Theorem 2 in Section 3.2.In Section 3.3, we generalize our results and propose a new RVFL network architecture for approximating continuous functions defined on smooth compact manifolds.We show that RVFL approximation results can be extended to this setting by proving an analog of Theorem 1 in Section 3.3.2and Theorem 3 in Section 3.5.5.Finally, in Section 3.4, we provide numerical evidence to illustrate the result of Theorem 3.

MATERIALS AND METHODS
In this section, we briefly introduce supporting material and theoretical results which we will need in later sections.This material is far from exhaustive, and is meant to be a survey of definitions, concepts, and key results.

A concentration bound for classic Monte-Carlo integration
A crucial piece of the proof technique employed in [14], which we will use repeatedly, is the use of the Monte-Carlo method to approximate high-dimensional integrals.As such, we start with the background on Monte-Carlo integration.The following introduction is adapted from the material in [8].
Let f : R N → R and S ⊂ R N a compact set.Suppose we want to estimate the integral I(f, S) := S f dµ, where µ is the uniform measure on S. The classic Monte Carlo method does this by an equal-weight cubature rule, where {x j } n j=1 are independent identically distributed uniform random samples from S and vol(S) := S dµ is the volume of S. In particular, note that EI n (f, S) = I(f, S) and Let us define the quantity It follows that the random variable I n (f ) has mean I(f, S) and variance vol 2 (S)σ(f, S) 2 /n.Hence, by the Central Limit Theorem, provided that 0 < vol 2 (S)σ(f, S) 2 < ∞, we have for any constant C > 0, where ε(f, S) := vol(S)σ(f, S).This yields the following well-known result: THEOREM 4. For any f ∈ L 2 (S, µ), the mean-square error of the Monte Carlo approximation I n (f, S) satisfies where the expectation is taken with respect to the random variables {x j } n j=1 and σ(f, S) is defined in (1).
In particular, Theorem 4 implies In the non-asymptotic setting, we are interested in obtaining a useful bound on the probability P(|I n (f, S) − I(f, S)| ≥ t) for all t > 0. The following lemma follows from a generalization of Bennett's inequality (Theorem 7.6 in [19]; see also [24,37]).LEMMA 1.For any f ∈ L 2 (S) and n ∈ N we have for all t > 0 and a universal constant C > 0, provided |vol(S)f (x)| ≤ K for almost every x ∈ S.

Smooth, compact manifolds in Euclidean space
In this section we review several concepts of smooth manifolds that will be useful to us later.Many of the definitions and results that follow can be found, for instance, in [33].Let M ⊂ R N be a smooth, compact d-dimensional manifold.A chart for M is a pair (U, ϕ) such that U ⊂ M is an open set and ϕ : U → R d is a homeomorphism.One way to interpret a chart is as a tangent space at some point x ∈ U ; in this way, a chart defines a Euclidean coordinate system on U via the map ϕ.A collection {(U j , ϕ j )} j∈J of charts defines an atlas for M if ∪ j∈J U j = M.We now define a special collection of functions on M called a partition of unity.DEFINITION 1.Let M ⊂ R N be a smooth manifold.A partition of unity of M with respect to an open cover {U j } j∈J of M is a family of nonnegative smooth functions {η j } j∈J such that for every x ∈ M we have 1 = j∈J η j (x) and, for every j ∈ J, supp(η j ) ⊂ U j .
It is known that if M is compact there exists a partition of unity of M such that supp(η j ) is compact for all j ∈ J [see 40].In particular, such a partition of unity exists for any open cover of M corresponding to an atlas.
Fix an atlas {(U j , ϕ j )} j∈J for M, as well as the corresponding, compactly supported partition of unity {η j } j∈J .Then we have the following, useful result [see 33, Lemma 4.8].
LEMMA 2. Let M ⊂ R N be a smooth, compact manifold with atlas {(U j , ϕ j )} j∈J and compactly supported partition of unity {η j } j∈J .For any f ∈ C(M) we have In later sections, we use the representation of Lemma 2 to integrate functions f ∈ C(M) over M. To this end, for each j ∈ J, let Dϕ j (y) denote the differential of ϕ j at y ∈ U j , which is a map from the tangent space T y M into R d .One may interpret Dϕ j (y) as the matrix representation of a basis for the cotangent space at y ∈ U j .As a result, Dϕ j (y) is necessarily invertible for each y ∈ U j , and so we know that | det(Dϕ j (y))| > 0 for each y ∈ U j .Hence, it follows by the change of variables theorem that

RESULTS
In this section, we prove our main results formulated in Section 1.2 and also use numerical simulations to illustrate the RVFL approximation performance in a low-dimensional submanifold setup.To improve readability of this section, we postpone the proofs of technical lemmas till Section 3.5.

Proof of Theorem 1
We split the proof of the theorem into two parts, the first handling the case ρ ∈ L 1 (R) ∩ L ∞ (R) and the second, addressing the case We begin by restating the theorem in a form that explicitly includes the distributions that we draw our random variables from.

THEOREM 5 ([14])
For any ε > 0, there exist constants α, Ω > 0 such that the following holds: If, for k ∈ N, the random variables are independently drawn from their associated distributions, and then there exist hidden-to-output layer weights {v k } n k=1 ⊂ R (that depend on the realization of the weights {w k } n k=1 and biases {b k } n k=1 ) such that the sequence of RVFL networks {f n } ∞ n=1 defined by as n → ∞.
PROOF.Our proof technique is based on that introduced by Igelnik and Pao, and can be divided into four steps.The first three steps essentially consist of Lemma 3, Lemma 4, and Lemma 5, and the final step combines them to obtain the desired result.First, the function f is approximated by a convolution, given in Lemma 3. The proof of this result can be found in Section 3.5.1.
Then we have uniformly for all x ∈ R N .
Next, we represent f as the limiting value of a multidimensional integral over the parameter space.In particular, we replace (f * h Ω )(x) in the convolution identity (4) with a function of the form K F (y)ρ(⟨w, x⟩ + b(y))dy, as this will introduce the RVFL structure we require.To achieve this, we first use a truncated cosine function in place of the activation function ρ and then switch back to a general activation function.
To that end, for each fixed Moreover, introduce the functions where y, w ∈ R N and u ∈ R and ϕ = A * A for any even function Then we have the following lemma, a detailed proof of which can be found in Section 3.5.2.
uniformly for every x ∈ K, where The next step in the proof of Theorem 5 is to approximate the integral in (7) using the Monte-Carlo method.Define v k := vol(K(Ω)) n F α,Ω y k , w k α , u k for k = 1, . . ., n, and the random variables {f n } ∞ n=1 by Then, we have the following lemma that is proven in Section 3.5.3. where To complete the proof of Theorem 5 we combine the limit representation (7) with the Monte-Carlo error guarantee (9) and show that, given any ε > 0, there exist α, Ω > 0 such that as n → ∞.To this end, let ε ′ > 0 be arbitrary and consider the integral I(x; p) given by for x ∈ K and p ∈ N. By (7), there exist α, Ω > 0 such that |f (x) − I(x; 1)| < ε ′ holds for every x ∈ K, and so it follows that for every x ∈ K. Jensen's inequality now yields that By ( 9), we know that the second term on the right-hand side of ( 11) is O(1/n).Therefore, we have and so the proof is completed by taking ε ′ = ε/2vol(K) and choosing α, Ω > 0 accordingly.
The full statement of the theorem is identical to that of Theorem 5 albeit now with , so we omit it for brevity.Its proof is also similar to the proof of the case where ρ ∈ L 1 (R) ∩ L ∞ (R) with some key modifications.Namely, one uses an integration by parts argument to modify the part of the proof corresponding to Lemma 4. The details of this argument are presented in the appendix, Section 3.5.4.

Proof of Theorem 2
In this section we prove the non-asymptotic result for RVFL networks in R N , and we begin with a more precise statement of the theorem that makes all the dimensional dependencies explicit.THEOREM 6.Consider the hypotheses of Theorem 5 and suppose further that ρ is κ-Lipschitz on R for some κ > 0. For any , where M := sup x∈K |f (x)|, c > 0 is a numerical constant, and Σ is a constant depending on f and ρ, and let parameters {w k } n k=1 , {b k } n k=1 , and {v k } n k=1 be as in Theorem 5. Then the RVFL network defined by with K := supp(f ) and suppose ε > 0, η ∈ (0, 1) are fixed.Take an arbitrarily κ-Lipschitz activation function ρ ∈ L 1 (R) ∩ L ∞ (R).We wish to show that there exists an RVFL network {f n } ∞ n=1 defined on K that satisfies the with probability at least 1 − η when n is chosen sufficiently large.The proof is obtained by modifying the proof of Theorem 5 for the asymptotic case.
We begin by repeating the first two steps in the proof of Theorem 5 from Sections 3.5.1 and 3.5.2.In particular, by Lemma 4 we have the representation (7), namely, F α,Ω (y, w, u)ρ α⟨w, x⟩ + b α (y, w, u) dydwdu holds uniformly for all x ∈ K. Hence, if we define the random variables f n and I n from Section 3.5.3as in ( 8) and (30), respectively, we seek a uniform bound on the quantity over the compact set K, where I(x; 1) is given by (10) for all x ∈ K. Since equation (7) allows us to fix α, Ω > 0 such that holds for every x ∈ K simultaneously, the result would follow if we show that, with high probability, |I n (x) − I(x; 1)| < ε/2vol(K) uniformly for all x ∈ K. Indeed, this would yield with high probability.To this end, for δ > 0 let C(δ, K) ⊂ K denote a minimal δ-net for K, with cardinality N (δ, K).Now, fix x ∈ K and consider the inequality where z ∈ C(δ, K) is such that ∥x − z∥ 2 < δ.We will obtain the desired bound on ( 12) by bounding each of the terms ( * ), ( * * ), and ( * * * ) separately.
First, we consider the term ( * ).Recalling the definition of I n , observe that we have where M := sup x∈K |f (x)| and we define Now, since ρ is assumed to be κ-Lipschitz, we have Hence, an application of the Cauchy-Schwarz inequality yields R α,Ω (x, z) ≤ καΩδ

√
N for all x ∈ K, from which it follows that holds for all x ∈ K.
Next, we bound ( * * * ) using a similar approach.Indeed, by the definition of I( • ; 1) we have Using the fact that R α,Ω (x, z) ≤ καΩδ holds for all x ∈ K, just like (13).
Combined with the bounds ( 13) and ( 14), it follows from (12) that simultaneously for all x ∈ K with probability at least 1 − η, provided δ and n satisfy ( 15) and ( 16), respectively.Since we require |I n (x) − I(x; 1)| < ε/2vol(K), the proof is then completed by setting ε ′ + t = ε/2vol(K) and choosing δ and n accordingly.In particular, it suffices to choose ε 15) and ( 16) become , as desired.REMARK 2. The implication of Theorem 6 is that, given a desired accuracy level ε > 0, one can construct a RVFL network f n that is ε-close to f with high probability, provided the number of nodes n in the neural network is sufficiently large.In fact, if we assume that the ambient dimension N is fixed here, then δ and n depend on the accuracy ε and probability η as ) for small values of x, the requirement on the number of nodes behaves like n ≳ log η −1 N √ ε, K ε whenever ε is sufficiently small.Using a simple bound on the covering number, this yields a coarse estimate of n ≳ ε −1 log(η −1 /ε).
REMARK 3. If we instead assume that N is variable, then, under the assumption that f is Hölder continuous with exponent β, one should expect that n = ω(N 2βN ) as N → ∞ (in light of Remark 10 and in conjunction with Theorem 6 with log(1 + 1/x) ≈ 1/x for large x).In other words, the number of nodes required in the hidden layer is superexponential in the dimension.This dependence of n on N may be improved by means of more refined proof techniques.As for α, if follows from Remark 12 that α = Θ(1) as REMARK 4. The κ-Lipschitz assumption on the activation function ρ may likely be removed.Indeed, since ( * * * ) in (12) can be bounded instead by leveraging continuity of the L 1 norm with respect to translation, the only term whose bound depends on the Lipschitz property of ρ is ( * ).However, the randomness in I n (that we did not use to obtain the bound (13)) may be enough to control ( * ) in most cases.Indeed, to bound ( * ) we require control over quantities of the form ρ α ⟨w For most practical realizations of ρ, this difference will be small with high probability (on the draws of y k , w k , u k ) whenever ∥x − z∥ 2 is sufficiently small.

Results on submanifolds of Euclidean space
The constructions of RVFL networks presented in Theorems 5 and 6 depend heavily on the dimension of the ambient space R N .Indeed, the random variables used to construct the input-to-hidden layer weights and biases for these neural networks are N -dimensional objects; moreover, it follows from ( 15) and ( 16) that the lower bound on the number n of nodes in the hidden layer depends superexponentially on the ambient dimension N .If the ambient dimension is small, these dependencies do not present much of a problem.However, many modern applications require the ambient dimension to be large.Fortunately, a common assumption in practice is that signals of interest have (e.g., manifold) structure that effectively reduces their complexity.Good theoretical results and algorithms in a number of settings typically depend on this induced smaller dimension rather than the ambient dimension.For this reason, it is desirable to obtain approximation results for RVFL networks that leverage the underlying structure of the signal class of interest, namely, the domain of f ∈ C c (R N ).
One way to introduce lower-dimensional structure in the context of RVFL networks is to assume that supp(f ) lies on a subspace of R N .More generally, and motivated by applications, we may consider the case where supp(f ) is actually a submanifold of R N .To this end, for the remainder of this section, we assume M ⊂ R N to be a smooth, compact d-dimensional manifold and consider the problem of approximating functions f ∈ C(M) using RVFL networks.As we are going to see, RVFL networks in this setting yield theoretical guarantees that replace the dependencies of Theorems 5 and 6 on the ambient dimension N with dependencies on the manifold dimension d.Indeed, one should expect that the random variables {w k } n k=1 , {b k } n k=1 are essentially d-dimensional objects (rather than N -dimensional) and that the lower bound on the number of network nodes in Theorem 6 scales as a (superexponential) function of d rather than N .

Adapting RVFL networks to d-manifolds
As in Section 2.2, let {(U j , ϕ j )} j∈J be an atlas for the smooth, compact d-dimensional manifold M ⊂ R N with the corresponding compactly supported partition of unity {η j } j∈J .Since M is compact, we assume without loss of generality that |J| < ∞.Indeed, if we additionally assume that M satisfies the property that there exists an r > 0 such that, for each x ∈ M, M ∩ B N 2 (x, r) is diffeomorphic to an ℓ 2 ball in R d with diffeomorphism close to the identity.Then one can choose an atlas {(U j , ϕ j )} j∈J with |J| ≲ 2 d T d vol(M)r −d by intersecting M with ℓ 2 balls in R N of radii r/2 [33].Here T d is the so-called thickness of the covering and there exist coverings such that T d ≲ d log(d).
for all x ∈ M, where As we will see, the fact that M is smooth and compact implies fj ∈ C c (R d ) for each j ∈ J, and so we may approximate each fj using RVFL networks on R d as in Theorems 5 and 6.In this way, it is reasonable to expect that f can be approximated on M using a linear combination of these low-dimensional RVFL networks.More precisely, we propose approximating f on M via the following process: 1.For each j ∈ J, approximate fj uniformly on ϕ j (U j ) ⊂ R d using a RVFL network fn j as in Theorems 5 and 6; 2.Approximate f uniformly on M by summing these RVFL networks over J, i.e., for all x ∈ M.

Main results on d-manifolds
We now prove approximation results for the manifold RVFL network architecture described in Section 3.3.1.For notational clarity, from here onward we use lim {n j } j∈J →∞ to denote the limit as each n j tends to infinity simultaneously.The first theorem that we prove is an asymptotic approximation result for continuous functions on manifolds using the RVFL network construction presented in Section 3.3.1.This theorem is the manifold-equivalent of Theorem 5. THEOREM 7. Let M ⊂ R N be a smooth, compact d-dimensional manifold with finite atlas {(U j , ϕ j )} j∈J and f ∈ C(M).Fix any activation function ρ ∈ L 1 (R) ∩ L ∞ (R) with R ρ(z)dz = 1.For any ε > 0, there exist constants α j , Ω j > 0 for each j ∈ J such that the following holds.If, for each j ∈ J and for k ∈ N, the random variables are independently drawn from their associated distributions, and b (j) then there exist hidden-to-output layer weights {v k=1 ⊂ R such that the sequences of RVFL networks { fn j } ∞ n j =1 defined by as {n j } j∈J → ∞.
PROOF.We wish to show that there exist sequences of RVFL networks { fn j } ∞ n j =1 defined on ϕ j (U j ) for each j ∈ J which together satisfy the asymptotic error bound as {n j } j∈J → ∞.We will do so by leveraging the result of Theorem 5 on each ϕ j (U j ) ⊂ R d .
To begin, recall that we may apply the representation (17) for f on each chart (U j , ϕ j ); the RVFL networks fn j we seek are approximations of the functions fj in this expansion.Now, as supp(η j ) ⊂ U j is compact for each j ∈ J, it follows that each set ϕ j (supp(η j )) is a compact subset of R d .Moreover, because fj (z) ̸ = 0 if and only if z ∈ ϕ j (U j ) and ϕ −1 j (z) ∈ supp(η j ) ⊂ U j , we have that fj = fj | ϕ j (supp(η j ) is supported on a compact set.Hence, fj ∈ C c (R d ) for each j ∈ J, and so we may apply Lemma 4 to obtain the uniform limit representation (7) on ϕ j (U j ), that is, where we define In this way, the asymptotic error bound that is the final result of Theorem 5, namely holds.With these results in hand, we may now continue with the main body of the proof.
Since the representation (17) for f on each chart (U j , ϕ j ) yields for all x ∈ M, Jensen's inequality allows us to bound the mean square error of our RVFL approximation by To bound ( * ), note that the change of variables (2) implies for each j ∈ J. Defining β j := inf y∈U j | det(Dϕ j (y))|, which is necessarily bounded away from zero for each j ∈ J by compactness of M, we therefore have Hence, applying (18) for each j ∈ J yields because j∈J 1/n j ≤ |J|/ min j∈J n j .With the bound (20) in hand, ( 19) becomes as {n j } j∈J → ∞, and so the proof is completed by taking each ε j > 0 in such a way that and choosing α j , Ω j > 0 accordingly for each j ∈ J.
REMARK 5. Note that the neural-network architecture obtained in Theorem 7 has the following form in the case of a generic atlas.To obtain the estimate of f (x), the input x is first "pre-processed" by computing ϕ j (x) for each j ∈ J such that x ∈ U j , and then put through the corresponding RVFL network.However, using the Geometric Multi-Resolution Analysis approach from [1] (as we do in Section 3.4), one can construct an approximation (in an appropriate sense) of the atlas, with maps ϕ j being linear.In this way, the pre-processing step can be replaced by the layer computing ϕ j (x), followed by the RVFL layer f j .We refer the reader to Section 3.4 for the details.
The biggest takeaway from Theorem 7 is that the same asymptotic mean-square error behavior we saw in the RVFL network architecture of Theorem 5 holds for our RVFL-like construction on manifolds, with the added benefit that the input-to-hidden layer weights and biases are now d-dimensional random variables rather than N -dimensional.Provided the size of the atlas |J| isn't too large, this significantly reduces the number of random variables that must be generated to produce a uniform approximation of f ∈ C(M).
One might expect to see a similar reduction in dimension dependence for the non-asymptotic case if the RVFL network construction of Section 3.3.1 is used.Indeed, our next theorem, which is the manifold-equivalent of Theorem 6, makes this explicit: THEOREM 8. Let M ⊂ R N be a smooth, compact d-dimensional manifold with finite atlas {(U j , ϕ j )} j∈J and f ∈ C(M).Fix any activation function ρ ∈ L 1 (R) ∩ L ∞ (R) such that ρ is κ-Lipschitz on R for some κ > 0 and R ρ(z)dz = 1.For any ε > 0, there exist constants α j , Ω j > 0 for each j ∈ J such that the following holds.Suppose, for each j ∈ J and for k = 1, ..., n j , the random variables are independently drawn from their associated distributions, and b (j) Then there exist hidden-to-output layer weights {v k=1 ⊂ R such that, for any , and , where M j := sup z∈ϕ j (U j ) | fj (z)|, c > 0 is a numerical constant, and C (j) := 2M j ∥ρ∥ ∞ vol(ϕ j (U j )), the sequences of RVFL networks { fn j } ∞ n j =1 defined by with probability at least 1 − η.
As alluded to earlier, an important implication of Theorems 7 and 8 is that the random variables {w k=1 are d-dimensional objects for each j ∈ J.Moreover, bounds for δ j and n j now have superexponential dependence on the manifold dimension d instead of the ambient dimension N .Thus, introducing the manifold structure removes the dependencies on the ambient dimension, replacing them instead with the intrinsic dimension of M and the complexity of the atlas {(U j , ϕ j )} j∈J .REMARK 6.The bounds on the covering radii δ j and hidden layer nodes n j needed for each chart in Theorem 8 are not optimal.Indeed, these bounds may be further improved if one uses the local structure of the manifold, through quantities such as its curvature and reach.In particular, the appearance of |J| in both bounds may be significantly improved upon if the manifold is locally well-behaved.

Numerical Simulations
In this section, we provide numerical evidence to support the result of Theorem 8. Let M ⊂ R N be a smooth, compact d-dimensional manifold.Since having access to an atlas for M is not necessarily practical, we assume instead that we have a suitable approximation to M. For our purposes, we will use a Geometric Multi-Resolution Analysis (GMRA) approximation of M (see [1]; and also, e.g., [15] for a complete definition).
A GMRA approximation of M provides a collection {(C j , P j )} j∈{1,...J} of centers C j = {c j,k } K j k=1 ⊂ R N and affine projections P j = {P j,k } K j k=1 on R N such that, for each j ∈ {1, . . ., J}, the pairs {(c j,k , P j,k )} K j k=1 define d-dimensional affine spaces that approximate M with increasing accuracy in the following sense.For every x ∈ M, there exists C x > 0 and k ′ ∈ {1, . . ., K j } such that holds whenever ∥x − c j,k ′ ∥ 2 is sufficiently small.In this way, a GMRA approximation of M essentially provides a collection of approximate tangent spaces to M. Hence, a GMRA approximation having fine enough resolution (i.e., large enough j) is a good substitution for an atlas.In practice, one must often first construct a GMRA from empirical data, assumed to be sampled from appropriate distributions on the manifold.Indeed, this is possible, and yields the so-called empirical GMRA, studied in [23], where finite-sample error bounds are provided.The main point is that given enough samples on the manifold, one can construct a good GMRA approximation of the manifold.
Let {(c j,k , P j,k )} K j k=1 be a GMRA approximation of M for refinement level j.Since the affine spaces defined by (c j,k , P j,k ) for each k ∈ {1, . . ., K j } are d-dimensional, we will approximate f on M by projecting it (in an appropriate sense) onto these affine spaces and approximating each projection using an RVFL network on R d .To make this more precise, observe that, since each affine projection acts on x ∈ M as P j,k x = c j,k + Φ j,k (x − c j,k ) for some othogonal projection Φ j,k : R N → R N , for each k ∈ {1, . . .K j } we have where Φ j,k = U j,k D j,k V T j,k is the compact singular value decomposition (SVD) of Φ j,k (i.e., only the left and right singular vectors corresponding to nonzero singular values are computed).In particular, the matrix of right-singular vectors V j,k : R d → R N enables us to define a function fj,k : R d → R, given by which satisfies fj,k (V T j,k x) = f (P j,k x) for all x ∈ M. By continuity of f and ( 21), this means that for any ε > 0 there exists j ∈ N such that |f (x) − fj,k (V T j,k x)| < ε for some k ∈ {1, . . ., K j }.For such k ∈ {1, . . ., K j }, we may therefore approximate f on the affine space associated with (c j,k , P j,k ) by approximating fj,k using a RFVL network fn j,k : R d → R of the form where {w ℓ=1 ⊂ R are random input-to-hidden layer weights and biases (resp.)and the hidden-to-output layer weights {v ℓ=1 ⊂ R are learned.Choosing the activation function ρ and random input-to-hidden layer weights and biases as in Theorem 8 then guarantees that |f (P j,k x) − fn j,k (V T j,k x)| is small with high probability whenever n j,k is sufficiently large.In light of the above discussion, we propose the following RVFL network construction for approximating functions f ∈ C(M): Given a GMRA approximation of M with sufficiently high resolution j, construct and train RVFL networks of the form (23) for each k ∈ {1, . . ., K j }.Then, given x ∈ M and ε > 0, choose k ′ ∈ {1, . . ., K j } such that and evaluate fn j,k ′ (x) to approximate f (x).We summarize this algorithm in Algorithm 1.Since the structure of the GMRA approximation implies ∥x − P j,k ′ x∥ 2 ≤ C x 2 −2j holds for our choice of k ′ ∈ {1, . . ., K j } [see 15], continuity of f and Lemma 5 imply that, for any ε > 0 and j large enough, holds with high probability, provided n j,k ′ satisfies the requirements of Theorem 8.

Algorithm 1 Approximation Algorithm
Step 1: For each k ∈ {1, . . ., K j }, construct and train1 a RVFL network fn j,k of the form (23) Step 2: For any x ∈ M, find c j,k ′ ∈ arg min c j,k ∈C j ∥x − c j,k ∥ 2 Step 3: Set y ♯ = fn j,k ′ (x) REMARK 7. In the RVFL network construction proposed above we require that the function f be defined in a sufficiently large region around the manifold.Essentially, we need to ensure that f is continuously defined on the set S := M ∪ M j , where M j is the scale-j GMRA approximation This ensures that f can be evaluated on the affine subspaces given by the GMRA.
To simulate Algorithm 1, we take M = S 2 embedded in R 20 and construct a GMRA up to level j max = 15 using 20,000 data points sampled uniformly from M. Given j ≤ j max , we generate RVFL networks fn j,k : R 2 → R as in (23) and train them on V T j,k (B N 2 (c j,k , r) ∩ T j,k ) using the training pairs {(V T k,j x ℓ , f (P j,k x ℓ ))} p ℓ=1 , where T k,j is the affine space generated by (c j,k , P j,k ).For simplicity, we fix n j,k = n to be constant for all k ∈ {1, . . ., K j } and use a single, fixed pair of parameters α, Ω > 0 when constructing all RVFL networks.We then randomly select a test set of 200 points x ∈ M for use throughout all experiments.In each experiment (i.e., point in Figure 1), we use Algorithm 1 to produce an approximation y ♯ = fn j,k ′ (x) of f (x). Figure 1 displays the mean relative error in these approximations for varying numbers of nodes n; to construct this plot, f is taken to be the exponential f (x) = exp( N k=1 x(k)) and ρ the hyperbolic secant function.Notice that for small numbers of nodes the RVFL networks are not very good at approximating f , regardless of the choice of α, Ω > 0. However, the error decays as the number of nodes increases until reaching a floor due to error inherent in the GMRA approximation.Hence, as suggested by Theorem 3, to achieve a desired error bound of ε > 0, one needs to only choose a GMRA scale j such that the inherent error in the GMRA (which scales like 2 −j ) is less than ε, then adjust the parameters α j , Ω j , and n j,k accordingly.REMARK 8.As we just mentioned, the error can only decay so far due to the resolution of the GMRA approximation.However, that is not the only floor in our simulation; indeed, the ε in Theorem 3 is determined by the α j 's and Ω j 's, which we keep fixed (see the caption of Figure 1).Consequently, the stagnating accuracy as n increases, as seen in Figure 1, is also predicted by Theorem 3. Since the solid and dashed lines seem to reach the same floor, the floor due to error inherent in the GMRA approximation seems to be the limiting error term for RVFL networks with large numbers of nodes.REMARK 9. Utilizing random inner weights and biases resulted in us needing to approximate the atlas to the manifold.To this end, knowing the computational complexity of the GMRA approximation would be useful in practice.As it turns out in [22], calculating the GMRA approximation has computational complexity O(C d N m log(m)), where m is the number of training data points and C > 0 is a numerical constant.

Proofs of technical lemmas 3.5.1 Proof of Lemma 3
Observe that h Ω defined in (3) may be viewed as a multidimensional bump function; indeed, the parameter Ω > 0 controls the width of the bump.In particular, if Ω is allowed to grow very large, then h Ω becomes very localized near the origin.Objects that behave in this way are known in the functional analysis literature as approximate δ-functions: For such functions, we write δ 0 (x) = lim t→∞ φ t (x) for all x ∈ R N , where δ 0 denotes the N -dimensional Dirac δ-function centered at the origin.
Given φ ∈ L 1 (R N ) with R N φ(x)dx = 1, one may construct approximate δ-functions for t > 0 by defining φ t (x) := t N φ(tx) for all x ∈ R N [34].Such sequences of approximate δ-functions are also called approximate identity sequences [31] since they satisfy a particularly nice identity with respect to convolution, namely, lim t→∞ ∥f * φ t − f ∥ 1 = 0 for all f ∈ C c (R N ) [see 31, Theorem 6.32].In fact, such an identity holds much more generally.LEMMA 6. [34, Theorem 1.18] Let φ ∈ L 1 (R N ) with R N φ(x)dx = 1 and for t > 0 define φ t (x) := ).For each j, we fix α j = 2 and vary Ω j = 10, 15 (solid and dashed lines, resp.).Reconstruction error decays as a function of n until reaching a floor due to error in the GMRA approximation of M. The code used to obtain these numerical results is available upon direct request sent to the corresponding author.
To prove (4), it would suffice to have lim Ω→∞ ∥f * h Ω − f ∥ ∞ = 0, which is really just Lemma 6 in case p = ∞.Nonetheless, we present a proof by mimicking [34] for completeness.Moreover, we will use a part of proof in Remark 10 below.
Then, for all f ∈ C 0 (R N ), we have PROOF.By symmetry of the convolution operator in its arguments, we have Since a simple substitution yields Finally, expanding the function h Ω , we obtain where we have used the substitution z = Ωy.Taking limits on both sides of this expression and observing that using the Dominated Convergence Theorem, we obtain So, it suffices to show that, for all z ∈ R N , To this end, let ε > 0 and z ∈ R N be arbitrary.Since f ∈ C 0 (R N ), there exists r > 0 sufficiently large such that |f (x)| < ε/2 for all x ∈ R N \ B(0, r), where B(0, r) ⊂ R N is the closed ball of radius r centered at the origin.Let B := B(0, r + ∥z/Ω∥ 2 ), so that for each x ∈ R N \ B we have both x and x − z/Ω in R N \ B(0, r).Thus, both |f (x)| < ε/2 and |f (x − z/Ω)| < ε/2, implying that Hence, we obtain Now, as B is a compact subset of R N , the continuous function f is uniformly continuous on B, and so the remaining limit and supremum may be freely interchanged, whereby continuity of f yields Since ε > 0 may be taken arbitrarily small, we have proved the result.
REMARK 10.While Lemma 7 does the approximation we aim for, it gives no indication of how fast decays in terms of Ω or the dimension N. Assuming h(z) = g(z(1)) • • • g(z(N )) for some nonnegative g (which is how we will choose h in Section 3.5.2) and f to be β-Hölder continuous for some fixed β ∈ (0, 1) yields that where the third inequality follows from Jensen's inequality.

Proof of Lemma 5: Monte-Carlo integral approximation
The next step in the proof of Theorem 5 is to approximate the integral in (7) using the Monte-Carlo method.To this end, let {y k } n k=1 , {w k } n k=1 , and {u k } n k=1 be independent samples drawn uniformly from K, [−Ω, Ω] N , and [− π 2 (2L + 1), π 2 (2L + 1)], respectively, and consider the sequence of random variables {I n (x)} ∞ n=1 defined by for each x ∈ K, where we note that vol(K(Ω)) = (2Ω) N π(2L + 1)vol(K).If we also define for x ∈ K and p ∈ N, then we want to show that as n → ∞, where the expectation is taken with respect to the joint distribution of the random samples {y k } n k=1 , {w k } n k=1 , and {u k } n k=1 .For this, it suffices to find a constant C(f, ρ, α, Ω, N ) < ∞ independent of n satisfying Indeed, an application of Fubini's theorem would then yield which implies (32).To determine such a constant, we first observe by Theorem 4 that where we define the variance term for all y, w ∈ R N and u ∈ R, where M := sup x∈K |f (x)| < ∞, we obtain the following simple bound on the variance term Since we assume ρ ∈ L ∞ (R), we then have Substituting the value of vol(K(Ω)), we obtain is a suitable choice for the desired constant.
because A is even.

Proof of Theorem 1 when ρ
) with K := supp(f ) and suppose ε > 0 is fixed.Take the activation function ρ : R → R to be differentiable with ρ ′ ∈ L 1 (R) ∩ L ∞ (R).We wish to show that there exists a sequence of RVFL networks {f n } ∞ n=1 defined on K which satisfy the asymptotic error bound as n → ∞.The proof of this result is a minor modification of second step in the proof of Theorem 5.
If we redefine h α (z) as αρ ′ (αz), then (27) plainly still holds and ( 29) reads Recalling the definition of cos Ω in (5) and integrating by parts, we obtain cos Ω (u)dρ(α(z − u)) = − cos Ω (u)ρ(α(z − u)) for y, w ∈ R N and u ∈ R, we again obtain the uniform representation (7) for all x ∈ K.The remainder of the proof proceeds from this point exactly as in the proof of Theorem 5.

DISCUSSION
The central topic of this paper is the study of the approximation properties of a randomized variation of shallow neural networks known as RVFL.In contrast with the classical single-layer neural networks, training of an RVFL involves only learning the output weights, while the input weights and biases of all the nodes are selected at random from an appropriate distribution and stay fixed throughout the training.The main motivation for studying the properties of such networks is as follows: 1.Random weights are often utilized as an initialization for a NN training procedure.Thus, establishing the properties of the RVFL networks is an important first step toward understanding how random weights are transformed during training.2.Due to their much more computationally efficient training process, the RVFL networks proved to be a valuable alternative to the classical SLFNs.They were successfully used in several modern applications, especially those that require frequent re-training of a neural network [6,38,7].
Despite their practical and theoretical importance, results providing rigorous mathematical analysis of the properties of RVFLs are rare.The work of Igelnik and Pao [14] showed that RVFL networks are universal approximators for the class of continuous, compactly supported functions and established the asymptotic convergence rate of the expected approximation error as a function of the number of nodes in the hidden layer.While this result served as a theoretical justification for using RVFL networks in practice, a close examination led us to the conclusion that the proofs in [14] contained several technical errors.
In this paper, we offer a revision and a modification of the proof methods from [14] that allow us to prove a corrected, slightly weaker version of the result announced by Igelnik and Pao.We further build upon their work and show a non-asymptotic probabilistic (instead of on average) approximation result, which gives an explicit bound on the number of hidden layer nodes that are required to achieve the desired approximation accuracy with the desired level of certainty (that is, with high enough probability).In addition to that, we extend the obtained result to the case when the function is supported on a compact, low-dimensional submanifold of the ambient space.
While our work closes some of the gaps in the study of the approximation properties of RVFL, we believe that it just starts the discussion and opens many directions for further research.We briefly outline some of them here.
In our results, the dependence of the required number n of the nodes in the hidden layer on the dimension N of the domain is superexponential, which is likely an artifact of the proof methods we use.We believe this dependence can be improved to be exponential by using a different, more refined approach to the construction of the limit-integral representation of a function.A related interesting direction for future research is to study how the bound on n changes for more restricted classes of (e.g., smooth) functions.
Another important direction that we did not discuss in this paper is learning the output weights and studying the robustness of the RVFL approximation to the noise in the training data.Obtaining provable robustness guarantees for an RVFL training procedure would be a step towards the robustness analysis of neural networks.

Figure 1 .
Figure1.Log-scale plot of average relative error for Algorithm 1 as a function of the number of nodes n in each RVFL network.Black (cross), blue (circle), and red (square) lines correspond to GMRA refinement levels j = 12, j = 9, and j = 6 (resp.).For each j, we fix α j = 2 and vary Ω j = 10, 15 (solid and dashed lines, resp.).Reconstruction error decays as a function of n until reaching a floor due to error in the GMRA approximation of M. The code used to obtain these numerical results is available upon direct request sent to the corresponding author.