# Integrable Abel equation and asymptotics of
symmetry solutions of Korteweg-de Vries equation

B.I. Suleimanov, A.M. Shavlukov
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We provide a general solution for a first order ordinary differential equation
with a rational right-hand side, which arises in constructing asymptotics for
large time of simultaneous solutions of the Korteweg-de Vries equation and the
stationary part of its higher non-autonomous symmetry. This symmetry is
determined by a linear combination of the first higher autonomous symmetry of
the Korteweg-de Vries equation and of its classical Galileo symmetry. This
general solution depends on an arbitrary parameter. By the implicit function
theorem, locally it is determined by the first integral explicitly written in
terms of hypergeometric functions. A particular case of the general solution
defines self-similar solutions of the Whitham equations, found earlier by G.V.
Potemin in 1988. In the well-known works by A.V. Gurevich and L.P. Pitaevsky
in early 1970s, it was established that these solutions of the Whitham
equations describe the origination in the leading term of non-damping
oscillating waves in a wide range of problems with a small dispersion. The
result of this article supports once again an empirical rule saying that under
various passages to the limits, integrable equations can produce only
integrable, in certain sense, equations. We propose a general conjecture:
integrable ordinary differential equations similar to that considered in the
present paper should also arise in describing the asymptotics at large times
for other symmetry solutions to evolution equations admitting the application
of the method of inverse scattering problem.

##  1 Limiting passages in integrable equations

According to a well-known empirical law, when describing the asymptotics of
solutions to integrable differential equations, as a result of reasonable
passage to the limit, there can arise only solutions of integrable equations,
but often in another sense different from the sense of integrability of
original equations.

For instance, all possible continuous limits of such integrable differential-
difference equations as Volterra chain

$\displaystyle(c_{n})^{\prime}_{t}=c_{n}(c_{n+1}-c_{n-1})$

and Toda chain

$\displaystyle(z_{n})^{\prime\prime}_{tt}=e^{z_{n+1}-z_{n}}-e^{z_{n+1}-z_{n}}$

can be only differential equations integrable in some sense. The same concerns
continuous limits of integrable chain discrete in two independent variables.
At present, some of this limiting passages are rigorously justified [1].

###### Remark 1.1

The scientific achievements of R.I. Yamilova made known first of all owing to
his papers devoted to various aspects of the integrability of chains of such
kind. In his articles [2]-[5], in particular, specific examples of transitions
from integrable chains to their continual limits were given.

One of the simplest examples of the limiting passage from one integrable
equation to another is the transition from the Korteweg-de Vries (KdV)
equation with a small dispersion integrable by the inverse scattering problem
method (ISPM) ($\varepsilon\to 0$)

$\displaystyle
u^{\prime}_{t}+uu^{\prime}_{x}+\varepsilon^{2}u^{\prime\prime\prime}_{xxx}=0,$

to the Hopf equation

$\displaystyle u^{\prime}_{t}+uu^{\prime}_{x}=0,$

the general solution of which is locally given by the formula

$\displaystyle x-tu=F(u).$

The same Hopf equation is a non-dissipative limit of the Burgers equation as
$\varepsilon\to 0$:

$\displaystyle
u^{\prime}_{t}+uu^{\prime}_{x}=\varepsilon^{2}u^{\prime\prime}_{xx}.$

This equation is reduced to the heat equation
$\Lambda^{\prime}_{t}=\varepsilon^{2}\Lambda^{\prime\prime}_{xx}$ by a
linearizing Cole-Hopf change
$u=-2\varepsilon^{2}\Lambda^{\prime}_{x}/\Lambda$.

A quasiclassical approximation of solutions to an integrable by ISPM nonlinear
Schrödinger equation (NSE) with a small dispersion

$\displaystyle-i\varepsilon
q^{\prime}_{t}=\varepsilon^{2}q^{\prime\prime}_{xx}+2\delta|q|^{2}q=0\qquad(\varepsilon\ll
1)$ (1.1)

provides a more complicated example of such passage: the substitution

$\displaystyle q=\sqrt{h}\exp{\left(i\frac{\varphi}{\varepsilon}\right)}$

reduces (1.1) to a system of two real evolution equations

$\displaystyle
h^{\prime}_{t}+2(h\varphi^{\prime}_{x})^{\prime}_{x}=0,\qquad\varphi^{\prime}_{t}+(\varphi^{\prime}_{x})^{2}-2\delta
h=\varepsilon^{2}\frac{(\sqrt{h})^{\prime\prime}_{xx}}{\sqrt{h}},$

and their non-dispersion limit is the system

$h^{\prime}_{t}+2(h\varphi^{\prime}_{x})^{\prime}_{x}=0,\qquad\varphi^{\prime}_{t}+(\varphi^{\prime}_{x})^{2}-2\delta
h=0.$ (1.2)

After differentiating the second equation in system (1.2) with respect to the
variable $x$ and the change $v=2\varphi^{\prime}_{x}$, this non-dispersion
limit becomes a classical hydrodynamical system

$\displaystyle h^{\prime}_{t}+(hv)^{\prime}_{x}=0,\quad
v^{\prime}_{t}+vv^{\prime}_{x}-4\delta h^{\prime}_{x}=0;$

the sign of the constant $\delta$ defines its hyperbolic and elliptic version.
Under the hodograph transformation

$\displaystyle v(t,x),h(t,x)\to t(h,v),x(h,v),$

the solutions of this equation are expressed via the solutions of the linear
system of equations

$\displaystyle x^{\prime}_{h}=vt^{\prime}_{h}+4\delta t^{\prime}_{v},\qquad
x^{\prime}_{v}=vt^{\prime}_{v}-ht^{\prime}_{h}.$

And a very curious transformation of the meaning of integrability occurs
[6]-[9] when passing to non-dispersion limits of spatially multidimensional
integrable ISM equations of the type, for example, the Kadomtsev-Petviashvili
equation

$\displaystyle\frac{\partial}{\partial
t}\left(u^{\prime}_{t}+uu^{\prime}_{x}+\varepsilon^{2}u^{\prime\prime\prime}_{xxx}\right)=u^{\prime\prime}_{yy}.$

The integrability of Whitham’s hydrodynamic equations by the generalized
hodograph method of S.P. Tsarev [11], [12] [10], obtained as a result of
averaging integrable ISM equations of the KdV, NSE and sine-Gordon type, is
also one of the confirmations of the empirical law formulated in the beginning
of the paper.

For example, using this method in combination with the algebraic-geometric
method of I.M. Krichever [13], G.V. Potemin [14], self-similar solutions of
the Whitham equations were found in explicit form, which, according to the
well-known results of A.V. Gurevich and L.P. Pitaevskii [15], [16] determine,
in the leading order, the behavior as $t\to\infty$ of a universal special
solution of the KdV equation

$\displaystyle u^{\prime}_{t}+uu^{\prime}_{x}+u^{\prime\prime\prime}_{xxx}=0,$
(1.3)

with the asymptotics $u=-x^{\frac{1}{3}}+o(1)$ as $x\pm\infty$.

###### Remark 1.2

In [16], A.V. Gurevich and L.P. Pitaevsky expressed an opinion that the system
of ODEs determining these self-similar solutions of Whitham equations cannot
be solved explicitly.

Ordinary differential equations (ODEs) of Painlevé type integrable by the
isomonodromic deformations method [17], often arising in describing the
asymptotics at large times for integrable by ISPM equations, are also examples
of such limiting equations.

###### Remark 1.3

Among the numerous important scientific results of A.B. Shabat, there is also
the following one: when describing the asymptotics at large times of a general
initial value problem for the KdV equation (1.3) in [18] and in Section 8.3 of
his Habilitation thesis [19], he gave the first example of such transition.

Along with the Painlevé equations themselves, in nonlinear problems with a
small parameter including those described by solutions of integrable
equations, their higher isomonodromic analogs play a universal role in the
description of various contrast transition regimes. Both the Painlevé equation
and their higher analogs can be regarded as a kind of nonlinear special
functions of wave catastrophes.

The general theory of such nonlinear special functions, the foundations of
which were laid in the work of A.V. Kitaev [20] (simultaneously and
independently one particular case, a higher analogue of the second Painlevé
equation, was considered by the first of the authors of this work in [21]) is
actively developed today and it finds numerous applications [22]–[50].

In particular, when describing solutions of non-dispersion equations in the
vicinity of the folding point, in the leading term with respect to a small
parameter $\varepsilon$ one uses a special universal solution of the Gurevich-
Pitaevsky KdV equation (1.3), which was discussed in the paragraph above
Remark 1.2. And as it was stated in [24], [26], this special solution at the
same time satisfies ODE

$\displaystyle
u^{\prime\prime\prime\prime}_{xxxx}+\frac{5uu^{\prime\prime}_{xx}}{3}+\frac{5(u^{\prime}_{x})^{2}}{6}+\frac{5(u^{3}-tu+x)}{18}=0,$
(1.4)

being the first higher representative of in the hierarchy of isomonodromic
analogues of the first Painlevé equation in the mentioned paper by A.V. Kitaev
[20].

###### Remark 1.4

The same solution of ODE (1.4) appears in describing of continuous limits of
isomonodromic chains [26], [37]. They are related with problems of quantum
gravitation theory and were considered in [51], [52].

In [28], V.R. Kudashev observed the following fact: after differentiating
(1.4) with respect to $x$, one obtains a stationary part of a non-autonomous
higher symmetry of the KdV equation (1.3), which is determined by the linear
combination of the stationary parts of its classical Galilean symmetry
$u^{\prime}_{\tau_{G}}=1-tu^{\prime}_{x}$ and by its first higher autonomous
symmetry

$\displaystyle
u^{\prime}_{\tau_{5}}=\left(u^{\prime\prime\prime\prime}_{xxxx}+\frac{5u^{\prime\prime}_{xx}u}{3}+\frac{5(u^{\prime}_{x})^{2}}{6}+\frac{5u^{3}}{18}\right)^{\prime}_{x}.$
(1.5)

The symmetry nature of ODE (1.4) and the results of paper [53] allow one to
obtain easily the form of aforementioned solution by G.V. Potemin of self-
similar Whitham equation [24], [26], [28]. In the last quarter of a century,
various properties of solutions to ODE (1.4), and mainly, of course, the
special solution of Gurevich-Pitaevsky, have been considered from various
points of view in many works, see, for example, [30], [32]–[35], [37],
[39]–[42], [54].

##  2 Kudashev integrable equation

The asymptotics of the universal Gurevich-Pitaevsky solution as $t\to\infty$
outside the zone of fast oscillation is given by two roots of the cubic fold
equation:

$\displaystyle u^{3}-tu+x=0,$ (2.1)

while in the area of non-decaying oscillations it depends on a slow variable
$z=xt^{-\frac{3}{2}}$ and fast phase

$\displaystyle\Phi=t^{\frac{7}{4}}f(z)+f_{0}(z)$ (2.2)

and it reads as

$\displaystyle
u=\sqrt{t}(v_{0}(z,\Phi)+t^{-\frac{7}{4}}v_{1}(z,\Phi)+t^{-\frac{7}{2}}v_{2}(z,\Phi)+\dots).$
(2.3)

This is almost obviously implied by the form

$\displaystyle u=\sqrt{t}v(z)$ (2.4)

of the solution to fold equation (2.1) and the form of linearizations by (2.4)
of KDV equation (1.3) and ODE (1.4).

In the end of the last century, V.R. Kudashev tried, without employing Whitham
averaging method, simply to look for a joint asymptotic solution of form (2.3)
to equations (1.3) and (1.4). While doing this, he discovered that the
function $f(z)$ determining the leading term of the fast phase (2.2) by the
formula

$\displaystyle R(z)=\frac{7f(z)}{4f^{\prime}_{z}}-\frac{3z}{2},$

is related with a solution of a rather simple first order ODE

$\displaystyle
R^{\prime}_{z}=\frac{486R^{4}-171R^{2}+9zR+5}{9(54R^{3}-9R+z)(2R+3z)}.$ (2.5)

This equation is obviously equivalent to the Abel equation of second kind:

$\displaystyle(486R^{4}-171R^{2}+5+9Rz)z^{\prime}_{R}=972R^{4}-162R^{2}+(1458R^{3}-225R)z+27z^{2}.$
(2.6)

###### Remark 2.1

In terms of solutions of equation (2.5), it is easy to describe the solution
of a more general ODE

$\displaystyle
r^{\prime}_{x}=\frac{486r^{4}-171tr^{2}+9xr+5t^{2}}{9(54r^{3}-9tr+x)(2tr+3x)}.$
(2.7)

Indeed, as $t=0$, ODE (2.7) is reduced to ODE $r^{\prime}_{x}=r/(3x)$ with the
general solution $r=const\cdot x^{\frac{1}{3}}$. As $t\neq 0$, the solutions
of (2.5) and (2.7) are related by the dilatations

$\displaystyle r(t,x)=t^{\frac{1}{2}}R(z),\qquad z=xt^{-\frac{3}{2}}.$

In [37], there were written formulae according to which in the case of a
special Gurevich - Pitaevsky solution, the solution $R(z)$ of the ODE (2.5) is
related to the self-similar solutions of the Whitham equations from the
article by G.V. Potemin [14]. V.R. Kudashev knew these formulae but never
published them during his lifetime.

However, after obtaining V.R. Kudashev equations (2.6) a natural question
immediately arose: according to the heuristic law outlined in Section 1 of
this article, this equation must undoubtedly be integrable. But how? Below we
provide an answer to this question found last year by using the computer
algebra system Maple.

Kudashev equation (2.6) turns out to possess the first integral, which is
written explicitly in terms of hypergeometric functions
${}_{2}\textrm{{F}}_{1}(\alpha,\beta;\gamma;y)\equiv
F(\alpha,\beta;\gamma;y)$. This integral reads as ($I$ is an arbitrary
constant)

$\displaystyle
I=\frac{997920i\sqrt{15}(F(-\frac{7}{6},-\frac{1}{6};\frac{2}{3};c_{1})(a_{1}b_{6}+b_{8})+\frac{12789F(-\frac{1}{6},\frac{5}{6};\frac{5}{3};c_{1})a_{1}c_{2}}{220})}{c_{1}^{\frac{1}{3}}(F(-\frac{5}{6},\frac{1}{6};\frac{4}{3};c_{1})b_{7}+77711400F(\frac{1}{6},\frac{7}{6};\frac{7}{3};c_{1})a_{1}(b_{5}+b_{9}))},$

where the functions $a_{i}$, $b_{i}$, $c_{i}$ are given by the formulae

$\displaystyle
a_{1}(R)=R^{2}-\frac{1}{3},a_{2}(R)=R^{4}-\frac{25}{111}R^{2}+\frac{1}{1332},a_{3}(R,z)=R^{3}-\frac{2R}{7}+\frac{z}{14},$
$\displaystyle
a_{4}(R)=R^{6}-\frac{55}{174}R^{4}+\frac{19}{1044}R^{2}-\frac{1}{6264},a_{5}(R)=R^{6}-\frac{503}{3198}R^{4}+\frac{25}{6396}R^{2}-\frac{1}{115128},$
$\displaystyle
a_{6}(R)=R^{2}(3R^{2}-1),a_{7}(R)=-\frac{2665}{1421}R^{6}+\frac{2515}{8526}R^{4}-\frac{125}{17052}R^{2}+\frac{5}{306936},$
$\displaystyle
a_{8}(R)=\sqrt{3R^{4}-R^{2}},a_{9}(R)=R^{2}-\frac{1}{24},a_{10}(R)=660R^{4}-76R^{2}+1,$
$\displaystyle a_{11}(R)=1980R^{6}-888R^{4}+79R^{2}-1,$ $\displaystyle
b_{1}(R,z)=2i\sqrt{15}a_{8}(R)(3R^{2}-1)+126R^{4}+9zR-36R^{2},$ $\displaystyle
b_{2}(R)=\sqrt{21}(31968R^{6}-8532R^{4}+324R^{2}-1)+a_{8}(R)(-83160R^{4}+9576R^{2}-126),$
$\displaystyle
b_{3}(R)=-\frac{148\sqrt{21}a_{2}(R)a_{8}(R)a_{9}(R)}{1155}+R^{8}-\frac{74}{165}R^{6}+\frac{79R^{4}-R^{2}}{1980},$
$\displaystyle
b_{4}(R)=-\frac{533i\sqrt{35}a_{1}(R)a_{5}(R)}{1421}+5a_{4}(R)R^{2},b_{5}(R)=a_{1}(R)R^{2}\left(3i\sqrt{35}a_{5}(R)+\frac{4263a_{4}(R)}{533}\right),$
$\displaystyle
b_{6}(R)=\frac{1}{220}\left(i\sqrt{35}a_{6}(R)^{\frac{3}{2}}a_{10}(R)+19188\sqrt{21}a_{3}(R)a_{5}(R)R-7308i\sqrt{15}a_{1}(R)a_{4}(R)R^{2}\right),$
$\displaystyle
b_{7}(R,z)=15120i\sqrt{15}a_{6}(R)^{\frac{3}{2}}a_{11}(R)-124338240R(i\sqrt{15}a_{1}(R)R(\sqrt{21}a_{1}(R)a_{5}(R)$
$\displaystyle\hphantom{b_{7}(R,z)=}-\frac{1036a_{2}(R)a_{8}(R)a_{9}(R)}{533})-\frac{29841a_{3}(R,z)a_{4}(R)R^{2}}{533}+7\sqrt{21}a_{3}(R)a_{5}(R)a_{8}(R)),$
$\displaystyle
b_{8}(R,z)=\frac{444a_{8}(R)}{385}\left(3i\sqrt{35}a_{1}(R)^{2}a_{2}(R)a_{9}(R)-\frac{29841a_{3}(R,z)a_{4}(R)R}{148}\right),$
$\displaystyle
b_{9}(R)=a_{8}(R)\left(-\sqrt{21}a_{1}(R)a_{5}(R)-\frac{1421i\sqrt{15}a_{4}(R)R^{2}}{533}\right),$
$\displaystyle
c_{1}(R,z)=-\frac{2993760i\sqrt{15}a_{1}(R)b_{3}(R)}{b_{1}(R,z)b_{2}(R)},$
$\displaystyle
c_{2}(R)=a_{8}(R)b_{4}(R)+a_{1}(R)R^{2}\left(\sqrt{21}a_{7}(R)+i\sqrt{15}a_{4}(R)\right).$

In the following table we provide approximate values of $I$ corresponding to
solutions of five different initial problems for ODE (2.6).

Initial value | Approximate value of $I$ | Interval of calculation in $R$
---|---|---
$z(0)=0$ | $0.0960605+0.1663816i$ | $R\in[-10,10]$
$z(1)=0$ | $0.0194046+0.0336097i$ | $R\in[-11,11]$
$z(5)=7$ | $0.0308202+0.0533821i$ | $R\in[-1,15]$
$z(50)=75$ | $0.03177193+0.05503056i$ | $R\in[-1,60]$
$z(-0.57735)=0.3849$ | $0.000591+0.0010237i$ | $R\in[-10,0]$
Table 1: Approximate values of $I$ corresponding to five different initial
problems.

##  3 Conclusion

In the opinion of the authors of the article, its result is potentially not
only of a particular importance.

For example, there arises a conjecture that the systems of two non-autonomous
ODEs

$\displaystyle\begin{aligned}
&A^{\prime}_{s}=\frac{(2A-1)(-288A^{3}+192A^{2}+24sA-27A-4s+4B)}{(A+2s)(-576A^{3}+504A^{2}-126A+48sA+8B-12s+9)},\\\
&B^{\prime}_{s}=\left(36A+3s-54A^{2}-\frac{9}{2}\right)A^{\prime}_{s}+\frac{108A^{3}-108A^{2}-6sA+6B+27A+4s}{2A+4s},\end{aligned}$
(3.1) $\displaystyle\begin{aligned}
&A^{\prime}_{s}=\frac{(2A-1)(288A^{3}-192A^{2}+24sA+27A-4s-4B)}{(A-2s)(-576A^{3}+504A^{2}-126A-48sA+8B+12s+9)},\\\
&B^{\prime}_{s}=(36A-3s-54A^{2}-\frac{9}{2})A^{\prime}_{s}-\frac{108A^{3}-108A^{2}+6sA+6B+27A-4s}{2A-4s},\end{aligned}$
(3.2)

considered in papers by R.N. Garifullin [44], [46] should also possess two
first integrals, which can be explicitly written in terms of hypergeometric
functions. By means of solutions to systems of ODEs (3.1) and (3.2), one
describes the asymptotics as $t\to\infty$ of form

$\displaystyle
u=t(v_{0}(z,\Phi)+t^{-\frac{5}{4}}v_{1}(z,\Phi)+t^{-\frac{5}{2}}v_{2}(z,\Phi)+\dots),\quad\Phi=t^{\frac{5}{2}}f(s)+n(s)\quad\left(s=\frac{x}{t^{2}}\right)$

for joint solutions of KdV equation (1.3) and fifth order ODE

$\displaystyle\left(u^{\prime\prime\prime\prime}_{xxxx}+\frac{5u^{\prime\prime}_{xx}u}{3}+\frac{5(u^{\prime}_{x})^{2}}{6}+\frac{5u^{3}}{18}\right)^{\prime}_{x}\pm\frac{2u+xu^{\prime}_{x}-3t(u^{\prime\prime\prime}_{xxx}+uu^{\prime}_{x})}{6}=0,$

determined by linear combinations of stationary parts of the first higher non-
autonomous symmetry (1.5) and a classical dilatation symmetry
$u_{\tau_{r}}=2u+xu^{\prime}_{x}-3t(u^{\prime\prime\prime}_{xxx}+uu^{\prime}_{x})$.

We propose a conjecture: a similar asymptotics at large times for joint
solutions of various nonlinear integrable by ISPM evolution equations with
ODEs of the higher Painlevé type determined by non-autonomous symmetries and,
more generally, by invariant manifolds of these evolution equations, should
also be described by ODEs integrable in a similar sense.

## References

  * [1] L.A. Kalyakin. Asymptotic transitions from discrete to continuous models // Teor. Matem. Fiz. 76:3, 323–327 (1988). [Theor. Math. Phys. 76:3, 891–894 (1988).]
  * [2] R.I. Yamilov. On classification of discrete evolution equations // Uspekih Matem. Nauk. 38:6(234), 155–156 (1983). (in Russian).
  * [3] D. Levi, P. Winternitz, R.I. Yamilov. Symmetries of the continuous and discrete Krichever-Novikov equation // SIGMA. 7, 097 (2011).
  * [4] R.N. Garifullin, R.I. Yamilov. On integrability of a discrete analogue of Kaup-Kupershmidt equation // Ufimskji Matem. Zhurn. 9:3, 158–164 (2017). [Ufa Math. J. 9:3, 158–164 (2017).]
  * [5] R.N. Garifullin, R.I. Yamilov. On the integrability of a lattice equation with two continuum limits // Itogi Nauki Tekh. Ser. Sovrem. Matem. Pril. Temat. Obz. VINITI RAN. 152, 159–164 (2018). (in Russian).
  * [6] V.E. Zakharov. Dispersionless limit of integrable systems in 2 + 1 dimensions // In: ‘‘Singular limits of dispersive waves’’, Ercolani, N.M. et al. (eds.). Plenum Press, New York, 165–174 (1994).
  * [7] E.V. Ferapontov, K.R. Khusnutdinova. Hydrodynamic reductions of multi- dimensional dispersionless PDEs: the test for integrability // J. Math. Phys. 45:6, 2365–2377 (2004).
  * [8] E.V. Ferapontov, K.R. Khusnutdinova. On integrability of (2+1)-dimensional quasilinear systems // Comm. Math. Phys. 48, 187–206 (2004).
  * [9] E.V. Ferapontov, B. Kruglikov. Dispersionless integrable systems in 3D and Einstein-Weyl geometry // J. Diff. Geom. 97, 215–254 (2014).
  * [10] B.A. Dubrovin, S.P. Novikov. Hydrodinamics of weakly deformed soliton lattices. Differential geometry and Hamiltonian theoy // Uspekhi Matem. Nauk. 44:6(270), 29–98 (1989). [ Russ. Math. Surv. 44:6, 35–124 (1989).]
  * [11] S.P. Tsarev. On Poisson brackets and one-dimensional Hamiltonian systems of hydrodynamic type // Dokl. AN SSSR. 282:3, 534–537 (1985). [Sov. Math. Dokl. 31, 488–491 (1985).]
  * [12] S.P. Tsarev. The geometry of hamiltonian systems of hydrodynamic type. The generalized hodograph method // Izv. AN SSSR. Ser. Matem. 54:5, 1048–1068 (1990). [Math. USSR-Izv. 37:2, 397–419 (1991).]
  * [13] I.M. Krichever. Method of averaging for two-dimensional ‘‘integrable’’ equations  // Funkts. Anal. Pril. 22:3, 37–52 (1988). [Funct. Anal. Appl. 22:3, 200–213 (1988).]
  * [14] G.V. Potemin. Algebro-geometric construction of self-similar solutions of the Whitham equations // Uspekhi Matem. Nauk. 43:5(263), 211–212 (1988). [Russ. Math. Surv. 43:5, 252–253 (1988).]
  * [15] A.V. Gurevich, L. P. Pitaevskii. Breaking a simple wave in the kinetics of a Rarefied Plasma // Zhurn. Eksper. Teor. Fiz. 60:6, 2155–2174 (1971). [Sov. Phys. JETP. 33:6, 1159–1167 (1971).]
  * [16] A.V. Gurevich, L.P. Pitaevskii. Non stationare structure of collisionless shock wave // Zhurn. Eksper. Teor. Fiz. 65:2, 590–604 (1973). [Sov. Phys. JETP. 38:2, 291–297 (1974).]
  * [17] A Fokas, A.R. Its, A.A. Kapaev, V.Yu. Novokshenov. Painlevé transcendents. The Riemann-Hilbert approach. Regul. Khaot. Dinam., Moscow-Izhevsk (2005). [Amer. Math. Soc., Providence, RI (2006).]
  * [18] A.B. Shabat. On the Korteweg-de Vries equation // Dokl. AN SSSR. 211:6, 1310–1313 (1973). [Sov. Math. Dokl. 14, 1266–1270 (1973).]
  * [19] A.B. Shabat. Transformation operators and nonlinear equations // Habiliation thesis. Bashkir State Univ., Ufa (1974). (in Russian).
  * [20] A.V. Kitaev. Turning points of linear systems and double asymptotics of the Painlevé transcendents // Zap. LOMI 187, 53–74 (1991). [J. Math. Sci. 73:4, 44–459, (1995).]
  * [21] B.I. Suleimanov. The second Painlevé equation in a problem concerning nonlinear effects near caustics // Zap. LOMI. 187, 110–128 (1991). [J. Math. Sci. 73:4, 482–493 (1995).]
  * [22] B.I. Suleimanov. A ‘‘nonlinear’’ generalization of special functions of wave catastrophes described by double integrals  // Matem. Zamet. 52:5, 102–106 (1992). [Math. Notes. 52:5,6, 1146–1149 (1992).]
  * [23] B.I. Suleimanov, I.T. Habibullin. Symmetries of the Kadomtsev-Petviashvili equation, the isomonodromic deformations, and ‘‘nonlinear’’ generalizations of special functions of wave catastrophes // Teor. Matem. Fiz. 97:2, 213–216 (1993). [Theor. Math. Phys. 97:2, 1250–1258, (1993).]
  * [24] B.I. Suleimanov.Solution of the Korteweg—de Vries equation which arises near the breaking point in problems with a slight dispersion // Pisma ZhETF. 58:11, 606–610 (1993). [JETP Lett. 58:11, 849–854 (1993).]
  * [25] B.I. Suleimanov. Influence of a weak nonlinearity on the high-frequency asymptotics in caustic rearrangements // Teor. Matem. Fiz. 98:2, 198–206 (1994). [Theor. Math. Phys. 98:2, 132–138 (1994).]
  * [26] B.I. Suleimanov. Onset of nondissipative shock waves and the ‘‘nonperturbative’’ quantum theory of gravitattion // Zhurn. Eskper. Teor. Fiz. 105:5, 1089–1099 (1994). [J. Exper.Theor. Phys. 78:5, 583–587 (1994).]
  * [27] A.V. Kitaev. Caustics in $1+1$ integrable systems // J. Math. Phys. 35:2, 2934–2954 (1994).
  * [28] V.R. Kudashev. KdV shock-like waves as invariant solutions of KdV equaton symmetry // Preprint: arXiv preprint patt-sol/9404002. (1994).
  * [29] V.R. Kudashev, B.I. Suleimanov. Characteristic features of some typical spontaneous intensivity collapse processes in unstable media // Pisma ZhETF. 62:4, 358–362 (1995). [JETP Lett. 62:4, 382–388, (1995).]
  * [30] V. Kudashev, B. Suleimanov. A soft mechanism for generation the dissipationless shock waves // Phys. Lett. A. 221:3, 4, 204–208 (1996). 1996.
  * [31] V.R. Kudashev, B.I. Suleimanov. Small-amplitude dispersion oscillations on the background of the nonlinear geometric optics approximation // Teor. Matem. Fiz. 118:3, 413–422 (1999). [Theor. Math. Phys. 118:3, 325–332 (1999).]
  * [32] V.R. Kudashev, B.I. Suleimanov. The effect of small dissipation on the onset of one-dimensional shock waves // Prikl. Matem. Mekh. 65:3, 456–466 (2001). [J. Appl. Math. Mech. 65:3, 441–451 (2001).]
  * [33] B. Dubrovin. On Hamiltonian perturbations of hyperbolic systems of conservation laws, II: Universality of critical behaviour // Comm. Math. Phys. 267, 117–139 (2006).
  * [34] T. Claeys, M. Vanlessen. The existence of a real pole-free solution of the fourth order analogue of the Painlevé I equation // Nonlinearity. 20:5, 1163–1184 (2007).
  * [35] T. Grava, C. Klein. Numerical study of a multiscale expansion of the Korteweg-de Vries equation and Painlevé-II equation // Proc. Royal Society A. 464:2091, 733–757 (2008).
  * [36] B. Dubrovin, T. Grava, C. Klein. On Universality of Critical Behavior in the Focusing Nonlinear Schrödinger Equation, Elliptic Umbilic Catastrophe and the Tritronquée Solution to the Painlevé-I Equation // J. Nonl. Sci. 19:1, 57–94 (2009).
  * [37] R. Garifullin, B. Suleimanov, N. Tarkhanov. Phase Shift in the Whitham Zone for the Gurevich-Pitaevskii Special Solution of the Korteweg–de Vries Equation // Phys. Lett. A. 374:13, 14, 1420–1424 (2010). 2010.
  * [38] R.N. Garifullin, B.I. Suleimanov. From weak discontinues to nondissipative shock waves  // Zh. Eksp. Teor. Fiz. 137:1, 149–165 (2010). [J. Exper. Theor. Phys. 110:1, 133–146 (2010).]
  * [39] T. Claeys. Asymptotics for a special solutions to the second member of the Painleve I hierarhy // J. Phys. A. 43:43, 434012. 18 pp. (2010).
  * [40] T. Claeys, T. Grava. Solitonic asymptotics for the Korteweg-de Vries equation in the small dispersion limit // SIAM J. Math. Anal. 42:5, 2132–2154 (2010).
  * [41] T. Claeys. Pole-free solutions of the first Painlevé Hierarchy and non-generic critical behavior for the KdV equation // Physica D: Nonl. Phen. 241:23, 2226–2236 (2011).
  * [42] T. Grava, A. Kapaev, C. Klein. On the Tritronquée Solutions of $P_{I}^{2}$ // Constr. Approx. 41:3, 425–466 (2015).
  * [43] B. Dubrovin, M. Elaeva. On the critical behavior in nonlinear evolutionary PDEs with small viscosity // Russ. J. Math. Phys. 19:4, 449–460 (2012).
  * [44] R.N. Garifullin. Phase shift for the common solution of KdV and fifth order differential equation // Ufimskij Matem. Zhurn. 4:2, 80–86 (2012). [Ufa Math. J. 4:2, 52-61 (2012).]
  * [45] M. Bertola, A. Tovbis. Universality for the Focusing Nonlinear Schrödinger Equation at the Gradient Catastrophe Point: Rational Breathers and Poles of the Tritronquée Solution to Painlevé I // Comm. Pure and Appl. Math. 66, 678–752 (2013).
  * [46] R.N. Garifullin. On simultaneous solution of the KdV equation and a fifth-order differential equation // Ufimskij Matem. Zhurn. 8:4, 53–62 (2016). [Ufa Math. J. 8:4, 80-86 (2016).]
  * [47] B. I. Suleimanov. Effect of a small dispersion on self-focusing in a spatially one-dimensional case // Pisma ZhETF. 106:6, 375–380 (2017). [JETP Lett. 6:6, 400–405 (2017).]
  * [48] B.I. Suleimanov. On analogs of wave catastrophe functions that are solutions of nonlinear integrable equations // Itogi Nauki Tekh. Ser. Sovrem. Matem. Pril. Temat. Obz. VINITI RAN. 163, 81–95 (2019). (in Russian).
  * [49] D. Bilman, L. Ling, P.D. Miller. Extreme superposition: rogue waves of infinite order and the Painlevé-III hierarchy // Duke Math. J. 169, 671–760 (2020).
  * [50] V.E. Adler. Nonautonomous symmetries of the KdV equation and step-like solutions// J. of Nonlinear Math. Phys. 27:3, 478–493 (2020).
  * [51] E. Bresin, E. Marinari, G. Parisi. A nonperturbative ambiguty free solution of a string model // Phys. Lett. B. 242:1, 35-38 (1990).
  * [52] M. Douglas, N. Seiberg, S. Shenker. Flow and unstability in quantum gravity // Phys. Lett. B. 244:3,4, 381–386 (1990).
  * [53] V.R. Kudashev, S.E. Sharapov. Inheritance of KdV symmetries under Whitham averaging and hydrodynamic symmetries of the Witham equations // Teor. Matem. Fiz. 87:1, 40–47 (1991). [Theor. Math. Phys. 87:1, 358–363 (1991).]
  * [54] B.I. Suleimanov. ‘‘Quantizations’’ of higher Hamiltonian analogues of the Painlevé I and Painlevé II equations with two degrees of freedom // Funkts. Anal. Pril. 48:3, 52–62 (2014). [Funct. Anal.Appl. 48:3, 198–207 (2014).]


# Expressivity and Approximation Properties of Deep
Neural Networks with ${\rm ReLU}^{k}$ Activation

Juncai He111Computer, Electrical and Mathematical Science and Engineering
Division, King Abdullah University of Science and Technology, Thuwal 23955,
Saudi Arabia. Tong Mao111Computer, Electrical and Mathematical Science and
Engineering Division, King Abdullah University of Science and Technology,
Thuwal 23955, Saudi Arabia. Jinchao Xu111Computer, Electrical and Mathematical
Science and Engineering Division, King Abdullah University of Science and
Technology, Thuwal 23955, Saudi Arabia. 222Department of Mathematics, The
Pennsylvania State University, University Park, PA 16802, USA.

###### Abstract

In this paper, we investigate the expressivity and approximation properties of
deep neural networks employing the ReLUk activation function for $k\geq 2$.
Although deep ReLU networks can approximate polynomials effectively, deep
ReLUk networks have the capability to represent higher-degree polynomials
precisely. Our initial contribution is a comprehensive constructive proof for
polynomial representation using deep ReLUk networks. This allows us to
establish an upper bound on both the size and count of network parameters.
Consequently, we are able to demonstrate a suboptimal approximation rate for
functions from Sobolev spaces as well as for analytic functions. Additionally,
through an exploration of the representation power of deep ReLUk networks for
shallow networks, we reveal that deep ReLUk networks can approximate functions
from a range of variation spaces, extending beyond those generated solely by
the ReLUk activation function. This finding demonstrates the adaptability of
deep ReLUk networks in approximating functions within various variation
spaces.

##  1 Introduction

Over the past decade, neural networks have achieved remarkable success in a
variety of fields, including image recognition [30, 29, 23, 21], speech
recognition [40, 24, 45], and natural language processing [8, 3, 26, 53],
among others. These achievements have sparked increased interest in the
theoretical understanding of neural networks. Theoretically, the
generalization error of a neural network can be decomposed into the sum of the
approximation error and the sampling error, as detailed in works such as [9].
In practice, the reduction of approximation error depends on carefully
designed neural network architectures, whereas the control of sample error
requires sufficient and reliable data resources. For instance, the
approximation capabilities of neural networks with sigmoid activation
functions were explored in various works [10, 25, 2, 39]. Additionally, the
universality of neural networks with non-polynomial activation functions was
established in [31, 5, 38]. In recent years, due to their simple computational
form and the potential to overcome the gradient vanishing problem, neural
networks employing the Rectified Linear Unit (ReLU, defined as
$\text{ReLU}(x)=\max\\{0,x\\}$) and ReLU${}^{k}=(\text{ReLU})^{k}$ activation
functions have gained widespread use. This paper focuses on the expressivity
and approximation properties of neural networks with ReLUk activation
functions.

Given the fact that any ReLU networks are continuous piecewise linear (CPwL)
functions, [1] initially demonstrated that any CPwL functions on $[0,1]^{d}$
can be represented by deep ReLU networks with $\lceil\log_{2}(d+1)\rceil$
hidden layers. Subsequently, [20] established that shallow ReLU networks
cannot represent any piecewise linear function on $[0,1]^{d}$ with $d\geq 2$.
In contrast, deep ReLU networks with more layers (but fewer parameters) can
recover any linear finite element functions. Recently, [17] explored the
optimal representation expressivity for continuous piecewise linear functions
using deep ReLU networks on $[0,1]$. By integrating this expressivity with the
Kolmogorov Superposition Theorem [28], a significantly improved approximation
rate was achieved. Most recently, through a detailed examination of the basis
functions and high-dimensional simplicial mesh characteristics of Lagrange
finite element functions, which are a special type of piecewise polynomial
functions, the authors in [22] demonstrated that Lagrange finite element
functions of any order in arbitrary dimensions can be represented by DNNs with
ReLU or ReLU2 activation functions.

On the other hand, the approximation theory of shallow and deep neural
networks with ReLU and ReLUk activation functions has been studied from
various aspects. For functions from the Barron space, it has been proved in
[27] that shallow ReLU networks uniformly with the rate
$\mathcal{O}\left(N^{-\frac{1}{2}-\frac{1}{d}}\right)$ with $\mathcal{O}(N)$
parameters. This result has been generalized to $L^{p}$-approximation,
spectral Barron spaces, and ReLUk activation functions in [49, 50, 51, 34]. In
particular, a optimal rate
$\tilde{\mathcal{O}}(N^{-\frac{1}{2}-\frac{2k+1}{2d}})$ for functions from the
so-called spectral Barron spaces, where $\tilde{\mathcal{O}}$ indicates up to
a constant power of $\log N$. Furthermore, for $\alpha$-Hölder functions, [37]
proved the approximation rate
$\tilde{\mathcal{O}}\left(N^{-\frac{\alpha}{d}\frac{d+2}{d+4}}\right)$ by
shallow ReLU networks. This result is generalized to shallow ReLUk networks
for $k=0,1,2,\dots$ in [56]. In particular, for $k=0,1$, the rate is improved
to the optimal rate $\mathcal{O}\left(N^{-\frac{\alpha}{d}}\right)$. Compared
to shallow neural networks, deep network architectures have significant
advantages as demonstrated in [6, 11, 44, 14]. Specifically, by using a bit-
extraction technique, it has been proved that fully-connected deep ReLU
networks with depth $L$ and width $N$ attain the optimal approximation rate
$\tilde{\mathcal{O}}\left(N^{2}L^{2}\right)$ [46, 47, 48].

For functions with specific compositional structures, deep ReLU networks can
achieve significantly better expressivity power or approximation accuracy than
their shallow ReLU counterparts [43, 52, 36]. In particular, it has been shown
that any function constructed by a shallow ReLU neural network can be
represented by a deep ReLU neural network with width $d+4$ [33] or a deep
convolutional neural network whenever the filter size $\geq 2$ [58, 18]. In
addition, it was demonstrated in [54] that sawtooth functions on $[0,1]$ with
$2^{N}$ oscillations can be represented by a deep ReLU network with
$\mathcal{O}(N)$ parameters. In contrast, a shallow ReLU network would require
$\mathcal{O}(2^{N})$ parameters to construct such a function. Subsequent
studies in [33, 57] showed that multivariate polynomials could be approximated
by deep ReLU networks exponentially, achieving an approximation rate of
$\tilde{\mathcal{O}}(N^{-\frac{\alpha}{d}})$ for functions from
$C^{\alpha}(\mathbb{B}^{d})$. Later, a further study of approximating
holomorphic functions by deep ReLU networks has been introduced in [42].
Additionally, [41] demonstrated that the sparse grid basis can be
exponentially approximated by deep ReLU networks, leading to a suboptimal
approximation rate of $\tilde{\mathcal{O}}(N^{-2})$ for functions from Korobov
spaces. However, a foundational aspect of these results is a special
expressivity property of deep ReLU networks. A more in-depth and systematic
analysis of this result, from the hierarchical basis perspective, is presented
in [19].

Given these findings, it is reasonable to anticipate that deep ReLUk networks
might also offer superior representation and approximation properties compared
to their shallow counterparts and deep ReLU networks. However, the
expressivity and approximation properties of neural networks exclusively using
the ReLUk activation function are not extensively documented in the
literature. For example, it has been noted in studies such as [4] and [55]
that ReLUk networks can represent any global polynomial, with existence proofs
provided to support this claim. Furthermore, the authors in [32] constructed a
deep ReLU2 network with $\mathcal{O}(n^{d})$ neurons, demonstrating its
capability to represent all polynomials of degree $\leq n$ on
$\mathbb{R}^{d}$. In this paper, we generalize this to all ReLUk networks,
providing a constructive proof for polynomial representation and estimating
the upper bound of parameter size and count. Specifically, we show that
shallow ReLUk networks with $\mathcal{O}(k^{d})$ parameters can represent
polynomials of degree at most $k$ in $\mathbb{R}^{d}$, and that deep ReLUk
networks of depth $L$ with $\mathcal{O}(k^{Ld})$ parameters can represent
polynomials of degree $\leq k^{L}$. Moreover, deep ReLUk networks can
approximate analytic functions and functions from Sobolev spaces with
suboptimal rates. Additionally, we prove that deep ReLUk networks with depth
$L$ and width $2(k+1)n$ can represent functions constructed by any shallow
ReLU${}^{k^{\ell}}$ networks of width $n$ for any $\ell\leq L$. Together with
recent results on variation spaces [51], we demonstrate an approximation error
of $\mathcal{O}(N^{-\frac{1}{2}-\frac{2k^{\ell}+1}{2d}})$ for functions from
the variation space $\mathcal{K}_{1}(\mathbb{P}_{k^{\ell}}^{d})$ for any
$\ell\leq L$. This underscores the adaptability of deep ReLUk networks in
approximating functions from variation spaces without precise knowledge of
their regularity, providing new insights into the advantages of deep
architectures with ReLUk activations.

The rest of this section will formally define our ReLUk network architectures,
including bounds on parameter size and count. Section 2 presents the
constructive proof for global polynomial representation by ReLUk networks and
discusses their approximation properties for analytic functions and functions
from Sobolev spaces. Section 3 explores the adaptive approximation properties
of deep ReLUk networks for functions in variation spaces. Section 4 offers
concluding remarks and reflections on our findings.

######  Definition 1.

Let $k\in\mathbb{N}$, then we denote the univariate function ReLUk as

$\sigma_{k}(x)=\left\\{\begin{array}[]{ll}x^{k},&\quad x\geq 0,\\\ 0,&\quad
x<0.\end{array}\right.$ (1.1)

In particular, the ReLU function is the ReLUk function with $k=1$.

For any vector $v=(v_{1},\dots,v_{n})$, we abuse the notation and denote the
vector mapping

$\sigma_{k}(v)=(\sigma_{k}(v_{1}),\dots,\sigma_{k}(v_{n}))\in\mathbb{R}^{n}.$
(1.2)

######  Definition 2 (ReLUk networks).

Let $d,k,n\in\mathbb{N}^{+}$, $B>0$, $\mathbb{B}^{d}$ be the unit ball of
$\mathbb{R}^{d}$ centered at the origin. Denote $\Sigma_{n}^{k}(B)$ be the
class of shallow ReLUk network on $\mathbb{B}^{d}$

$\Sigma_{n}^{k}(B)=\left\\{x\mapsto c\cdot\sigma_{k}(Ax+b):\
A\in\mathbb{R}^{n\times d},\ b\in\mathbb{R}^{n},\ c\in\mathbb{R}^{n},\
\|A\|_{\max}\leq B,\ \|b\|_{\infty}\leq B\right\\}.$ (1.3)

where

$\|A\|_{\max}:=\max_{\begin{subarray}{c}1\leq i\leq n\\\ 1\leq j\leq
d\end{subarray}}|(A)_{ij}|.$

In addition, if we suppose further in (1.3) that $\|c\|_{\infty}<M$ for some
$M>0$, we will denote the set as $\Sigma_{n,M}^{k}(B)$. This coincides with
the notation in the previous works [49, 50, 51].

The class of fully-connected deep ReLUk networks with depth $L$ and width
$\\{n_{i}\\}_{i=1}^{L}$ is defined inductively as

$\begin{split}\Sigma_{n_{1:L}}^{k}(B)=\\{x\mapsto c\cdot h_{L}(x):&\
h_{i+1}(x)=\sigma_{k}\left(A_{i}h_{i}(x)+b_{i}\right),\
A_{i}\in\mathbb{R}^{n_{i}\times n_{i-1}},\ b_{i}\in\mathbb{R}^{n_{i}},\
c\in\mathbb{R}^{n_{L}},\\\ &\ \|A_{i}\|_{\max}\leq B,\ \|b_{i}\|_{\infty}\leq
B,i=1,\dots,L\\},\end{split}$ (1.4)

where $n_{0}=d$, $h_{0}(x)=x$.

Similarly, if we further require the weights in the output layers to be
bounded as $\|c\|_{\infty}\leq M$, then we denote the class as
$\Sigma_{n_{1:L},M}^{k}(B)$.

Before introducing our results, we specify that the number of parameters in
the fully connected deep ReLUk network is
$n_{L}+\sum\limits_{i=1}^{L}n_{i}(n_{i-1}+1)$. However, deep neural network
structures are often considered to be sparse (e.g.,[7, 41, 35, 15]), i.e.,
most of the parameters are fixed to be $0$. Such a network is called a sparse
deep network. The sparsity remarkably decreases the number of parameters and,
thus, the computation cost. The deep ReLUk networks we constructed in this
paper are sparse neural networks.

##  2 Representing Polynomials

In this section, we construct the ReLUk neural networks that represent
polynomials. We will first show how shallow ReLUk networks represent
polynomials of degree bounded by $k$, and extend the result to deep networks.
As a consequence, we obtain suboptimal approximation rates for analytic
functions and functions from Sobolev spaces.

The following lemma shows how monomials with the form $x^{\alpha}$ of degree
$n$ are constructed by $n$-th powers of linear combinations of
$(x_{1},\dots,x_{d})$ (see, e.g., [5, 38]). It has been known for many years
that the homogeneous polynomial space of degree $n$ is a span of such linear
combinations. However, our following lemma construct it explicitly, which
allows an estimation of the coefficients.

######  Lemma 1.

Let $n,d\in\mathbb{N}^{+}$, and
$\alpha=(\alpha_{1},\dots,\alpha_{d})\in\mathbb{N}^{d}$ with

$\alpha_{1}+\dots+\alpha_{d}=n,$

then the monomial $x^{\alpha}$ can be written as a linear combination

$x^{\alpha}=\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{n}{2}\right\rfloor}^{n-\left\lfloor\frac{n}{2}\right\rfloor}c_{n_{2},\dots,n_{d}}\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}\right)^{n}$
(2.1)

with each
$c_{n_{2},\dots,n_{d}}\in\left[-\left(\frac{n}{2}+1\right)^{2d},\left(\frac{n}{2}+1\right)^{2d}\right]$.

###### Proof.

We prove this inductively on $d$. Clearly, this is true for $d=1$. We assume
the result is true for $d-1$, thus

$x_{1}^{\alpha_{1}}\dots
x_{d-1}^{\alpha_{d-1}}=\sum\limits_{n_{2},\dots,n_{d-1}=-\left\lfloor\frac{j}{2}\right\rfloor}^{j-\left\lfloor\frac{j}{2}\right\rfloor}c_{n_{2},\dots,n_{d-1}}\left(x_{1}+n_{2}x_{2}+\dots+n_{d-1}x_{d-1}\right)^{j},$
(2.2)

where $j=n-\alpha_{d}$ and
$c_{n_{2},\dots,n_{d-1}}\in\left[-\left(\frac{n}{2}+1\right)^{2(d-1)},\left(\frac{n}{2}+1\right)^{2(d-1)}\right]$.

Let $\xi(x)$ be any function of $x$ and consider $(\xi(x)+n_{d}x_{d})^{n}$ for
$n_{d}=-\left\lfloor\frac{n}{2}\right\rfloor,\dots,n-\left\lfloor\frac{n}{2}\right\rfloor$,
then we can write

$(\xi(x)+n_{d}x_{d})^{n}=\sum\limits_{i=0}^{n}\binom{n}{i}n_{d}^{n-i}\xi(x)^{i}x_{d}^{n-i}.$

Consequently, for any $b=(b_{0},\dots,b_{n})^{T}\in\mathbb{R}^{n+1}$, we have

$\sum\limits_{s=0}^{n}b_{s}\left[\xi(x)+\left(s-\left\lfloor\frac{n}{2}\right\rfloor\right)x_{d}\right]^{n}=b^{T}B_{n}\lambda$
(2.3)

where
$\lambda=\left(\xi(x)^{n},x_{d}\xi(x)^{n-1},\dots,x_{d}^{n-1}\xi(x),x_{d}^{n}\right)^{\top}$,
and

$\begin{split}B_{n}=&\left[\begin{array}[]{cccc}\binom{n}{0}1&\binom{n}{1}\left(-\left\lfloor\frac{n}{2}\right\rfloor\right)&\dots&\binom{n}{n}\left(-\left\lfloor\frac{n}{2}\right\rfloor\right)^{n}\\\
\binom{n}{0}1&\binom{n}{1}\left(1-\left\lfloor\frac{n}{2}\right\rfloor\right)&\dots&\binom{n}{n}\left(1-\left\lfloor\frac{n}{2}\right\rfloor\right)^{n}\\\
\vdots&\vdots&\ddots&\vdots\\\
\binom{n}{0}1&\binom{n}{1}\left(n-\left\lfloor\frac{n}{2}\right\rfloor\right)&\dots&\binom{n}{n}\left(n-\left\lfloor\frac{n}{2}\right\rfloor\right)^{n}\end{array}\right]\\\
=&\left[\begin{array}[]{cccc}1&\left(-\left\lfloor\frac{n}{2}\right\rfloor\right)&\dots&\left(-\left\lfloor\frac{n}{2}\right\rfloor\right)^{n}\\\
1&\left(1-\left\lfloor\frac{n}{2}\right\rfloor\right)&\dots&\left(1-\left\lfloor\frac{n}{2}\right\rfloor\right)^{n}\\\
\vdots&\vdots&\ddots&\vdots\\\
1&\left(n-\left\lfloor\frac{n}{2}\right\rfloor\right)&\dots&\left(n-\left\lfloor\frac{n}{2}\right\rfloor\right)^{n}\end{array}\right]\left[\begin{array}[]{cccc}\binom{n}{0}&0&\dots&0\\\
0&\binom{n}{1}&\dots&0\\\ \vdots&\vdots&\ddots&\vdots\\\
0&0&\dots&\binom{n}{n}\end{array}\right].\end{split}$ (2.4)

Clearly, $B_{n}$ is invertible. Let $\hat{b}=e_{n-j+1}^{T}B_{n}^{-1}$ with

$e_{n-j+1}=\left(\underbrace{0,\dots,0}_{n-j},1,\underbrace{0,\dots,0}_{j}\right)^{T}$

and $\xi(x)=x_{1}+n_{2}x_{2}+\dots+n_{d-1}x_{d-1}$, then we have

$\begin{split}&\sum\limits_{s=0}^{n}\hat{b}_{s}^{T}\left[x_{1}+n_{2}x_{2}+\dots+n_{d-1}x_{d-1}+\left(s-\left\lfloor\frac{n}{2}\right\rfloor\right)x_{d}\right]^{n}=\hat{b}^{T}B_{n}\lambda=e_{n-j+1}\lambda=\xi(x)^{j}x_{d}^{n-j}\\\
=&(x_{1}+n_{2}x_{2}+\dots+n_{d-1}x_{d-1})^{j}x_{d}^{\alpha_{d}}.\end{split}$

Substituting (2.2), we get

$\begin{split}&\sum\limits_{n_{2},\dots,n_{d-1}=-\left\lfloor\frac{j}{2}\right\rfloor}^{j-\left\lfloor\frac{j}{2}\right\rfloor}c_{n_{2},\dots,n_{d-1}}\sum\limits_{s=0}^{n}\hat{b}_{s}\left[x_{1}+n_{2}x_{2}+\dots+n_{d-1}x_{d-1}+\left(s-\left\lfloor\frac{n}{2}\right\rfloor\right)x_{d}\right]^{n}\\\
=&\sum\limits_{n_{2},\dots,n_{d-1}=-\left\lfloor\frac{j}{2}\right\rfloor}^{j-\left\lfloor\frac{j}{2}\right\rfloor}c_{n_{2},\dots,n_{d-1}}(x_{1}+n_{2}x_{2}+\dots+n_{d-1}x_{d-1})^{j}x_{d}^{\alpha_{d}}\\\
=&x_{1}^{\alpha_{1}}\dots
x_{d-1}^{\alpha_{d-1}}x_{d}^{\alpha_{d}}.\end{split}$

For
$n_{2},\dots,n_{d}\in\left\\{-\left\lfloor\frac{n}{2}\right\rfloor,\dots,n-\left\lfloor\frac{n}{2}\right\rfloor\right\\}$,
denote

$c_{n_{2},\dots,n_{d}}:=\left\\{\begin{array}[]{ll}c_{n_{2},\dots,n_{d-1}}\hat{b}_{n_{d}+\left\lfloor\frac{n}{2}\right\rfloor},&\quad
n_{2},\dots,n_{d-1}\in\left\\{-\left\lfloor\frac{j}{2}\right\rfloor,\dots,j-\left\lfloor\frac{j}{2}\right\rfloor\right\\}\\\
0,&\quad\hbox{otherwise,}\end{array}\right.$

then

$\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{n}{2}\right\rfloor}^{n-\left\lfloor\frac{n}{2}\right\rfloor}c_{n_{2},\dots,n_{d}}\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}\right)^{n}=x_{1}^{\alpha_{1}}\dots
x_{d}^{\alpha_{d}}.$

Finally, let us consider the bounds of $c_{n_{2},\dots,n_{d}}$. Here, we only
need to estimate the bounds of $\hat{b}_{s}$. Since
$\hat{b}=e_{n-j+1}^{T}B_{n}^{-1}$, for $s=0,\dots,n$, we have
$|b_{s}|\leq\|B_{n}^{-1}\|_{\max}=\max\limits_{1\leq i,k\leq
n+1}|(B_{n}^{-1})_{ik}|$.

We can bound $\|B_{n}^{-1}\|_{\max}$ by the norm of the inverse of the
Vandermonde matrix in (2.4). Then by [16, Theorem 1], we have

$\begin{split}\|B_{n}^{-1}\|_{\max}\leq&\frac{\prod\limits_{s=0}^{n}\left(1+\left|s-\left\lfloor\frac{n}{2}\right\rfloor\right|\right)}{\min\limits_{0\leq
i\leq
n}\left\\{(1+\left|i-\left\lfloor\frac{n}{2}\right\rfloor\right|)\prod\limits_{s\neq
i}|s-i|\right\\}}=\frac{\prod\limits_{s=0}^{n}\left(1+\left|s-\left\lfloor\frac{n}{2}\right\rfloor\right|\right)}{(1+\left|0\right|)\prod\limits_{s\neq\left\lfloor\frac{n}{2}\right\rfloor}\left|s-\left\lfloor\frac{n}{2}\right\rfloor\right|}\\\
=&\left(\prod\limits_{m=1}^{n-\left\lfloor\frac{n}{2}\right\rfloor}\frac{1+m}{m}\right)\times
1\times\left(\prod\limits_{m=1}^{\left\lfloor\frac{n}{2}\right\rfloor}\frac{1+m}{m}\right)\leq\left(\frac{n}{2}+1\right)^{2}.\end{split}$
(2.5)

Thus

$|c_{n_{2},\dots,n_{d}}|\leq|b_{s}||c_{n_{2},\dots,n_{d-1}}|\leq\left(\frac{n}{2}+1\right)^{2}\left(\frac{n}{2}+1\right)^{2(d-1)}\leq\left(\frac{n}{2}+1\right)^{2d}.$

This completes the proof of Lemma 1. ∎

By noticing that

$t^{k}=\sigma_{k}(t)+(-1)^{k}\sigma_{k}(-t),\quad t\in\mathbb{R},$

Lemma 1 enables the representation of the polynomials of degree $\leq k$ by
shallow ReLUk neural networks, which is also the key to showing the
expressively of deep ReLUk neural networks.

######  Remark 1.

In fact, the choice of $n_{2},\dots,n_{d}$ in the proof of Lemma (1) is not
the only choice. From (2.4), we can observe that our proof relies on writing
$B_{n}$ as a product Vandermonde matrix and a diagonal matrix. The attempt to
lessen the bound (2.5) resulted in our choice of $n_{2},\dots,n_{d}$. However,
we specify that our choice does not necessarily lead to the optimal upper
bound $M$ in the next lemma.

######  Lemma 2.

Let $d,k\in\mathbb{N}^{+}$, $k\geq 2$, $N=2(k+1)^{d}$, and $B>0$, then any
polynomial

$P(x)=\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\ \|\alpha\|_{1}\leq
k\\}}a_{\alpha}x^{\alpha},\ x\in\mathbb{B}^{d}$

with degree $\leq k$ can be represented by the shallow ReLUk network
$\Sigma_{N,M}^{k}(B)$, where

$M=B^{-k}\left(\frac{k}{2}+1\right)^{2(d+1)+k}\left(\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k\\}}|a_{\alpha}|\right).$

###### Proof.

By applying Lemma 1 with $d\leftarrow d+1$ and $n\leftarrow k$, we conclude
any monomial of degree $k$ can be represented as

$x_{1}^{\alpha_{1}}\dots
x_{d}^{\alpha_{d}}x_{d+1}^{\alpha_{d+1}}=\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{k}{2}\right\rfloor}^{k-\left\lfloor\frac{k}{2}\right\rfloor}c_{n_{2},\dots,n_{d+1}}\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1}x_{d+1}\right)^{k}.$

Replacing $x_{d+1}$ by $1$, then all monomials of degree $\leq k$ on variables
$(x_{1},\dots,x_{d})$ can be represented as

$x_{1}^{\alpha_{1}}\dots
x_{d}^{\alpha_{d}}=\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{k}{2}\right\rfloor}^{k-\left\lfloor\frac{k}{2}\right\rfloor}c_{n_{2},\dots,n_{d}}\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1}\right)^{k}$
(2.6)

with
$c_{n_{2},\dots,n_{d+1}}=c_{n_{2},\dots,n_{d+1}}(\alpha)\in\left[-\left(\frac{k}{2}+1\right)^{2(d+1)},\left(\frac{k}{2}+1\right)^{2(d+1)}\right]$.
Then any polynomial

$P(x)=\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\ \|\alpha\|_{1}\leq
k\\}}a_{\alpha}x^{\alpha}$

can be written as

$\begin{split}P(x)=&\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq
k\\}}a_{\alpha}\sum\limits_{n_{2},\dots,n_{d+1}=-\left\lfloor\frac{k}{2}\right\rfloor}^{k-\left\lfloor\frac{k}{2}\right\rfloor}c_{n_{2},\dots,n_{d+1}}(\alpha)\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1}\right)^{k}\\\
=&\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{k}{2}\right\rfloor}^{k-\left\lfloor\frac{k}{2}\right\rfloor}\beta_{n_{2},\dots,n_{d+1}}\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1}\right)^{k}\end{split}$

where

$\begin{split}\beta_{n_{2},\dots,n_{d+1}}=&\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq
k\\}}a_{\alpha}c_{n_{2},\dots,n_{d+1}}(\alpha)\leq\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k\\}}|a_{\alpha}|\max\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k\\}}|c_{n_{2},\dots,n_{d+1}}(\alpha)|\\\
\leq&\left(\frac{k}{2}+1\right)^{2(d+1)}\left(\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k\\}}|a_{\alpha}|\right).\end{split}$

Thus, the polynomial $P$ can be written as

$\begin{split}P(x)=\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{k}{2}\right\rfloor}^{k-\left\lfloor\frac{k}{2}\right\rfloor}\beta_{n_{2},\dots,n_{d+1}}&\left[\sigma_{k}\left(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1}\right)\right.\\\
&\ +\left.(-1)^{k}\sigma_{k}\left(-x_{1}-n_{2}x_{2}-\dots-
n_{d}x_{d}-n_{d+1}\right)\right].\end{split}$

By noticing the inequality $\sigma_{k}(ay)=a^{k}\sigma_{k}(y)$ for all $a>0$,
$y\in\mathbb{R}$, we can write $P$ as

$\begin{split}P(x)=\sum\limits_{n_{2},\dots,n_{d}=-\left\lfloor\frac{k}{2}\right\rfloor}^{k-\left\lfloor\frac{k}{2}\right\rfloor}&B^{-k}\left\lceil\frac{k}{2}\right\rceil^{k}\beta_{n_{2},\dots,n_{d+1}}\left[\sigma_{k}\left(B\left\lceil\frac{k}{2}\right\rceil^{-1}(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1})\right)\right.\\\
&\
+\left.(-1)^{k}\sigma_{k}\left(-B\left\lceil\frac{k}{2}\right\rceil^{-1}(x_{1}+n_{2}x_{2}+\dots+n_{d}x_{d}+n_{d+1})\right)\right].\end{split}$

Clearly, the parameters satisfy

$\left|B\left\lceil\frac{k}{2}\right\rceil^{-1}\right|\leq
B\quad\text{and}\quad\left|B\left\lceil\frac{k}{2}\right\rceil^{-1}n_{j}\right|\leq
B$

for all $j=2,\dots,d+1$, and

$\left|B^{-k}\left\lceil\frac{k}{2}\right\rceil^{k}\beta_{n_{2},\dots,n_{d+1}}\right|\leq
B^{-k}\left(\frac{k}{2}+1\right)^{2(d+1)+k}\left(\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k\\}}|a_{\alpha}|\right)$

for all $-\left\lfloor\frac{k}{2}\right\rfloor\leq n_{2},\dots,n_{d}\leq
k-\left\lfloor\frac{k}{2}\right\rfloor$. This proves Lemma 2. ∎

Now, we are ready to prove our first main result, which reveals the way to
construct polynomials by deep ReLUk neural networks, as well as an estimation
of the upper bound of the parameters.

######  Theorem 1.

Let $d,k,L\in\mathbb{N}^{+}$, $B>0$, and $n_{1}=\dots=n_{L}=2(k^{L}+1)^{d}$,
then any polynomial

$P(x)=\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\ \|\alpha\|_{1}\leq
k^{L}\\}}a_{\alpha}x^{\alpha},\ x\in\mathbb{B}^{d}$

with degree $\leq k^{L}$ can be represented by a deep ReLUk network
architecture $\Sigma_{n_{1:L},M}^{k}(B)$ with

$M=B^{-k\frac{k^{L}-1}{k-1}}\left(\frac{k^{L}}{2}+1\right)^{2(d+1)+k^{L}}\left(\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k^{L}\\}}|a_{\alpha}|\right).$

Furthermore, there are $2(2L+d)(k^{L}+1)^{d}$ nonzero parameters in the neural
network architecture.

###### Proof.

By Lemma 2, $P$ can be represented as

$P(x)=\sum\limits_{j=1}^{2(k^{L}+1)^{d}}c_{j}\sigma_{k^{L}}(w_{j}\cdot
x+b_{j}),$ (2.7)

where $w_{j}\in[-1,1]^{d}$, $b_{j}\in[-1,1]$ and
$c_{j}\in[-M^{\prime},M^{\prime}]$ for $j=1,\dots,2(k^{L}+1)^{d}$, where

$M^{\prime}=\left(\frac{k^{L}}{2}+1\right)^{2(d+1)+k^{L}}\left(\sum\limits_{\\{\alpha\in\mathbb{N}^{d}:\
\|\alpha\|_{1}\leq k^{L}\\}}|a_{\alpha}|\right).$

Given $h_{0}(x)=x$, there exists a matrix $A_{1}$ with $2d(k^{L}+1)^{d}$
parameters and bias term $b=(b_{1},\dots,b_{2d(k^{L}+1)})$ such that

$h_{1,j}(x)=\sigma_{k}(Bw_{j}\cdot x+Bb_{j})=B^{k}\sigma_{k}(w_{j}\cdot
x+b_{j}),$

$j=1,\dots,2(k^{L}+1)^{d}$.

Given

$h_{i}(x)=B^{k\frac{k^{i}-1}{k-1}}\left(\sigma_{k^{i}}(w_{1}\cdot
x+b_{1}),\dots,\sigma_{k^{i}}\left(w_{2(k^{L}+1)^{d}}\cdot
x+b_{2(k^{L}+1)^{d}}\right)\right)^{\top},$

the identity matrix $A_{i}=BI_{2(k^{L}+1)^{d}\times 2(k^{L}+1)^{d}}$ with
$2(k^{L}+1)^{d}$ parameters and bias term $b=(0,\dots,0)$ give

$\begin{split}h_{i+1}(x)=&\sigma_{k}(A_{i}h_{i}(x))=\sigma_{k}(B\times
B^{k\frac{k^{i}-1}{k-1}}h_{i}(x)),\\\
=&\left(\sigma_{k}\left(B^{\frac{k^{i+1}-1}{k-1}}\sigma_{k^{i}}(w_{1}\cdot
x+b_{1})\right),\dots,\sigma_{k}\left(B^{\frac{k^{i+1}-1}{k-1}}\sigma_{k^{i}}\left(w_{2(k^{L}+1)^{d}}\cdot
x+b_{2(k^{L}+1)^{d}}\right)\right)\right)^{\top},\\\
=&B^{k\frac{k^{i+1}-1}{k-1}}\left(\sigma_{k^{i+1}}(w_{1}\cdot
x+b_{1}),\dots,\sigma_{k^{i+1}}\left(w_{2(k^{L}+1)^{d}}\cdot
x+b_{2(k^{L}+1)^{d}}\right)\right)^{\top}.\end{split}$

Inductively, we have
$h_{L}(x)=B^{k\frac{k^{L}-1}{k-1}}\left(\sigma_{k^{L}}(w_{1}\cdot
x+b_{1}),\dots,\sigma_{k^{L}}\left(w_{2(k^{L}+1)^{d}}\cdot
x+b_{2(k^{L}+1)^{d}}\right)\right)^{\top}$. Thus, with
$\tilde{c}=B^{-k\frac{k^{L}-1}{k-1}}(c_{1},\dots,c_{2(k^{L}+1)^{d}})$, we have

$\tilde{c}\cdot
h_{L}(x)=\sum\limits_{j=1}^{2(k^{L}+1)^{d}}c_{j}\sigma_{k^{L}}(w_{j}\cdot
x+u_{j})=P(x).$

Here, we have

$\|\tilde{c}\|_{\infty}\leq B^{-k\frac{k^{L}-1}{k-1}}M^{\prime}=M.$

This implies $P\in\Sigma_{n_{1:L}}^{k}(M)$ with
$n_{1}=\dots=n_{L}=2(k^{L}+1)^{d}$.

There are $2d(k^{L}+1)^{d}$ parameters in $A_{1}$ and $2(k^{L}+1)^{d}$
parameters in $A_{2},\dots,A_{L}$. The number of parameters in bias terms is
$L\times 2(k^{L}+1)^{d}$, and $2(k^{L}+1)^{d}$ parameters in $\tilde{c}$.

∎

######  Remark 2.

As we mentioned before, the existence of such representation of the
polynomials has been mentioned in previous literature such as [55, 4]. Those
results rely on a generalized Vandermonde matrix theory for high-dimensional
cases, while our proof in Lemma 1 reduces the problem into $1$-dimension
Vandermonde matrices. This improvement allows an estimation of the
coefficients in the linear combination (2.1), and hence the parameters in the
neural networks.

As a consequence of the classical polynomial approximation theory (e.g. [12,
Chapter 7]), we can deduce the rate of approximation for analytic functions
and for functions from Sobolev spaces.

######  Definition 3 (Sobolev spaces).

Let $d\in\mathbb{N}^{+}$ and $1\leq p\leq\infty$. For $r\in\mathbb{N}$, the
Sobolev space $W^{r,p}(\mathbb{B}^{d})$ is defined as

$W^{r,p}(\mathbb{B}^{d})=\left\\{f\in L^{p}(\mathbb{B}^{d}):\
\max\limits_{\\{\alpha\in\mathbb{N}^{d},\|\alpha\|_{1}=r\\}}\|D^{\alpha}f\|_{L^{p}(\mathbb{B}^{d})}<\infty\right\\}$
(2.8)

with the norm given by

$\|f\|_{W^{r,p}(\mathbb{B}^{d})}^{p}=\|f\|_{L^{p}(\mathbb{B}^{d})}^{p}+\sum\limits_{\\{\alpha\in\mathbb{N}^{d},\|\alpha\|_{1}=r\\}}\|D^{\alpha}f\|_{L^{p}(\mathbb{B}^{d})}^{p},\quad
f\in W^{r,p}(\mathbb{B}^{d}).$

For $r\notin\mathbb{N}$, let $\theta=r-\lfloor r\rfloor$. The Sobolev semi-
norm is

$|f|_{W^{r,p}}(\mathbb{B}^{d})=\left\\{\begin{array}[]{ll}\max\limits_{\\{\alpha\in\mathbb{N}^{d},\|\alpha\|_{1}=r\\}}\left(\displaystyle\int_{\mathbb{B}^{d}\times\mathbb{B}^{d}}\frac{|D^{\alpha}f(x)-D^{\alpha}f(y)|^{p}}{|x-y|^{d+\theta
p}}dxdy\right)^{1/p},&\quad 1\leq p<\infty,\\\
\max\limits_{\\{\alpha\in\mathbb{N}^{d},\|\alpha\|_{1}=r\\}}\sup\limits_{x,y\in\mathbb{B}^{d},x\neq
y}\displaystyle\frac{|D^{\alpha}f(x)-D^{\alpha}f(y)|}{|x-y|^{\theta}},&\quad
p=\infty.\end{array}\right.$

The space $W^{r,p}(\mathbb{B}^{d})$ is defined as

$W^{r,p}(\mathbb{B}^{d})=\left\\{f\in W^{\lfloor r\rfloor,p}(\mathbb{B}^{d}):\
|f|_{W^{r,p}}(\mathbb{B}^{d})<\infty\right\\}$ (2.9)

with the norm

$\|f\|_{W^{r,p}(\mathbb{B}^{d})}=\|f\|_{W^{\lfloor
r\rfloor,p}(\mathbb{B}^{d})}+|f|_{W^{r,p}}(\mathbb{B}^{d}).$

######  Definition 4 (Analytic functions).

Let $d\in\mathbb{N}^{+}$, $\rho\in(0,1)$, $f$ is said to be an analytic
function on $U_{\rho}:=\left\\{z\in\mathbb{C}^{d}:\
\left|z_{j}+\sqrt{z_{j}-1}\right|<\rho^{-1},\ j=1,\dots,d\right\\}$ if it is
complex differentiable at each $z\in U_{\rho}$.

######  Corollary 1.

Let $d,k,L\in\mathbb{N}^{+}$, $B>0$, and

$n_{1}=\dots=n_{L}=2(k^{L}+1)^{d}.$

There exists a sparse deep ReLUk network $\Sigma_{n_{1:L}}^{k}(B)$ with
$2(2L+d)(k^{L}+1)^{d}$ parameters such that

  1. (a)

if $f\in W^{r,p}(\mathbb{B}^{d})$ for some $r>0$, then

$\inf\limits_{g\in\Sigma_{n_{1:L}}^{k}(B)}\|f-g\|_{L^{p}(\mathbb{B}^{d})}\leq\inf\limits_{g\in\mathcal{P}_{d}^{k^{L}}}\|f-g\|_{L^{p}(\mathbb{B}^{d})}\sim\|f\|_{W^{r,p}(\mathbb{B}^{d})}k^{-Lr};$
(2.10)

  2. (b)

if $f$ is analytic on $U_{\rho}:=\left\\{z\in\mathbb{C}^{d}:\
\left|z_{j}+\sqrt{z_{j}-1}\right|<\rho^{-1},\ j=1,\dots,d\right\\}$, then for
any $\epsilon>0$ with $\rho+\epsilon<1$,

$\inf\limits_{g\in\Sigma_{n_{1:L}}^{k}(B)}\|f-g\|_{L^{p}(\mathbb{B}^{d})}\leq\inf\limits_{g\in\mathcal{P}_{d}^{k^{L}}}\|f-g\|_{L^{p}(\mathbb{B}^{d})}\lesssim\left(\max\limits_{w\in\partial\Gamma_{\rho+\epsilon}}\left|f\left(\frac{w+w^{-1}}{2}\right)\right|\right)k^{Ld}(\rho+\epsilon)^{k^{L}},$
(2.11)

where $\Gamma_{\rho+\epsilon}:=\left\\{w\in\mathbb{C}^{d}:\
|w_{j}|\leq(\rho+\epsilon)^{-1},\ j=1,\dots,d\right\\}$. The corresponding
constant depends only on $d,\rho$, and $\epsilon$.

###### Proof.

  1. (a)

The first inequality is exactly Theorem 1. The relation

$\inf\limits_{g\in\mathcal{P}_{d}^{k^{L}}}\|f-g\|_{L^{p}(\mathbb{B}^{d})}\sim\|f\|_{W^{r,p}(\mathbb{B}^{d})}k^{-Lr}$

is the classical polynomial approximation theory (see, e.g., [12, Chapter 7]).

  2. (b)

Again, the first inequality is Theorem 1. For the second inequality, we apply
the classical estimation for analytic functions (see, e.g., [12, Chapter 7]
and [13, Lemma A.1] for the multivariate version) and conclude

$\begin{split}\inf\limits_{g\in\mathcal{P}_{d}^{k^{L}}}\|f-g\|_{L^{2}(\mathbb{B}^{d})}^{2}\leq&\left(\max\limits_{w\in\partial
U_{\rho+\epsilon}}|f(w+w^{-1})|\right)^{2}\sum\limits_{n=k^{L}+1}^{\infty}\binom{n+d-1}{d-1}(\rho+\epsilon)^{2n}\\\
\leq&C^{2}\left(\max\limits_{w\in\partial
U_{\rho+\epsilon}}|f(w+w^{-1})|\right)^{2}k^{2L}(\rho+\epsilon)^{2k^{L}},\end{split}$

where $C=C(d,\rho,\epsilon)$ is independent of $f$ and $L$.

∎

##  3 Adaptive Approximation for Functions from Variation Spaces

In this section, we show the adaptive approximation properties of deep ${\rm
ReLU}^{k}$ neural networks for functions in variation spaces. This result is
obtained by studying the representation power of deep ${\rm ReLU}^{k}$ neural
networks for shallow networks.

######  Lemma 3.

Let $d,k,L,n\in\mathbb{N}^{+}$, $M>0$, and $n_{1}=\dots=n_{L}=2(k+1)n$. Then
there exists a deep ReLUk neural network architecture
$\Sigma_{n_{1:L},M^{\prime}}^{k}(B^{\prime})$ with

$N=[(4L-2)(k+1)+d]n$

nonzero parameters, such that all functions $f\in\Sigma^{K}_{n,M}(B)$ with
$K\in\\{k^{1},\dots,k^{L}\\}$ can be represented this neural network.
Furthermore, the parameters in $\Sigma_{n_{1:L},M^{\prime}}^{k}(B^{\prime})$
are bounded as

$M^{\prime}=\left(\frac{k}{2}+1\right)^{4}M,\quad
B^{\prime}=B\vee\left(\frac{k}{2}+1\right)^{4}.$

###### Proof.

Let $1\leq\ell\leq L$, $K=k^{\ell}$, and $f\in\Sigma^{K}_{n,M}(B)$. Then

$f_{n}(x)=\sum\limits_{m=1}^{n}c_{m}\sigma_{K}(w_{m}\cdot x+u_{m}),\quad
x\in\mathbb{B}^{d},$

where $|c_{m}|\leq M$, $|u_{m}|\leq B$, and $\|w_{m}\|_{\infty}\leq B$ holds
for $m=1,\dots,n$.

Given $h_{0}(x)=x$, there exists a sparse matrix
$A_{0}\in\mathbb{R}^{n_{1}\times d}$ with $(A_{0})_{i,j}=0$ for $i>n$ and
$b_{0}=(u_{1},\dots,u_{n})^{\top}$ such that

$h_{1,j}(x)=\left\\{\begin{array}[]{ll}\sigma_{k}\left(w_{j}\cdot
x+u_{j}\right),&\quad j\leq n,\\\ 0,&\quad\hbox{otherwise.}\end{array}\right.$
(3.1)

Clearly, all the parameters in $A_{0}$ and $b_{0}$ are bounded by
$B^{\prime}$.

Given

$h_{i,j}(x)=\left\\{\begin{array}[]{ll}\sigma_{k^{i}}(w_{j}\cdot
x+u_{j}),&\quad j\leq n,\\\ 0,&\quad\hbox{otherwise,}\end{array}\right.$

let $A_{i}=B^{\prime}I_{2(k+1)n\times 2(k+1)n}$, and
$b_{i}=(0,\dots,0)^{\top}$,

$h_{i+1,j}(x)=\left\\{\begin{array}[]{ll}\sigma_{k^{i+1}}(w_{j}\cdot
x+u_{j}),&\quad j\leq n,\\\ 0,&\quad\hbox{otherwise.}\end{array}\right.$

So we get

$h_{\ell,j}(x)=\left\\{\begin{array}[]{ll}\sigma_{K}(w_{j}\cdot
x+u_{j}),&\quad j\leq n,\\\ 0,&\quad\hbox{otherwise.}\end{array}\right.$

By (2.6), the function $y\mapsto y$ can be represented as a linear combination

$y=\sum\limits_{t=0}^{k}a_{t}\left(y+t-\left\lfloor\frac{k}{2}\right\rfloor\right)^{k}=\sum\limits_{t=0}^{k}a_{t}\left[\sigma_{k}\left(y+t-\left\lfloor\frac{k}{2}\right\rfloor\right)+(-1)^{k}\sigma_{k}\left(-y-t+\left\lfloor\frac{k}{2}\right\rfloor\right)\right]$
(3.2)

with $|a_{t}|\leq\left(\frac{k}{2}+1\right)^{4}$.

At the $\ell$-th layer, there exists a sparse matrix
$A_{\ell}\in\mathbb{R}^{n_{\ell+1}\times n_{\ell}}$ with

$\begin{split}&(A_{\ell})_{2(k+1)(m-1)+2t+1,m}=1,\\\
&(A_{\ell})_{2(k+1)(m-1)+2t+2,m}=-1,\end{split}$

for $m=1,\dots,n,\ t=0,\dots,k$ and proper bias term $b_{\ell}$ with
$\|b_{\ell}\|_{\infty}\leq\frac{k}{2}+1$ such that

$\begin{split}&h_{\ell+1,2(k+1)(m-1)+2t+1}(x)=\sigma_{k}\left(\sigma_{L}(w_{m}\cdot
x+u_{m})+t-\left\lfloor\frac{k}{2}\right\rfloor\right)\\\
&h_{\ell+1,2(k+1)(m-1)+2t+2}(x)=\sigma_{k}\left(-\sigma_{L}(w_{m}\cdot
x+u_{m})-t+\left\lfloor\frac{k}{2}\right\rfloor\right).\end{split}$

for $m=1,\dots,n,\ t=0,\dots,k$.

At the $(\ell+1)$-th layer, there exists a sparse matrix
$A_{\ell+1}\in\mathbb{R}^{n_{\ell+2}\times n_{\ell+1}}$ with

$\begin{split}&(A_{\ell+1})_{m,2(k+1)(m-1)+2t+1}=a_{t},\\\
&(A_{\ell+1})_{m,2(k+1)(m-1)+2t+2}=(-1)^{k}a_{t}\end{split}$

and $b_{\ell+1}=(0,\dots,0)^{\top}$ such that

$\begin{split}h_{\ell+2,j}(x)=&\sum\limits_{t=0}^{k}a_{t}\left[\sigma_{k}\left(\sigma_{L}(w_{j}\cdot
x+u_{j})+t-\left\lfloor\frac{k}{2}\right\rfloor\right)+(-1)^{k}\sigma_{k}\left(-\sigma_{L}(w_{j}\cdot
x+u_{j})-t+\left\lfloor\frac{k}{2}\right\rfloor\right)\right]\\\
=&\sigma_{L}(w_{j}\cdot x+u_{j}),\quad j\leq n,\\\
h_{\ell+2,j}(x)=&0,\quad\hbox{otherwise.}\end{split}$

We repeat these two steps up to the $L$-th layer, then the last hidden layer
$h_{L}(x)$ is given as either

$\begin{split}h_{L,j}(x)=&\sigma_{K}(w_{j}\cdot x+u_{j}),\quad j\leq n,\\\
h_{L,j}(x)=&0,\quad\hbox{otherwise}\end{split}$

or

$\begin{split}&h_{L,2(k+1)(m-1)+2t+1}(x)=\sigma_{k}\left(\sigma_{K}(w_{m}\cdot
x+u_{m})+t\right)\\\
&h_{L,2(k+1)(m-1)+2t+2}(x)=\sigma_{k}\left(-\sigma_{K}(w_{m}\cdot
x+u_{m})-t\right).\end{split}$

for $m=1,\dots,n,\ t=0,\dots,k$.

In either case, by (3.2), there exists
$c\in\left[-\left(\frac{k}{2}+1\right)^{4},\left(\frac{k}{2}+1\right)^{4}\right]^{2(k+1)m}$
such that

$c\cdot h_{L}(x)=\sum\limits_{m=1}^{n}c_{m}\sigma_{K}(w_{m}\cdot
x+u_{m})=f_{n}(x).$

By the form of the matrices $A_{0},\dots,A_{L-1}$, the number of parameters is

$dn+2(k+1)n+\sum\limits_{j=2}^{L}\left[2(k+1)n+2(k+1)n\right]+2(k+1)n=N.$

This completes the proof of Lemma 3.

∎

This lemma enables us to consider the rate of approximation for functions from
the variation spaces generated by ReLUK activation functions with
$K=k,k^{2},\dots,k^{L}$. The approximation theory of such spaces has been
studied in [51, 34, 27], etc.

######  Definition 5 (Variation space).

Let $K,d\in\mathbb{N}^{+}$, denote the dictionary $\mathbb{P}_{K}^{d}$ as

$\mathbb{P}_{K}^{d}=\left\\{x\mapsto\sigma_{K}(\omega\cdot x+b):\
\omega\in\mathbb{S}^{d-1},\ b\in[-1,1]\right\\}.$

Consider the closure of the convex, symmetric hull of $\mathbb{P}_{K}^{d}$

$\overline{\mathrm{Conv}(\mathbb{P}_{K}^{d})}=\overline{\left\\{\sum\limits_{j=1}^{n}a_{j}h_{j}:\
n\in\mathbb{N}^{+},\ \sum\limits_{j=1}^{n}|a_{j}|\leq 1,\
h_{j}\in\mathbb{P}_{K}^{d}\right\\}},$ (3.3)

we define the variation space with the dictionary $\mathbb{P}_{L}^{d}$ as

$\mathcal{K}_{1}(\mathbb{P}_{K}^{d}):=\left\\{f\in
R\times\overline{\mathrm{Conv}(\mathbb{P}_{K}^{d})}:\ R>0\right\\}.$ (3.4)

with the norm

$\|f\|_{\mathcal{K}_{1}(\mathbb{P}_{K}^{d})}=\inf\left\\{R>0:\ f\in
R\times\overline{\mathrm{Conv}(\mathbb{P}_{K}^{d})}\right\\}.$ (3.5)

We will compare the result with the rate of approximation of shallow ReLUk
networks given in [51].

######  Theorem 2 (Sigel and Xu, 2022).

Let $k,d\in\mathbb{N}^{+}$, then

$\sup\limits_{f\in\mathcal{K}_{1}(\mathbb{P}_{k}^{d})}\inf\limits_{f_{n}\in\Sigma_{n}^{k}(2)}\|f-f_{n}\|_{L^{2}(\mathbb{B}^{d})}\lesssim\|f\|_{\mathcal{K}_{1}(\mathbb{P}_{k}^{d})}n^{-\frac{1}{2}-\frac{2k+1}{2d}}.$
(3.6)

Combining this theorem with Lemma 3, we have the following theorem.

######  Theorem 3.

Let $k,n,L\in\mathbb{N}^{+}$, and $n_{1}=\dots=n_{L}=2(k+1)n$, then there
exists a deep ReLUk network architecture $\Sigma_{n_{1:L}}^{k}(B)$ with
$B=\left(\frac{k}{2}+1\right)^{4}$ such that for any $K=k^{\ell}$ with
$1\leq\ell\leq L$,

$\sup\limits_{f\in\mathcal{K}_{1}(\mathbb{P}_{K}^{d})}\inf\limits_{f_{n}\in\Sigma_{n_{1:L}}^{k}(B)}\|f-f_{n}\|_{L^{2}(\mathbb{B}^{d})}\leq
C\|f\|_{\mathcal{K}_{1}(\mathbb{P}_{K}^{d})}n^{-\frac{1}{2}-\frac{2K+1}{2d}},$

where $N=[(4L-2)(k+1)+d]n$, $C$ is a constant independent of $n$.

######  Remark 3.

The index $K$ in the variation spaces defined can be correlated to the
regularity of the functions in the spaces (cf. [51]). Suppose the target
function has a large but unknown regularity, saying
$f\in\cap_{j=1}^{K}\mathcal{K}_{1}(\mathbb{P}_{j}^{d})$, then the
approximation property of shallow ReLUk networks would be bad if the
activation function is not chosen properly: If $k\ll K$, then the
approximation rate that we can achieve is at most

$\sup\limits_{f\in\mathcal{K}_{1}(\mathbb{P}_{k}^{d})}\inf\limits_{f_{n}\in\Sigma_{n}^{k}}\|f-f_{n}\|_{L^{2}(\mathbb{B}^{d})}\lesssim\|f\|_{\mathcal{K}_{1}(\mathbb{P}_{k}^{d})}n^{-\frac{1}{2}-\frac{2k+1}{2d}}.$

If $k>K$,
$\cap_{j=1}^{K}\mathcal{K}_{1}(\mathbb{P}_{j}^{d})\not\subseteq\mathcal{K}_{1}(\mathbb{P}_{k}^{d})$,
the approximation rate is unknown.

However, by taking the depth $L\geq\lceil\log_{k}K\rceil$, the deep ReLUk
networks we constructed in this section can automatically approximate
functions from $\cap_{j=1}^{K}\mathcal{K}_{1}(\mathbb{P}_{j}^{d})$ with this
unknown $K$. The rate of approximation is at least

$\sup\limits_{f\in\mathcal{K}_{1}(\mathbb{P}_{K}^{d})}\inf\limits_{f_{n}\in\Sigma_{n_{1:L}}^{k}}\|f-f_{n}\|_{L^{2}(\mathbb{B}^{d})}\leq
C\|f\|_{\mathcal{K}_{1}(\mathbb{P}_{K}^{d})}n^{-\frac{1}{2}-\frac{2K+1}{2d}},$

since $k^{\lceil\log_{k}K\rceil}\geq K$.

##  4 Conclusions

In conclusion, this paper has significantly advanced the understanding of
${\rm ReLU}^{k}$ neural networks and their capabilities in function
approximation and representation. Our constructive proofs have not only
demonstrated the representational power of shallow ${\rm ReLU}^{k}$ networks
in polynomial approximation but also paved the way for constructing deep ${\rm
ReLU}^{k}$ networks with an explicit parameterization for any polynomial of
degree less than $k^{L}$, scaling efficiently with network depth and
dimensions. Furthermore, we established a connection between the coefficients
of polynomials and the bounds of neural network parameters, thereby enabling
precise estimations of network parameters that uphold the approximation
integrity. The elucidation of deep ${\rm ReLU}^{k}$ networks’ ability to
approximate functions from Sobolev spaces, albeit with a suboptimal rate,
opens new avenues for exploring their application in complex functional
spaces. Most notably, our construction that deep ${\rm ReLU}^{k}$ networks can
emulate the function representation of shallower networks with enhanced
efficiency underscores the hierarchical strength of depth in neural
architectures. Coupled with recent findings on variation spaces, our results
underscore the adaptivity and potential of deep ${\rm ReLU}^{k}$ networks to
achieve significant approximation accuracy, contributing to the theoretical
foundations that may stimulate further breakthroughs in the application of
deep learning across varied domains.

## References

  * [1] R. Arora, A. Basu, P. Mianjy, and A. Mukherjee.  Understanding deep neural networks with rectified linear units.  In International Conference on Learning Representations, 2018.
  * [2] A. R. Barron.  Universal approximation bounds for superpositions of a sigmoidal function.  IEEE Transactions on Information theory, 39(3):930–945, 1993.
  * [3] A. Bordes, S. Chopra, and J. Weston.  Question answering with subgraph embeddings.  arXiv preprint arXiv:1406.3676, 2014.
  * [4] Q. Chen, W. Hao, and J. He.  Power series expansion neural network.  Journal of Computational Science, 59:101552, 2022.
  * [5] C. K. Chui and X. Li.  Approximation by ridge functions and neural networks with one hidden layer.  Journal of Approximation Theory, 70(2):131–141, 1992.
  * [6] C. K. Chui, X. Li, and H. N. Mhaskar.  Limitations of the approximation capabilities of neural networks with one hidden layer.  Advances in Computational Mathematics, 5:233–243, 1996.
  * [7] C. K. Chui, S.-B. Lin, and D.-X. Zhou.  Deep neural networks for rotation-invariance approximation and learning.  Analysis and Applications, 17(05):737–772, 2019.
  * [8] R. Collobert, J. Weston, L. Bottou, M. Karlen, K. Kavukcuoglu, and P. Kuksa.  Natural language processing (almost) from scratch.  Journal of machine learning research, 12(ARTICLE):2493–2537, 2011\.
  * [9] F. Cucker and D. X. Zhou.  Learning theory: an approximation theory viewpoint, volume 24.  Cambridge University Press, 2007.
  * [10] G. Cybenko.  Approximation by superpositions of a sigmoidal function.  Mathematics of control, signals and systems, 2(4):303–314, 1989\.
  * [11] O. Delalleau and Y. Bengio.  Shallow vs. deep sum-product networks.  Advances in neural information processing systems, 24, 2011.
  * [12] R. A. DeVore and G. G. Lorentz.  Constructive approximation, volume 303.  Springer Science & Business Media, 1993.
  * [13] K. Doctor, T. Mao, and H. Mhaskar.  Encoding of data sets and algorithms.  arXiv preprint arXiv:2303.00984, 2023.
  * [14] R. Eldan and O. Shamir.  The power of depth for feedforward neural networks.  In Conference on learning theory, pages 907–940. PMLR, 2016.
  * [15] Z. Fang, H. Feng, S. Huang, and D.-X. Zhou.  Theory of deep convolutional neural networks ii: Spherical analysis.  Neural Networks, 131:154–162, 2020.
  * [16] W. Gautschi.  On inverses of vandermonde and confluent vandermonde matrices iii.  Numerische Mathematik, 29:445–450, 1978.
  * [17] J. He.  On the optimal expressive power of relu dnns and its application in approximation with kolmogorov superposition theorem.  arXiv preprint arXiv:2308.05509, 2023.
  * [18] J. He, L. Li, and J. Xu.  Approximation properties of deep relu cnns.  Research in the mathematical sciences, 9(3):38, 2022.
  * [19] J. He, L. Li, and J. Xu.  Relu deep neural networks from the hierarchical basis perspective.  Computers & Mathematics with Applications, 120:105–114, 2022.
  * [20] J. He, L. Li, J. Xu, and C. Zheng.  Relu deep neural networks and linear finite elements.  Journal of Computational Mathematics, 38(3):502–527, 2020.
  * [21] J. He and J. Xu.  Mgnet: A unified framework of multigrid and convolutional neural network.  Science china mathematics, 62:1331–1354, 2019.
  * [22] J. He and J. Xu.  Deep neural networks and finite elements of any order on arbitrary dimensions.  arXiv preprint arXiv:2312.14276, 2023.
  * [23] K. He, X. Zhang, S. Ren, and J. Sun.  Deep residual learning for image recognition.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 770–778, 2016.
  * [24] G. Hinton, L. Deng, D. Yu, G. E. Dahl, A.-r. Mohamed, N. Jaitly, A. Senior, V. Vanhoucke, P. Nguyen, T. N. Sainath, et al.  Deep neural networks for acoustic modeling in speech recognition: The shared views of four research groups.  IEEE Signal processing magazine, 29(6):82–97, 2012.
  * [25] K. Hornik, M. Stinchcombe, and H. White.  Multilayer feedforward networks are universal approximators.  Neural networks, 2(5):359–366, 1989.
  * [26] S. Jean, K. Cho, R. Memisevic, and Y. Bengio.  On using very large target vocabulary for neural machine translation.  arXiv preprint arXiv:1412.2007, 2014.
  * [27] J. M. Klusowski and A. R. Barron.  Approximation by combinations of ReLU and squared ReLU ridge functions with $\ell^{1}$ and $\ell^{0}$ controls.  IEEE Transactions on Information Theory, 64(12):7649–7656, 2018\.
  * [28] A. N. Kolmogorov.  On the representation of continuous functions of many variables by superposition of continuous functions of one variable and addition.  In Doklady Akademii Nauk, volume 114, pages 953–956. Russian Academy of Sciences, 1957.
  * [29] A. Krizhevsky, I. Sutskever, and G. E. Hinton.  Imagenet classification with deep convolutional neural networks.  Advances in neural information processing systems, 25, 2012.
  * [30] Y. LeCun, L. Bottou, Y. Bengio, and P. Haffner.  Gradient-based learning applied to document recognition.  Proceedings of the IEEE, 86(11):2278–2324, 1998.
  * [31] M. Leshno, V. Y. Lin, A. Pinkus, and S. Schocken.  Multilayer feedforward networks with a nonpolynomial activation function can approximate any function.  Neural networks, 6(6):861–867, 1993.
  * [32] B. Li, S. Tang, and H. Yu.  Better approximations of high dimensional smooth functions by deep neural networks with rectified power units.  Communications in Computational Physics, 27(2):379–411, 2019.
  * [33] Z. Lu, H. Pu, F. Wang, Z. Hu, and L. Wang.  The expressive power of neural networks: A view from the width.  Advances in neural information processing systems, 30, 2017.
  * [34] L. Ma, J. W. Siegel, and J. Xu.  Uniform approximation rates and metric entropy of shallow neural networks.  Research in the Mathematical Sciences, 9(3):46, 2022.
  * [35] T. Mao, Z. Shi, and D.-X. Zhou.  Theory of deep convolutional neural networks iii: Approximating radial functions.  Neural Networks, 144:778–790, 2021.
  * [36] T. Mao, Z. Shi, and D.-X. Zhou.  Approximating functions with multi-features by deep convolutional neural networks.  Analysis and Applications, 21(01):93–125, 2023.
  * [37] T. Mao and D.-X. Zhou.  Rates of approximation by ReLU shallow neural networks.  Journal of Complexity, 79:101784, 2023.
  * [38] H. N. Mhaskar.  Approximation properties of a multilayered feedforward artificial neural network.  Advances in Computational Mathematics, 1:61–80, 1993.
  * [39] H. N. Mhaskar and C. A. Micchelli.  Degree of approximation by neural and translation networks with a single hidden layer.  Advances in applied mathematics, 16(2):151–183, 1995.
  * [40] T. Mikolov, A. Deoras, D. Povey, L. Burget, and J. Černockỳ.  Strategies for training large scale neural network language models.  In 2011 IEEE Workshop on Automatic Speech Recognition & Understanding, pages 196–201. IEEE, 2011.
  * [41] H. Montanelli and Q. Du.  New error bounds for deep ReLU networks using sparse grids.  SIAM Journal on Mathematics of Data Science, 1(1):78–92, 2019.
  * [42] J. A. Opschoor, C. Schwab, and J. Zech.  Exponential relu dnn expression of holomorphic maps in high dimension.  Constructive Approximation, 55(1):537–582, 2022.
  * [43] T. Poggio, H. Mhaskar, L. Rosasco, B. Miranda, and Q. Liao.  Why and when can deep-but not shallow-networks avoid the curse of dimensionality: a review.  International Journal of Automation and Computing, 14(5):503–519, 2017.
  * [44] I. Safran and O. Shamir.  Depth separation in ReLU networks for approximating smooth non-linear functions.  arXiv preprint arXiv:1610.09887, 14, 2016.
  * [45] T. N. Sainath, B. Kingsbury, A.-r. Mohamed, G. E. Dahl, G. Saon, H. Soltau, T. Beran, A. Y. Aravkin, and B. Ramabhadran.  Improvements to deep convolutional neural networks for lvcsr.  In 2013 IEEE workshop on automatic speech recognition and understanding, pages 315–320. IEEE, 2013.
  * [46] Z. Shen, H. Yang, and S. Zhang.  Deep network approximation characterized by number of neurons.  arXiv preprint arXiv:1906.05497, 2019.
  * [47] Z. Shen, H. Yang, and S. Zhang.  Optimal approximation rate of ReLU networks in terms of width and depth.  Journal de Mathématiques Pures et Appliquées, 157:101–135, 2022.
  * [48] J. W. Siegel.  Optimal approximation rates for deep ReLU neural networks on sobolev spaces.  arXiv preprint arXiv:2211.14400, 2022.
  * [49] J. W. Siegel and J. Xu.  Optimal approximation rates and metric entropy of ReLUk and cosine networks.  arXiv preprint arXiv:2101.12365, 2021.
  * [50] J. W. Siegel and J. Xu.  High-order approximation rates for shallow neural networks with cosine and ReLUk activation functions.  Applied and Computational Harmonic Analysis, 58:1–26, 2022.
  * [51] J. W. Siegel and J. Xu.  Sharp bounds on the approximation rates, metric entropy, and n-widths of shallow neural networks.  Foundations of Computational Mathematics, pages 1–57, 2022.
  * [52] L. Song, J. Fan, D.-R. Chen, and D.-X. Zhou.  Approximation of nonlinear functionals using deep ReLU networks.  arXiv preprint arXiv:2304.04443, 2023.
  * [53] I. Sutskever, O. Vinyals, and Q. V. Le.  Sequence to sequence learning with neural networks.  Advances in neural information processing systems, 27, 2014.
  * [54] M. Telgarsky.  Representation benefits of deep feedforward networks.  arXiv preprint arXiv:1509.08101, 2015.
  * [55] J. Xu.  The finite neuron method and convergence analysis.  arXiv preprint arXiv:2010.01458, 2020.
  * [56] Y. Yang and D.-X. Zhou.  Optimal rates of approximation by shallow ReLUk neural networks and applications to nonparametric regression.  arXiv preprint arXiv:2304.01561, 2023.
  * [57] D. Yarotsky.  Error bounds for approximations with deep ReLU networks.  Neural Networks, 94:103–114, 2017.
  * [58] D.-X. Zhou.  Universality of deep convolutional neural networks.  Applied and computational harmonic analysis, 48(2):787–794, 2020\.