Nonparametric Quantile Regression with ReLU-Activated Recurrent Neural Networks
Authors: Hang Yu
Lyumin Wu
Wenxin Zhou
Zhao Ren
Venue: NeurIPS 2025 poster ()
Keywords: nonparametric regression; quantile regression; neural networks; deep learning; non-asymptotic bounds; stationarity
Primary area: deep_learning

Abstract: This paper investigates nonparametric quantile regression using recurrent neural networks (RNNs) and sparse recurrent neural networks (SRNNs) to approximate the conditional quantile function, which is assumed to follow a compositional hierarchical interaction model. We show that RNN- and SRNN-based estimators with rectified linear unit (ReLU) activation and appropriately designed architectures achieve the optimal nonparametric convergence rate, up to a logarithmic factor, under stationary, exponentially $\boldsymbol{\beta}$-mixing processes. To establish this result, we derive sharp approximation error bounds for functions in the hierarchical interaction model using RNNs and SRNNs, exploiting their close connection to sparse feedforward neural networks (SFNNs). Numerical experiments and an empirical study on the Dow Jones Industrial Average (DJIA) further support our theoretical findings.

Nonparametric Quantile Regression with ReLU-Activated Recurrent Neural Networks
Hang Yu?.2, Lyumin Wu Wen-Xin Zhou", Zhao Ren 5.1
National Key Laboratory for Novel Software Technology, Nanjing University, China School of Artificial Intelligence, Nanjing University, China Department of Al and Data Science, The University of Hong Kong. China Department of Information and Decision Sciences, University of Illinois Chicago, USA 5 Department of Statistics, University of Pittsburgh, USA
Abstract
This paper investigates nonparametric quantile regression using recurrent neural networks (RNNs) and sparse recurrent neural networks (SRNNs) to approximate the conditional quantile function, which is assumed to follow a compositional hierarchical interaction model. We show that RNN- and SRNN-based estimators with rectified linear unit (ReLU) activation and appropriately designed architectures achieve the optimal nonparametrie convergence rate, up to logarithmic factor, under stationary, exponentially 3-mixing processes. To establish this result, we derive sharp approximation error bounds for functions in the hierarchical interac- tion model using RNNs and SRNNs, exploiting their close connection to sparse feedforward neural networks (SFNNs). Numerical experiments and empirical study on the Dow Jones Industrial Average (DJIA) further support our theoretical findings.
Introduction
Quantile regression (QR) [Koenker and Bassett, 1978] provides a flexible framework for estimating conditional quantiles of response variable, offering more comprehensive characterization of the conditional distribution than least squares regression, which focuses solely on the conditional mean. By modeling different quantiles, QR reveals how covariate effects vary across the response distribu- tion, making it particularly valuable when errors are non-normal, heteroscedastic, or when the focus is on tail behavior. Since its introduction, QR has evolved through wide range of methodological developments, including linear QR [Koenker and Bassett, 1978), quantile autoregression [Koenker and Xiao, 2006), quantile regression forests [Meinshausen, 2006], and quantile boosting [Zheng, 2012], among others. An early step toward integrating neural networks into QR was taken by White [1992], who established theoretical guarantees for single-layer feedforward networks. Building on this foundation, subsequent research has focused on multi-layer feedforward neural networks (FNNs) with ReLU activation functions [Nair and Hinton, 2010]. Assuming that the true conditional quantile function admits a compositional structure consisting of diver-dimensional component functions, Shen et [2021] derived a sub-optimal convergence rate for Resee FNN estimators in the presence of heavy-tailed response distributions, which was later refined in a subsequent work [Shen al., 2025]. Extending this line of work, Padilla al. [2022] investigated nonparametric QR using Rel.U-activated sparse feedforward neural networks (SFNNs) with bounded parameters, achieving optimal convergence rates under Hölder smoothness or Besov-space assumptions on the quantile
"Equal contribution.
"ooreeponeeece: Zhao Ren <zren @zren@pit.edu>>> pitt.edu>
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

function. More recently. Feng et al. [2024] addressed the problem of covariate shift in nonparametric QR via ReLU FNNs, obtaining minimax-ovn -onimaa- rates under an adaptive calibratton condition.
However, FNNs often struggle to capture temporal dynamics and sequential dependencies, thereby limiting their effectiveness in modeling dependent data. In contrast, recurrent neural networks (RNNs) [Rumelhart et al., 1986], specifically designed for sequential data, are naturally better suited to such settings. Their capacity retain and integrate information across time steps has proven effective in applications such as time series forecasting, sequence-to-sequeee and modeling long-term dependencies [Hochreiter and Schmidhuber, 1997, Chung et al., 2014]. Motivated by these strengths, we investigate RNNs and their sparse variants (SRNNs) predictive function classes and establish their theoretical properties within the framework of nonparametric QR.
learning [Sutskever et al., 2014],
Recent studies have investigated the theoretical properties of RNNs in the many-to-one setting. Several works have established PAC-style generalization guarantees for RNNs [Chen al., 2020, Tu et al., 2020, Cheng et al., 2025]. More closely related to our setting are results on their approximation properties under dependent data. In particular, Jíao al. [2024] showed that RNNs can attain the optimal convergence rate in nonparametric regression when the data are stationary and B- mixing. Building upon these developments, the present study establishes statistical guarantees for the performance of RNNs in QR, which requires distinct technical tools from those employed in least squares estimation. In what follows, we formally define the model, outline the main contributions, and introduce the notation.
Model. Consider sequentially stationary observations {(x, (1),))))111 where any consecutive N observations share the same joint distribution as Z = ((X1, Ki),...,
Yx)). Here, Y, E R dy-dimensional covariate vector for i [N]. Motivated by the recurrent structure of RNNs, we assume that Y1 depends on the sequence of covariates (Xt-N-1:.... ...........))
denotes the random outcome of interest, and X, E Rd
Xt). The formal regularity conditions governing this dependence are specified in Section 3.2. Given quantile level T E (0, 1) of interest, we define the conditional 7-th quantile of y/ (or YN) given N X₁ (or XI...., XN)
(Ue (ovlxx-An X+- xx) = N+1. .NI...,..,,,.. X,-N+1, X E R2
NSIS n,
where fo R + R is the unknown conditional quantile function. Equivalently, the relationship between YN and (X1,..., Xv) can be expressed in additive form
YN fo(,,....,, Xx) +
(1.1)
where e denotes the regression error satisfying P(e S ,I,
A key motivating example is the broad class of nonlinear autoregressive (AR) models with time- varying conditional variance.
Example (Nonlinear AR Model with Time-varying Conditional Variances). Consider the het- croscedastic nonlinear AR model
We ((W, Wr-p) + (oi Wr-g).
where σb IPP IR and σ RR + R are unknown functions, { is sequence of independently and identically distributed (i.i.d.) random variables satisfying yr (cx) = The integers p and denote the AR order, with p S q. This model conforms to the general formulation in (1.1) by setting N q- YN = Wi, and X, = We- 11 for i [q]. Under this representation, the conditional quantile function fo(X1, Xx) ) is given by ø (W,,I...., (..p.)...,.-p))..
Contributions. We summarize the main contributions of this work as follows.
To the best of our knowledge, this work provides the first approximation error bounds for functions within a hierarchical interaction model [Kohler and Langer, 2021] using both RNNs and SRNNs. These bounds highlight the ability of such architectures to effectively capture the complexity inherent in hierarchical interaction structures. Our analysis builds upon the close connections between SFNNs and SRNNs, as established in Lemma 4 and Lemma 14 of the supplementary material. particular, we show that any SRNN can be represented by an SFNN with a slightly larger but less sparse architecture, and vice versa, indicating that their respective function classes are comparable expressive power. This result complements the equivalence between FNNs and RNNs previously demonstrated by Jiao et al. [2024].

---

Built upon the established approximation error bounds, we conduct comprehensive error analysis for nonparametric QR with weakly dependent data using RNNs and SRNNs. Specif- ically, we estimate the true conditional quantile function fo within hierarchical interaction model characterized by intrinsic smoothness 7² We show that, for a stationary exponentially B-mixing sequence of observations, the empirical risk minimizers based on both RNNs and SRNNs achieve convergence rate of n 22 /(27*+11 under the squared L2 norm, up to logarithmic factor. This rate coincides with the minimax-optimal rate established by Schmidt- Hieber [2020]. Furthermore, we derive slower convergence rate for the case of algebraically --mixing dependence.
Notations. For two sequences {an and {bn], we write On Z bn if there exists constant C' > independent of n, such that an 2 Cb,, We write X by if bn Z On Additionally, we use the notation On bn when both an S bn and on x bn hold. For any E R, let |a denote the largest integer that is strictly smaller than a, and [a] denote the smallest integer that is strictly larger than a. We denote N {1.2.... .,, No NU {0], and İ f....-2 -1, 1,2..... For any N E N, we define [N] as (1, .I.....N..
NJ. For any set S, we denote its cardinality by [S]. For x (I), (I1,..., v.) E Rd we define ||x||e = ( x. p)//// ||x||oo max; |za]]. ||x||o 1 * 0). For A [ar,] E mm we define ||Al|o = E, 1 * 0) and vec(A) (aI.I...., al........ am. am.....mmmm)). am.n) E IR Moreover, for any real-valued function h defined on domain X, we define |//||b = supxex [1(x)]. Let Pxx be probability measure on X. The Lp norm of h s < 2o) with respect to Px is defined as ||bl|p (Ex-px (P))))
2 Methodologies
In this section, we formally define the architectures RNNs and SRNNs, both employing the ReLU activation function, (zx) max(2, these architectures, we then introduce nonparametric QR estimators based on RNNs and SRNNs, constructed through empirical risk minimization (ERM) over overlapping subsequences of the data.
[z, 0}, applied elementwise vector inputs. Building
2.1 RNNs and SRNNs
An RNN is characterized by the following parameters: the input dimension dx the output dimension dy. the width W, and the depth L. Given a time horizon N, an RNN processes an input sequence X = (x(1). (NN) R sequentially through three types of layers: input layer sequence of recurrent layers {ri)/e and an output layer q. The architecture generates an output sequence : (1) , E Rd according to
Vo -p( (Px Px(N))
Vi= m (V-1) = (r/ (V-1), 1 (Vi-))), lE [L],
==q(Vi.) (Qv()) ......QQ... (N)). (.....Qv[.))
where P E Wxxx Q Rdyxx xW and V12 = (v(")... ......y?
XN for E [L]U {O). For each time step E [N] and each layer E [Z], the recurrent operation 7((9) is defined by
vi (N) R
(e) (V-1) == σ (Am) (t" (Vi-1) Biv +c),
where A,, Bi RRW., C W and r (0) By composing all L + 2 layers, the overall RNN function ry can be expressed as
0 O.R.
To = -oror. ........ ori op.
where (vec(P) vec(AA) vee BB) T. C with de = (2W2 + W + W(dx + dy). In the many- one RNN setting considered here, the prediction function corresponds to the final element of the output sequence, denoted by r (N)
TT vec(Bz) c7,, vec(Q) e Rd
We define RNN as the class of RNN prediction functions with bounded outputs:
RNN dx dy (W, K)
(N) r8
R
A + B3
sup X XERdxxx ME[N
1| (X) S K
3

---

For simplicity, we assume K 2 1 throughout and omit it when boundedness is either understood or required. Likewise, when the input and output dimensions dx and dy are clear from context, we denote the class more compactly (.N.NWW L.).
Such RNN architectures, particularly those with large width depth, often substantially parameterized in practice, which can lead to severe overfitting during training, Introducing sparsity provides an effective mechanism to mitigate this overparameterizztion.. From a theoretical per- spective, sparsity reduces the effective hypothesis space, thereby improving generalization bounds. We now formally define RNN sparsity. For recurrent layer rı, define its sparsity as ((r)) ||Aillo + |B + ||c/ For the input and output layers, we set T(p) = Pllu and T(g) [Q] respectively. The total sparsity of an RNN rg is then given by T(ro) T(p) (()) Ei= T(r). Using this notation, we define the class of SRNNs as
SRAN de dW,
(N) E RNN d2 (W,,,, with T(re) S ss). (1om) 10,-0 (N)
where denotes the sparsity budget. SRNNS combine the sequential representation power of RNNs with the computational and statistical advantages of sparse architectures, making them partieu- larly attractive for both theoretical analysis and practical implementation in resource-constrained environments.
Nonparametric quantile regression
As is standard in the QR literature, we begin by imposing regularity conditions on the noise variable e in (1.1), conditioned on the covariates XI,...; N, as formalized the following assumption.
Assumption 1. The conditional density of given X1, Xi...., X denoted by Pe|x, Xxx exists and is continuous over its support. Moreover, it satisfies, almost surely over XI...., XN.
p S Pe|Xi .X.......... Xx (0) S sup P (u) p
ueR
for constants p 2 p >>
O. In addition, there exists constant lo such that, almost surely over X1, X XN. PP.j, (u1) - Pelxi (u2) S 40/u1 u2l for all un ug E I.
This assumption is standard in the literature and has been adopted in prior works such Belloni and Chernozhukov [2011], Belloni et al. [2019], and Padilla et al. [2022]. plays a crucial role in linking the excess risk to the squared L2 error of the estimators, established in Theorem 3.
Under Assumption 1, the target function fo is the unique minimizer of the population check loss
R+( E(x,..... Xx)- II, YN) (Yx
(xX1,..,xNN)..
where II denotes the joint distribution of X1, .I,..., loss. Hereafter, assume that the joint distribution II has compact support [0, 1|dxxN, Given dataset {(x,, 8/1 M we form overlapping subsequences to construct the training sample
.,. XN and (+( (T 1(u < 0))u is the check
S {((x)..... xN), yN). ((xg.... XN+1), (1++1),n+1),.... yN+1),.... ((xn-N N+1.-NN+......... x.), vn)).
We then estimate the target function fo by ERM over function class F, yiclding the estimator
f E arg min Ra(s)) = N Σ
-
(2.1)
fEF
t=
In the sections that follow, we derive error bounds for the cases where F RNN and J SRNN 4 (W, L, K, s), respectively.
3 Statistical Theory
In this section, we begin by introducing the hierarchical interaction model and derive error bounds for function approximation under this framework using both RNNs and SRNNs. Building on these results, we first establish oracle-type inequalities and subsequently use them to obtain separate upper bounds on the To2 error of QR estimators based on RNNs and SRNNs, under the assumption that the data-generating process satisfies stationary A-mixing condition.

---

3.1 Approximation error bounds
The theoretical performance of neural networks critically depends on the properties of the underlying function class. A commonly adopted assumption in nonparametric statistics is that the true regression function belongs to a Hölder class. We recall its definition below, following [Stone, 1982].
Definition 1 (Hölder Class of Functions C. (X, K)). Given domain X' C IR positive Hölder
smoothness parameter 83 and a constant K > 0, the 3-Hölder function class is defined as
000h [(x) 8 [(y)
C. (x)
f - R
Σ
a f Allc0
+
Σ
sup X.E'X'
K
X y 12
ailaalI <B
ac
x-y
where [3]. 3 -r, 0a = .g....00 gad with (ai,...,, aa) E N6 and ||a||| Ei= az.
Moreover, we refer to 7 = 3/d as the addension-adjuused degree of smoothness of C. (X, K).
Specifically, any function f C (X, CX, K) is bounded in magnitude by K. Without loss of generality, we assume throughout the paper that K 2 1. Stone [1982] established that the minimax convergence rate for estimating a regression function under the L2 norm over the Hölder class C² (X' K) is 1) However, in neural network applications, the input dimension d is often large, resulting in a small value of the dimension-adjusted smoothness y, which in turn leads to slow convergence rates. To mitigate the impact of high dimensionality, we consider a class of functions with compositional structure, known as the hierarchical interaction model, which captures intrinsic low-dimensional structures and allows for more favorable approximation and estimation properties.
Definition 2 (Hierarchical Interaction Model). Let I, d E N be positive integers, and let K > 1. Suppose P S [1, oc) x N is a parameter set such that ()P().)) max (3, [) < co. The hierarchical interaction model (P, K) is defined recursively as follows:
(i A function h IRd + IR belongs to the first-level model H,(P,K) K) if there exist (B,t) E P. a function h1 C ( K), and a set of indices
(ji....,je) S {1....,d] d} such that h(x)
by (1)) -.pp.-... for (5is-... 3r) E Rdd
(ii) For / >> 1, a function h belongs to the hierarchical interaction model HP, K) if there exist (B,t) E P, a function h4 C₂2 (IR K), and functions u1,..., Ue H1'- (P, K) such that h(x) hı(u1(x) ..., ue (x)) for Rdd
In analogy to the Hölder class, we define the intrinsic smoothness of H, K) by >" B where (3* t*) arg min
This quantity can be interpreted as the effective smoothness of the least regular component in the composition. Crucially, it does not depend on the ambient input dimension, thereby mitigating the curse of dimensionality, Extensive research has established the minimax-optimal convergence rates for models related to the hierarchical interaction model in nonparametric regression, demonstrating that these models retain favorable statistical performance even in high-dimensional settings [Bauer 2020, Kohler and Langer, 2021]. Moreover, the hierarchical
and Kohler, 2019, Schmidt- interaction model encompasses a broad class of structured functions, including classical models such additive models, single-index models, and other compositionally structured function classes [Kohler and Langer, 2021].
Based the preceding definitions, impose the following assumption on the true quantile function fo introduced in (I.1). Assumption 2. Let K 2 1, N, and P C [1, oo) x N satisfy (.P((.)))P max < cc. The true quantile function fo belongs to the hierarchical interaction model Hi x (P,K)
The following theorem provides an error bound for approximating functions within hierarchical interaction model using RNNs.
Theorem 1. Under Assumption 2, for any Wo. Lu 2 3. and probability measure μ on [0. 1]dxx that is absolutely continuous with respect to the Lebesgue measure, the following inequality holds
{ 1f(x) - fo(x))2((ds))
1/2
inf S.R.NN,..WW...E
< (WoL))
5

---

where
W = (d + W) log Wol + 1 and L= 20g Lo log Lol + 2N.
(3.1)
Here, the positive constants c1-03 depend on (1, (K,
Theorem establishes the ability of RNNs approximate hierarchical functions. In contrast, while Jiao et al. [2024] also explored the approximation capabilities of RNNs, their analysis relied on the restrictive assumption that the true regression function belongs Hölder class. As a result, their approximation error bounds exhibit a strong dependence on the input dimension. By leveraging the benefits of the bierarchical interaction model high dimensions, a more thorough exploration of the theoretical properties of SRNNs becomes possible. This is particularly relevant as SRNNs are commonly used in nigh-diimnssinal data scenarios. We then present the following theorem, which provides a result for SRNNs that is not achievable by Jiao et al. [2024].
Theorem 2. Under Assumption 2, for any Wo 2 (uP(S.(jeP max (3 + 1),, (K + 1)c"), Lg 2 1, and probability measure μ [0, 1]d,>* that is absolutely continuous with respect the Lebesgue measure, the following inequality holds
1/2
Sarnaaa 1s(x) - [0())
Sea Wo22 Lo + Wo'')
inf JESRNN az KW,,
where
W (dx Wo +1, T= = Tog 2N, e7 Lo W + 3(dx + 1)es NoN++6(coLa + N) + 1)csWo
(3.2)
Here, the positive constants Cτ-c7 depend on (,,
The approximation error bound established Theorem demonstrates the effectiveness SRNNs in approximating hierarchical functions. The key distinction between the error bounds in Theorem 1 and Theorem lies the construction of the neural networks. Specifically, Theorem 1 based on RNNs constructed using the methodology in Theorem 3.3 of Jiao et al. [2023], while Theorem 2 employs SRNNs constructed through local Taylor approximations, following the approach Yarotsky [2017]. Importantly, the difference in the resulting error bounds becomes negligible when applied to the derivation of convergence rates for QR estimators. As demonstrated in Theorem and Theorem both approximation bounds can be incorporated into oracle inequalities, leading to optimal rates for estimators based on RNNs and SRNNs, respectively.
3.2 Error bounds for QR estimators
The statistical performance of neural network estimators in regression tasks critically depends on the distribution the observed data. Classical statistical theory typically assumes that observations are i..d. In contrast, consider more general setting where the observations {(x,, yt) form stationary A-mixing sequence. To proceed, we introduce the following definitions.
Definition 3 (Stationarity). A sequence of random vectors {2 is said to be stationary if, for any given E Z and k E No. the distribution of the random matrix (Zz, Zeth..... Ze+m) is identical to that of (ZI+k. Zu+k+1..... Zt+m+k).
Definition (B-mixing [Bradley, 1983]). Let [2 be sequence of random vectors. For any i.j E oo,, +oo), define of = (z,,, Zitl...., Zi) as the o-algebra generated by Zı, i S k S j. For any a N, the B-mixing coefficient of the stochastic process o defined as
((a) = sup A
sup ((AI))P(AA)) AEATE
We say that {z, 1,0 is algebraically B-mixing if there exist positive constants By and r > 1 such that ((aa < So/a" a" for all Similarly, it is said to exponentially B3mixxng there exist positive constants 3g 31 and such that ((a) < 3o exp(-3,a") for all a.
The concept of --mixing has been extensively studied in the literature [Yu, 1994, Mohri and Rostamizadeh, 2010, Phandoidaen and Richter, 2020, Jiao et al., 2024]. In particular, number of

---

procedures have been developed to estimate the mixing rate [McDonald et al., 2011, 2015]. We emphasize that our data assumption includes the nonlinear AR model with time-varying conditional variances, which was not considered in Jiao et al. [2024] due to their differing stationarity assumption on the sequential observations.
The following theorem presents oracle-type inequalities for the QR estimator when the function class specified as RNNs and SRNNs.
Theorem 3. Assume Assumption holds, that {(x,, yt) is stationary G-mixing sequence, f x0 I K for some K 2 1. For any positive integer l such that n 2 4CN, we define Mmax 2m N + (Ne) and Mmin = (n N + 1)/(2N)).
and
(i) Let F RNN dx, (W, L, K) with W, L 2 3.
8
inf fl foll2 and 86 W (og L}) log (Mmax) /mmin- feF
Then, there exists constant og > 0, such that for any u 2 1, the ERM estimator fin (2.1) satisfies P (If - foll2 2 es (8e- + 6ú + v u/ lu/mmaxx) ) S Nie " + n (2k_)ogg22//6)]](((-1)) +1). (3.3)
(ii) Let F SRNN de, (W, L, K, s) with W, L 3,
5, = fer
inf || - foll2 and 5. = v VL. log( (.?) log(m
M
Then, there exists a constant Cg > O, such that for any u 2 1, the ERM estimator f in (2.1) satisfies P (IIf folla 2 co (8 + 6g + /ul mmax)) - Nie " n (1/6 + 1). (3.4)
Remark 1. To the best of our knowledge. the established oracle inequalities for QR under dependence represent a new contribution and are technically nontrivial. To place these contributions in context, we begin with a concise overview of the analytical framework, highlighting the key methodological innovations that arise at each stage of the analysis.
Novel decomposition. First, we introduce a donut-shaped decomposition that has not appeared in the context of nonparametric regression with neural networks. Specifically, we introduce donut- shaped sets, which allow us to decompose the probability bound IP(II] foll2 5.) into manageable components by bounding each term separately, where 6, = cs (6. + 6, + ww/mmmmm. for RNNs and 6. = co (8a + 8, + //u/mmax) for SRNNs. This decomposition differs from the direct argument used in Eq. (16) of Jiao et al. [2024], whose analysis of the excess risk critically relies the squared loss. To handle the check loss function, we develop the novel decomposition described above.
Refined blocking technique. Next, we introduce a refined blocking technique to relate the mixing sequence to its i.i.d. counterparts, which differs from the approach in Jiao et al. [2024]. In Step 2 of their proof of Theorem 13, part of the data within each partition is discarded to facilitate analysis under the squared loss. In comparison, our new partitioning procedure, i.e., Figure 2, retains all observations, thereby enabling more comprehensive analysis under dependent data. Moreover, while Lemma 16 in Jiao et al. [2024] plays central role their argument, it cannot be directly applied to our setting. Instead, we employ probabilistic counterpart, i.c., Lemma 5, to carry out analysis.
Sharper inequality. Finally, we develop a novel and sharp empirical process inequality. i.e., Lemma 7, that underpins the tightness of our oracle inequality. One reason the result of Shen et al. [2021] lacks tightness is that their analysis relies on non-sharp application of the Bernstein inequality. To address this limitation, we derive Lemma 7 by building on Theorem 7.3 of Bousquet [2003] and Corollary 5.1 of Chernozhukov et al. [2014].
For both function classes, the non-asymptotic bound in Theorem comprises three components: the approximation error éu. the stochastic error 8g and a dependence-adjustmee for the discrepancy between dependent and independent sequences. The integer serves as a key parameter, commonly introduced in time series analysis [Nobel and Dembo, 1993, Yu, 1994], to bridge the behavior of mixing sequences and their i.i.d counterparts. By appropriately tuning the network parameters and selecting to balance the trade-off's among these components, the ensuing theorems establish the convergence rates of the corresponding estimators when the target function So exhibits a hierarchical interaction structure.
term that accounts
7

---

Theorem 4. Let RNN (W,, K) be the hypothesis class F and assume that the probability measure 1I on [0, 1|dxXN is absolutely continuous with respect to the Lebesgue measure.
(i) Suppose Assumption 1 and Assumption 2 hold and {(x,, ye) is stationary exponentially A-mixing sequence. Let Wo, Lo 2 3 satisfy Lo = (log m) (G+1/1)) /244*++2 and define W and L according to (3.1). Then, there exists a constant C1o > such that the ERM estimator f in (2.1) satisfies
((oggn)oou +1/(2+++1)
u u(log /)
(logn)
log
P
|f foll2 2 ciu
n
ea
An immediate consequence is that Ilf foll2 = Op (n~> /22"" +1) (log n) (0+1/r/> 100--11////"""0/0/ +12*+11).
(ii) Suppose Assumption I and Assumption 2 hold and {(x, V/ | is a stationary algebraically B B-mixing sequence. Let Wo, Lo 2 3 satisfy WoLo x (n (n(1 (logn)") 1/(4"++3) +2) and define W
and L according to (3.1). Then, there exists a constant C11 such that the ERM estimator f in (2.1) satisfies
( (1--////)
-12(2++11
u
n //r log
P
folla 2
2011 ci
(log (1-1))
11-1//
e
An immediate consequence is that I1f- follz Op (n- (1-(111/0)) (1-1///1//22""+4) +1) (log (7 (277/(22++))))
Theorem shows that the RNN-based QR estimator f. with suitably chosen network architecture, achieves the minimax-optimal convergence rate n n--2² /(2-" +1) under the squared T-2 norm in stationary exponentially G-mixing setting, up to logarithmic factor. Under stationary algebraically G-mixing conditions, the convergence rate under squared L2 norm slows to n (1-1//r)2 /(2y" +1), Nevertheless, as r - oo, the exponent -(1 1/r)2>* /(2>" + 1) approaches -2 /22" + 1) for fixed y", indicating that the convergence rate becomes arbitrarily close to the optimal rate. Our results differ from prior QR studies in two key aspects. First, unlike the FNNs used in previous works [Shen et al., 2021, 2024], we employ RNNs as the approximating function class. Second, rather than assuming independent covariates as in Sangnier al. [2016] and Padilla et al. [2022], we consider a more general dependence structure, where the sequence {(x, 3t) 1 is stationary and G-mixing. Furthermore, compared to the convergence rate of RNNs in nonparametric regression [Jiao et al., 2024], our results capture intrinsic low- low-dimensional structures of the target function, thereby achieving optimal performance while cireumventing the curse of dimensionality.
The following theorem establishes the convergence rate of the SRNN-based estimator. Theorem 5. Let SRWWW (W, L., K, s) be the hypothesis class F and assume that the probability measure II on [0, 1] 1]2x*N is absolutely continuous with respect to the Lebesgue measure.
(i) Suppose Assumption and Assumption 2 hold and {(x,, w/ ) t=1 forms stationary exponentially --mixing sequence. Let Wo 2 sup(3.0)P max [(3 + 1)r, (K + 1)e') and Lo 2 3 satisfy Wo x
(n/ (log m) exists a constant 612 > such that the ERM estimator f in (2.1) satisfies
and Lo x logn. Define W. L, and according to (3.2). Then, there
[ī Solla 2 c12 (logn)(442//)) x² 72x(2""+>) +1)
u u(llonn) 2/
(logn) log
P
n
n
An immediate consequence is that - Op (n 00
/(2"" +1) (log n) ((1+1/1)""/(1) +1)).
(ii) Suppose Assumption and Assumption 2 hold and ((x,, yo ) forms stationary algebraically B a-mxiing sequence. Let Wo 2 suP(3.0)P max + 1),, (K + I)e'} and Lo 2 3 satisfy Wo x
(n (1(1-117) /log n) )))))))) 1/(22"++1) and Lp x logn. Define W, L. and s according to (3.2). Then, there exists constant C13 0 such that the ERM estimator f in (2.1) satisfies
((
(logn) n) 11-1//-
+1)
uu
1/" logn
1 folla 2 C13
+ \
+ (logn)
1-1
ea
An immediate consequence is that Iff foll2 Op (n- (1-11//)) /00""//22220000000000 (22"+1) ()0 (22"+(22+++1
8

---

Similar to Theorem 4, under stationary exponentially B-mixing conditions, f attains the minimax- optimal convergence rate, up to logarithmic factor. Under algebraic B-mixing, the convergence rate is slower but approaches the optimal rate - oc for fixed -- Compared to Theorem 4, Theorem 5 imposes more restrictions on the choice of width W and length I. of SRNNs, owing to technical reasons. However, when fixing the length I. to be the same for both theorems, Theorem provides guarantees for wider SRNNs. This aligns with the practical setting better because sparsity is commonly considered for deep and wide neural networks.
Numerical Study
In this section, we conduct numerical experiments to evaluate the finite-sample performance of RNN- and SRNN-based QR estimators in comparison to quantile random forest (QRF) and FNN-based Python. The QRF estimator is trained using the
estimators. All experiments are implemented scikit-garden package, and the number of trees is set 100. For all NN-based estimators, we employ early stopping to mitigate overfitting, as proposed by [Raskutti et al., 2014]. Specifically, the dataset is split into training (80%) and validation (20%) sets. At the end of each epoch, the model's s validation performance is evaluated using the mean check loss, and training is terminated if the validation loss does not improve for 20 consecutive epochs. For the FNN-based estimator, we adopt L = hidden layers with widths W1 = W 200, following the architecture in Padilla et al. [2022]. For the RNN- and SRNN-based estimators,
L = 3 layers with hidden width W = 100 and ReLU activation. To induce sparsity. we prune the smallest 40% of parameters and finetune the remaining weights post-training. The hyperparameters above are fixed throughout the main experiments, and we further investigate the impact of varying these parameters in Appendix D. All neural networks implemented in PyTorch.
We conduct experiments on the following two nonlinear AR models at quantile levels and 0.5 to validate the theoretical results established in our analysis.
Model I(SIM:-1):
= - - -00 6
-.6))-
(4.1)
5 2 = Σ (-1)) (Y, + --22++1))
where (-) is the standard normal distribution function.
Model 2:
Y; = -..hh (Y.-1,...,, + h2 Y1-10)
(4.2)
where
(Y cos(22
It (K-- Yo-4 + Y
Y
+ Y te 7h
+ cos( cos( cos(TYy-10),
h2(Yi- . sin (12=-++1-0
In In [1+ 11
1+ L+ Y23
sin sin
In the experiments, we examine both light-tailed and heavy-tailed noise distributions: the standard normal distribution (0, 1) and a scaled Student's !-distribution with 2.25 degrees freedom 22.25- respectively. Theorem in Chen and Chen [2000] ensures that the data generated from (4.1) 1) and (4.2) are strictly stationary and exponentially B-mixing. For each model, we generate 110,000 data points, discard the first 100 burn-in period, and split the remaining data by assigning the last 100,000 points to the test set and the rest to the training set. For each trained estimator, we compute the mean squared error (MSE) experiment. As shown in Figure 1, the RNN- and SRNN- estimators exhibit significantly lower variance and superior overall performance compared to the QRF- and FNN-based estimators.
the test set. We perform 500 Monte Carlo repetitions for each

---

0.6
Model 1
0.2-
Model
Model
2.0
Model 1
0.4
1.0-
0.2
0.0
ORF
FNN
RNN SRNN
ORF
FNN
RNN SRNN
ORF
FNN
RNN
SRNN
ORF
FNN
RNN
SRNN
(8)
[0
1),
[00 1).r
(c) 22.25
1g2 20.0
2.2
2.0
0.15
1.8-
0.30
1.4
0.05
(e)
[0,, 1). 0.1
( NOO, 1). 0.3
(g) 12.25. 0.1
(b) 12:25. 0.5
Figure I: MSE comparison Model Model with different noise distributions and quantile levels.
Application
In this section, we conduct experiments on the DJIA dataset' to empirically validate our approach under the stationarity assumption. To further assess robustness in nonstationary settings, we provide a complementary case study on GDP forecasting in Appendix C. Since log-returns of stock prices widely regarded approximately stationary, the DJTA dataset is well-suited for evaluation in the stationary environment. We partition the data chronologically. allocating the first 19 years for training and the final year for evaluation. We assess the out-of-sample predictive performance of RNN- SRNN-, FNN-, and QRF-based estimators in a bass-nnss-sdy--ahall forecasting task, with implementations following the specifications in Section 4. All models trained by minimizing the empirical check loss and evaluated on the test set using the same criterion.
Table presents the empirical check loss five quantile levels. The results show that the RNN and SRNN estimators consistently outperform FNN and QRF methods at all quantile levels. Moreover, the SRNNs achieve predictive accuracy comparable to RNNs while inducing sparsity. These empirical findings are consistent with, and provide validation for, our theoretical results.
Model
7 = 0.1
T= == 0.25
== 0.5
==
0.75
T == 0.9
QRF
0.456
0.698
0.817
0.616
0.365
FNN
0.538
0.735
0.810
0.662
0.404
RNN
0.410
0.640
0.760
0.562
0.306
SRNN
0.406
0.647
0.759
0.561
0.305
Table 1: Out-of-sample prediction errors at different quantiles for DJIA growth analysis.
Conclusion
This study investigates the convergence properties nonparametric quantile regression using RNNs and SRNNs. Error bounds are derived for the approximation of functions within hierarchical interac- tion model using RNNs and SRNNs, respectively. Based on these error bounds, we demonstrate that, for stationary, exponentially B-mixing sequence of observations, the empirical risk minimizers of both RNN- and SRNN-based methods achieve the optimal convergence rate.
Future research could explore the approximation capabilities of RNNs and SRNNs for nonparametric quantile process regression. Instead of focusing on fixed quantile level т, useful extension would involve approximating the entire quantile process indexed by E [TL, Tu] S (0,1), 1), which complements the results of Shen al. [2024] using FNNs.
We use DJLA data from Jan 1, 2000. Dec 31, 2020, obtained from https: http://twwwwww..tt.////
10

---

Acknowledgment
Zhao Ren was supported in part by the National Science Foundation (DMS-21 1A4S-2113566 and the National Institutes of Health (NIGMS ROI R01GM157600) Wen-Xin Zhou was supported in part by the National Science Foundation (DMS-2401268) and the Australian Research Council Discovery Project Grant (DP230100147).
References
T. Adrian, N. Boyarchenko, and D. Giannone. Vulnerable growth. American Economic Review, 109 (4):1263-1289, 2019.
M. Anthony and P. L. Bartlett. Neural Network Learning: Theoretical Foundations. Cambridge University Press, 2009.
P. L. Bartlett, N. Harvey, C. Liaw, and A. Mehrabian. Nearly-tight VC-dimension and pseudodimen- sion bounds for piecewise linear neural networks. Journal of Machine Learning Research, 20(63); 1-17, 2019.
B. Bauer and M. Kohler. deep learning as a remedy for the curse of dimensionality in nonpara- metric regression. The Annals of Stafistics, 2264)22261- 2019.
A. Belloni and V. Chernozhukov. & 11-penalized quantile regression high-dimensional sparse models. The Annals of Statistics, 39(1))82-130, 2011.
A. Belloni, V. Chernozhukov, D. Chetverikov, and Fernández-Ve Val. Conditional quantile processes based on series or many regressors. Journal of Econometrics, 213(1):4-29, 2019.
O. Bousquet. Concentration inequalities for sub-additive functions using the entropy method. In Stochastic Inequalities and Applications, pages 213-247, 2003.
R. C. Bradley. Absolute regularity and functions of Markov chains. Stochastic Processes and Their Applications, 100in,441111 14(1):67-77, 1983.
M. Chen and G. Chen. Geometric ergodicity of nonlinear autoregressive models with changing conditional variances. Canadian Journal of Statistics, 28(3):605-614, 2000.
M. Chen, X. Li, and T. Zhao. On generalization bounds of family of recurrent neural networks. In Proceedings of the 23rd International Conference on Artificial Intelligence and Statistics (AISTATS), pages 1233-1243, 2020.
X. Cheng, Huang, and S. Ma. Generalization and risk bounds for recurrent neural networks. Neurocomputing, 316:128825, 2025.
V. Cheroozuukvv, D. Chetverikov, and K. Kato. Gaussian approximation of suprema of empirical processes. The Annals of Statistics, 12644:156411597 2014.
Chung. C. Gulcehre, K. Cho, and Bengio. Empirical evaluation of gated recurrent neural networks on sequence modeling. arXiv preprint, arXiv:1412.3555, 2014.
J. Fan, Y. Gu, and W.-X. Zhou. How do noise tails impact deep ReLU networks? The Annals of Statistics, 18455)88858771 2024.
X. Feng, X. He, Y. Jiao, L. Kang, and C. Wang. Deep nonparametric quantile regression under covariate shift. Journal of Machine Learning Research, 25(385):1-50, 2024.
S. Hochreiter and J. Schmidhuber. Long short-term memory. Neural Computation, 98): 1997.
Y. Jiao, G. Shen, Y. Lin, and J. Huang. Deep nonparametric regression on approximate manifolds: Nonasymptotic error bounds with polynomial prefactors. The Annals of Statistics, 51(2):691-716. 2023.

---

Y. Jiao, Y. Wang, and B. Yan. Approximation bounds for recurrent neural networks with application to regression. arXiv preprint, arXivv2409.05577, 2024.
R. Koenker and G. Bassett. Regression quantiles. Econometrica, 46(1):33-50, 1978.
R. Koenker and Z. Xiao. Quantile autoregression. Journal of the American Statistical Association, 101(475):980 990, 2006.
M. Kohler and S. Langer. On the rate of convergence of fully connected deep neural network regression estimates. The Annals of Statistics, 49(4):2231-2249) 2021.
D. J. McDonald, C. R. Shalizi, and M. Schervish. Estimating beta-mixing coefficients. In Proceedings of the 14th International Conference on Artificial Intelligence and Statistics (AISTATS), pages 516-524, 2011.
D. J. McDonald, C. R. Shalizi, and M. Schervish. Estimating beta-mixing coefficients via histograms. Electronic Journal of Statistics, 9(2):2855 2883, 2015.
D. J. McDonald, C. R. Shalizi, and M. Schervish. Nonparametric risk bounds for time-series forecasting. Journal of Machine Learning Research, 18(32):1-40, 2017.
N. Meinshausen. Quantile regression forests. Journal of Machine Learning Research, 7(6):983-999, 2006.
M. Mohri and A. Rostamizadeh. Stability bounds for stationary e-mixing and f--mxiig processes. Journal of Machine Learning Research, 11(2):789-8 814, 2010.
V. Nair and G. E. Hinton. Rectified linear units improve restricted Boltzmann machines. In Proceedings of the 27th International Conference on Machine Learning (ICML), pages 807-814, 2010.
A. Nobel and A. Dembo. A note uniform laws of averages for dependent processes. Statistics & Probability Letters, 17(3):169- 172, 1993.
O. H. M. Padilla, W. Tansey, and Y. Chen. Quantile regression with ReLU networks: Estimators and minimax rates. Journal of Machine Learning Research, 23(247):1-42, 2022.
N. Phandoidaen and S. Richter. Forecasting time series with encoder-decoder neural networks. arXiv preprint, arXiv:2009.0885 08844, 2020.
G. Raskutti, M. Wainwright, and B. Yu. Early stopping and non-parametric regression: An optimal data-dependent stopping rule. Journal of Machine Learning Research, 15(1):335-366 2014.
D. E. Rumelhart, G. E. Hinton, and R. J. Williams. Learning representations by back-propagating errors. Nature, 323(6088):533-536 1986.
M. Sangnier, O. Fercoq, and F. d' Alché Buc. Joint quantile regression in vector-valued RKHSs. In Advances in Neural Information Processing Systems 29 (NIPS), pages 3693-3701, 2016.
Nonparametric regression using deep neural networks with ReLU activation J. function. The Annals of Statistics, 48(4):1875-189, 2020.
G. Shen, Y. Jiao, Y. Lin, L. Horowitz, and J. Huang. Deep quantile regression: Mitigating the curse of dimensionality through composition. arXiv preprint, arXiv:2107.04907 2021.
G. Shen, Y. Jiao, Lin, J. L.. Horowitz, and J. Huang. Nonparametric estimation of non-crossing quantile regression process with deep ReQU neural networks. Journal of Machine Learning Research, 258 2024.
G. Shen, R. Dai, G. Wu, S. Luo, C. Shi, and H. Zhu. Deep distributional learning with non-crossing quantile network. arXiv preprint, arXiv.2504.08215, 2025.
C. Song, G. Hwang, J. H. Lee, and M. Kang. Minimal width for universal property deep RNN. Journal of Machine Learning Research, 24(121):1-41,
2023, 2023,
12

---

C. J. Stone. Optimal global rates of convergence for nonparametric regression. The Annals of Statistics, 1044)1040-10533, 1982.
I. Sutskever, O. Vinyals, and Q V. Le. Sequence to sequence learning with neural networks. In Advances in Neural Information Processing Systems 27 (NIPS), pages 3104-3112, 2014.
Z. Tu, F. He, and D. Tao. Understanding generalization in recurrent neural networks. In Proceedings of the 8th International Conference on Learning Representations (ICLR), 2020.
M. Vidyasagar. Learning and generalisation: with applications to neural networks. Springer Science
& Business Media, 2013.
H. White. Nonparametric estimation of conditional quantiles using neural networks. In Computing Science and Statistics, pages 190-199, 1992.
D. Yarotsky. Error bounds for approximations with deep ReLU networks. Neural Networks, 94: 103-114, 2017.
B. Yu. Rates of convergence for empirical processes of stationary mixing sequences. The Annals of Probability, 22(1):94 116, 1994.
S. Zheng. Qboost: Predicting quantiles with boosting for regression and binary classification. Expert Systems with Applications, 39(2):1687-169) 2012.
13

---

Appendix
The appendix is organized as follows: Appendix A provides the proofs of the main theorems presented in Section 3. Appendix B contains the proofs of all auxiliary lemmas referenced in Appendix A. Appendix presents an additional case study on GDP forecasting that further supports our theoretical findings. Appendix D provides a sensitivity analysis of the hyperparameters for the SRNNs.
A Proofs of Theorems
We first provide the definition of FNNs for our analysis. An FNN fo :R - IR can be represented as
fo (x) := fL ....o fi(x),
x XE RER"
where each feedforward layer fi(x) defined as
[A(Ai bi), fi(x) = AAx + bi,
1E L-1-11 =
For each layer l, the weight matrix is A E RWixWi-1 and the bias vector is E WRW, The initial di- mension Wo dx and the final dimension WT = dy correspond to the input and output sizes, respec- tively. The complete parameter vector is defined as 0 = (vec(A1) b, ..... vec(AL) b[]" E REL +1) The width W of the network is defined max (Wi, Wi-13. We denote FNNo (W, L, K) as a class of FNNS with width W, depth I. and bounded output, defined as
FNN dx dy (W, W,L.K, K)-
Jo:Re
Rdy
sup XER
[|o(())l S K sK}
Similar to SRNNs, define T(fi) | + |Ibillo. for l E [L] and define T(fo) - E T(f) as the sparsity of the SFNN fo. Moreover, we define the class of SFNNS as follows
SFNN de dy (W, L. K) = [ffe fo E FNN de (W, L, K) with T(fe) S ss).
To maintain consistent input formatting, if the input to an FNN or SFNN a sequence X E R first stack its columns into vector and then use this vector as the input to the neural network, that is. (X) fe( fe(((x())) ...... ) )T)) For notational simplicity, the distinction between representing the input as a sequence or as its vectorized form is not explicitly made where the context is unambiguous.
For simplicity, we assume K > 1 throughout and omit when boundedness is understood or not required in the notation of the function classes FNN, SFNN, RNN and SRAN.
A.1 Proof of Theorem 1
We first introduce two lemmas: the first establishes an error bound for using FNNs to approximate functions within a hierarchical interaction model, and the second demonstrates that an FNN can be represented by RNN prediction function.
Lemma 1 (Proposition 3.4 in Fan et al. [2024]). Given hierarchical interaction model HP K). for any Wo. Ly 2 3, and probability measure μ on [0, 1]42xN that is absolutely continuous with respect to the Lebesgue measure, the following inequality holds
Vorieann 1f(xx) -f(())
1/2
inf JoEH)((P.K) JEFNN
sup
I <loooo)) 27"
where W cıs Wo log Wol and L C16 Lo log Lo]. Here, the positive constants C14 depend on (1, (.P.
Lemma (Proposition 2 in Jiao et al. [2024]). For any FNN f E FNN de N.d, (W, L), there exists an RNN prediction function f E RNN de,dy ((dx + 1)W + +1,2L + 2N) such that
f(x) f(X), X E [0, 1
14

---

Proof of Theorem 1. Applying Lemma to Xx (P, K), there exist positive constants C14-C16
such that
1/2
JEFNN (W,K) (OLAIAAA - - 50(())
inf
I ccu (cau(WLLo
where W = c15 Wo log Wo], L' = Ciel log Lo]- By Lemma 2, we have that for any f FNN NWW T/, K) there exists an RNN prediction
function f E RRN dR. ((dx + 1)W" 1, 2D' + +2N, K) such that f(X) = f(x). X [0, 1]
Combining the two results above, complete the proof of the theorem.
A.2 Proof of Theorem
Following similar line of reasoning in the proof of Theorem 1, we introduce two key lemmas that extend Lemma and Lemma 2 to the sparse setting. The first lemma establishes an error bound for approximating functions under hierarchical interaction model using SFNNs. The second lemma shows that an SFNN can be equivalently represented by an SRNN prediction function. Proofs of these lemmas deferred to Appendix B.1 and Appendix B.2, respectively.
Lemma 3. Given hierarchical interaction model (P,, K), for any Wo 2 SUP(8.))EP max
1)², (K + 1)e'), La 1, and probability measure μ on [0, 1]d,xN that absolutely continuous with respect to the Lebesgue measure, there exist W. L, > such that the following inequality holds
inf (oE PES.N (OOIEANNN 1f(x)) - (2 ((dx)
1/2
sup
S C17 Wo2 Lo + "
where W = C1s Wis = L and s = Loo Here C17-C20 are positive constants depending on K).
Lemma 4. For any SFNN f E SFNN dg x (W, L, there exists an SRNN prediction function
d ((dx + 1)W + 1,2L + 2N, + 3(dx + 1)WN +6L () + 1)W) such that
f(x) S(X), X [0, 1J<xN
Proof of Theorem 2. Applying Lemma 3 to H/, (P, K), there exist positive constants C17-C20
such that
1/2
TESFNNOXXXIIH KK', ( 1f(x) [6(x)
inf
S Sei7 Wo2 Io + +Wo ")
where W = = L' C19Lo and s' s' c20 Furthermore, Lemma establishes that for any SFNN f E SFNN de (WW (W','', L', K, s'), there exists an SRNN prediction function of E SRA M.NNN + 1)W" +1,2L' + 2N, K, 3(xx 1)W'N + 6(L' + (d + 1)W') such that f(X) = f(X), X e [0,1]xxN
Combining the two results above, we complete the proof of the theorem.
A.3 Proof of Theorem 3
As preparatory step for the proof of Theorem 3, we present the following three lemmas:
Lemma (Lemma 2 in Nobel and Dembo [1993]). Let real-valued Borel measurable function, and let { 1.00 be stationary G-mixing sequence of random vectors. We define block of length
15

---

L as a sequence of L consecutive vectors from {zilt= Consider sequence C consisting of m such blocks, where the gap between any two consecutive blocks is a 1. That is, for any two consecutive blocks in C, the index of the first vector in the second block exceeds the index of the last vector in the first block by exactly a + 1. Then for any constant I, the following bound holds
[P()CC) 2 t) P(h(C) 2 +) S mm + 1).
Here Č comprises m independent blocks, each drawn from the same distribution as in C.
The lemma above enables to extend concentration results for i.id. data the setting of f-mixing data. Similar results can be found in Yu [1994], Vidyasagar [2013], and McDonald et al. [2017].
Lemma 6. Under Assumption 1, for any function f: [0, 1jdxN |- K, KJ. the population check loss function satisfies
2 - & R-(f) R.(f)) S c2all - foll2.
where 02 = min (2) and 022 - 2/2.
In the analysis of Theorem 3, we consider a generie function class F. The main difference between (i) and (ii) lies in the properties og. which are determined by whether F az. (W, L, K) or SR.N 4_ (W, L, &). We present key concentration result for each of these function classes with their corresponding ò6. For any F. we define F(6) = ={ʃ f S 6}. For convenience, we denote 2, ((x,~ N+1, Mt-NN+1) -N11), ... (xt, v() and S =[+) Additionally, for any S. we introduce the following notation:
g(f, 2) Pr (Yn f(X XN)) pr (Yn - fo(X)....,N))).
(A.1)
With the above definition, have the following lemma.
Lemma 7. Let 21...., Zn be i..dd random vectors drawn from the distribution of 7, and let be defined as in (A.1). We choose ós = WL./log(mmx(WW, log n/n for RNNs, and 6 = WI, log(( W 1.2) log n/n for SRNNs. Then there exists universal constant C23 > 0 such that for any 5 > 8, and 0 < x S no², the following incquality holds uniformly over both F RNN d (W, L, K) and SRNN d... LW,,
P sup (EF())
g(S.z,) Ez|o(s. 2)] 2 023d
8b + V
I
SE
Proof of Theorem 3. First, let m n N + 1 and, for given 1, define 8, C(6a + 86 /2Néu/m), where is given by
C max /822//21 22223/c22 c21 1.
(A.2)
Here, the constants 02 and 022 are specified in Lemma and 623 is specified Lemma 7. We emphasize that in Theorem 3, both ca and Co can be set to 4C. Next, for an integer E N, we define the donut-shaped sets as
D: = \ - 2 < 11 - fullz s 2'5.}.
With the definition of D,, write
[log (2K/8.)
P (15-Joll foll2 > 6.)
Σ P (fEED))
(A.3)
Therefore, it reduces to bounding each probability P(f E DD) separately. Following Lemma 6, for any f E D,, we have
C21 22--282 < c2ill - Joll < R2 R+( - R-)(o)..
(A.4)
We next derive an upper bound of the right-hand side of (A.4). By the definition of 6a. there exists fm F such that Il folla _ 25,. Now, if f E D; we have
R, ( R, (fo)
16

---

-R,(f) - (f) + (R. - Rm((m)) + (m (fm)) + (. - R-(fo) SR-(I (π) + R. (Jm) (mm) + R. R.((o))
where the last incquality follows from the definition of f. By Lemma 6, it follows that R, (fm)
R-( v 46226" For any set S, we denote
ss( 1
ZES
Then the earlier inequality be further bounded as
1
R-( - (f6)
+ 4e220
(A.5)
m
Combining (A.4) and (A.5), we obtain an upper bound of the probability IP (f E D:) as
p ((EED) SP 3f E D; such that (fmm As:( 2 22 - 402203 (A.6)
S' = ={ [[z where = UI-N (x,, Mt)))
Partition
B, (by, where by. ZiN-j. 5n jE (0) IN 1].
Partition
B3 (bjjate) mn3-1 a E [6) jE -1].
Figure 2: Illustration of the Partitioning Procedure. Here, m, = |By| and m = | We first split S' into N sequences. Each sequence Bj is then partitioned into equidistant seq-sequences By.
We next relate the dependent observations independent observations. We partition S' into Ne sequences step by step. We first partition it into N equidistant sequences {B, For j E
(0) [N we define B,= {bya)==i. where bj, = Z j E (0) U|N 1] and m, [B)]. Let rı m mod N denote the remainder. Then, for < m, [m/N], and for j 2 Fi. my = [m/N] - 1. We continue partitioning each B, into equidistant sub-sequences { \
For [e]. define B {b (bja+te) m=-- where mß |8% Let mj mod & denote the remainder. Then, for T2, we have my [m,/E], and for r2, we have mf. / 1. The detail is shown in Figure 2. With the partition, we have
N As = Σ
N-1
Σ g(f.b) - m3
Eey (J).
1=
BEB
J= a=1
Plugging it into (A.6) and noting that fm E F(2'8.) for all i (since 26, < 2'8, for all i > 1). we obtain
P (fED,
N-1 8
Agg (fm) - Ag. (f)) 2 22 2212 - 402282
<P Bf E D, such that
m
4
1= a=1
17

---

3f Dj such that max j.0
( (fm) ABy 100(f)) 2 Ne
(21 20 6² -402263)
4
m
C21 22 202208)
max sup Ja (EF(216.)
Agg
(f)
N 8
m C21 2288?
max
sup
A (f)
jia fEJ(2)6.)
N2 16
N-1 8
m C21 2 52
P
sup SEF(2%A.)
A B (f)
(A.7)
Nė 16
where the fourth inequality follows from the choice of C' in (A.2) and the last inequality follows from the union bound.
For each By. we leverage the property of G-mixing to relate it to By. The key here is that the blocks comprising B3 are mutually independent, with each block containing exactly N elements. Furthermore, the elements within each block b By are identically distributed those in the corresponding block b E Bj.
22. By Lemma 5, we have For any j E {0} U IN - 1] and E [e], we denote c NY 16
sup SEF(2)6.)
Ag (f) 2 - IP
sup (EF(225.)
A (1)
4i
S ( + 1). (A.8)
For any j [0] U [N 1] and e [e]. we next bound the probability
m C21 22882
P
sup fEF(2'8.)
A (f)
Nł 16
C21 2 82
<P
sup (2.(2
A Rg (f) 2 m
32
C23 228 82
<P
sup FEF(2'6.)
A B (f) 2 m
C
where the first inequality follows that m/(NE) 2 m /2 and the last inequality follows that c2a/C S c21/32. /32.
To this end, we choose 5 = 2'8, and x = 2 Since C 2 1, we have d 2 5, and 0 S S mz02. Then, Lemma 7 yields
PP
sup fEF(2'6.)
A Ra (f)
223 23 221 82 m
C
Σ g(S.b) Rz|o(,, 7)]
C23 28 C
sup fEFF(2AA.)
m
bEB²
2N Cx
SP
sup
Σ
(({bb) - Ezla(s 2)]
> 0236
80 -
fEF(2'5.)
m
m
bEB
sup fET(256.)
Σ g(S,b) Ezly(S2
2 C236
80 +
m
m
bE
< exp(-x)= = exp( (xp(-222u)) u),
(A.9)
where the second inequality follows that m/ (NE)= 2m..
18

---

Combining (A.9) with (A.3), (A.7) and (A.8) implies
[1. Soll2 > 6x)
(2Kg2(2K/6.]]
SNE
exp(-22ku) + 22/ 1)N + 1)
N
exp( (i))
2m logg 201goo////////)(((111111
<NLe loga (2)))))))//))))(((((((
where the last inequality uses the fact that u 2 1. This proves the claim.
A.4
Proof of Theorem 4
According to Theorem obtain
8a = inf |I. -oll2 2 (WoLn)-27 JEF
(A.10)
Recall that
2N (max [, L}) log 2mm((NN)
5s WWL V
2N l log Wo log Wo, Lo log Lo}) log m
Wo log Wo log Lo
(A.11)
m
We now analyze how these bounds behave under different 3-mixing conditions.
Case Exponentially 8-mixing. Recall that for exponentially 3-mixing, there exist positive constants So. 31, and such that (aa) < Su exp( -3aa) for all a. By plugging (A.10) and (A.11) into (3.3), and é x (log such that (((E - 1)N + 1) x 1/n², together with setting WoLo x (n/(log (n/(lon)((0t)) (0+1)))))) +2), we obtain that there exists a constant c10 0 such that
(Logn))(+1)) 122(2++1)
(log n) /ru
(logn)
log
P | foll2 2 c10
T
n
Case Algebraically B-mixing. Recall that for algebraically B-mixing, there exist positive constants 9o and > such that B(a) Bo/a" for all k. By plugging (A.10) and (A.11) into (3.3), and choosing f x n 1/r log together with setting WoLo x (n(l-1/t/(lognn) 1/(47" +2) we obtain that there exists constant 6n such that
(log n)² 11/r
71/(2* +1)
log 1//r
n 1/r l/gloogn
P 1lf
folla en1
(log n) (1-r)
A.5 Proof of Theorem
According to Theorem 2, we obtain
5, inf || fo foll S Wo2 Lo + W."
(A.12)
FEF
Recall that
2N log MoLLllogWWI.. (WL2) log(2m/(Ne)
2. 2N log W, I.)) log
6, =
I
(A.13)
m
19

---

We now analyze how these bounds behave under different G-mixing conditions.
Case 1: Exponentially A-mixing. By plugging (A.12) and (A.13) into (3.4), and choos- ing l x (logn)4/r such that B((e 1)N + 1) S 1/n², together with setting Wo x (n/ (log (log))(+))) (0)g)))((/)))1113) (/2/))+)+/+++))) Lo x logn, and x WoLo, obtain that there exists constant C1 > such that
(logn)(411/m ((Log)))
(
+1/(2+++
(log n) 12
(log U/
log
P 1 foll2 2 612
ea
"
. Case 2: Algebraically f-mixing, By plugging (A.12) and (A.13) into (3.4), and choosing e x
n 1/1 log n, together with setting Wo x (n2011/)) (logg) 1/12"* +1) we obtain that there exists constant C13
Lo x log n. and x WoLo.
such that >" >2/(2"" +1)
(log (logn)' n
lolg n n 1/r
n 1/r log n
PP [ foll2 2 c1a
(log n) (1-r)
1-/r/
6
B Proofs of Lemmas
B.1 Proof of Lemma 3
This proof is adapted from Schmidt-Hieke [2020]. In preparation for the proof of Lemma 3, we first present two auxiliary lemmas that will be useful the analysis:
Schmidt-Hiele [2020]). For any function fo E (10, Wo 2 max((B + 1)S, (K + } and Lo 2 1, the following inequality holds
Lemma 8 (Theorem
feSFNNa.. (W, |lf - foll & C21 Wo2 Lo + c25 Wo
inf
where W = C26 Wo L = Lo and = C2s LoW Wo. Here, the positive constants C24-C2s depend on (B, d, K).
Lemma 9. Let fi E SFNN d, (W1,L= $1) and f2 E SFNN d (W2, L2, 82). Their aggrega- tion and composition satisfy the following properties
Aggregation rule: If L1 Le and di di, then (2) e SFNN di ddatd', (W1 W2, L1, $1 + 82).
Composition rule: d2 = di, then f = f2 ((o()( E SFNN dy (max (max (W1, W2), Li + L2 1, + 82).
We omit the proof for Lemma 9 since it is straightforward to obtain the result. We refer readers to Section 7 Schmidt- Schmidt-Hicber [2020] for more details.
Proof of Lemma 3. In this lemma, we consider a fixed input domain within [0, 1]d. By the definition of hierarchical interaction model, for any fo E H/ (P, (P,K), there exists a sequence functions {gi){ such that
fo = gi ...00...
where 91 1J + |- K and |- KJ + |-K, Kj+++1 for i 2,..., While this formulation appears slightly different from our original definition of the hierarchical inter- action model, they are indeed equivalent. Specifically, for each level i = 2,...,,, we have h T where each component function hijj E Cf, ([- K, K],, K) for j E [[++]]- Similarly, 91 (In,,,...,,... T with hn E Ci,, ((0, 11a, K) for j [12]. Here, (Bi, E for i [4], t < and t+ 1.
=
To facilitate the application of Lemma 8, the sequence of functions (oi) is transformed into another sequence {h+)1=1 as follows
1
h (x)
2K 91 (x)
Vx E [0, 1
2'
20

---

1,(x) = 2K 91 12 - K) +)+
Vx E 10,
2'
hx(x) 91 (2Kx KK,, v [0,
Building upon this transformation, follows that
fo hy o...o h1 where his E c (10, 11 1) for j E [12]. his E c (10, 1] (2K) for = 2,..., - and j E [[++1], and h E C (10, 1] ( jor [[+1]].
By Lemma 8, for each [1] and [4++1], there exist positive numbers C24, C25, C26, G27 and c2s and function fis E SFFNN LG, 8:), such that
2: 0
Ilfe, - h a S c2a Lo + C25 W
where L, = Lo. W; = c2,Wo. Wo. and = C28 Lo Wo. However, fü,j may not map into [0, 1]. For analytical convenience, we further transform into fr,, defined by
(x) min maxx (x). Vx.
This transformation can be implemented by adding two additional layers with four parameters to each fuz- Under this transformation, for each i E - and [t+11], there exists function f. E SFNN(WW + 2, I,, 4), taking values in [0. 1]. such that
IP. hl S C24 + W 0
We next construct f. which provides a consistent approximation to fo via the aggregation and composition rules in Lemma 9. For i
1], we aggregate fi,,,..,, into single function
fl defined by (f......... fiss.))) By composition, we construct f = f(o o .... .... = fi fi o ..... ..... f1 Then J (W,NNWW,,L s), with 3 = L; = 029 W max; f,+1 (W; + 2) = Wo and 8 8 = (*: + + = LoWoo The approximation error between fo and f can then be bounded as follows
then f (SFNN(te) + 2), L, 2, 1+ (s; + 4)).
Il fl
11 h-11 .... o h fi f( fillx
h1 1 .... o h - hy fi fill h filloc
(( 1| 0
oh fi 6 fi + [ fillo
<K(2K)
(1/11110 o pi o fill + 11 - f1lloo
1-1 <K
1-1
II (2
Σ
114 Pl + l14 filloo
1=1
where the second inequality follows from the definition of Hölder smoothness.
Since K II (2K) is constant, it follows that
Je S C32 Wi2 Lo + Wo w,r')²
||f folk2 < C32 max sEll
W,2 La + Wo
15
To ensure that the output lies within [K, K]. we define f by I min min 1100ll//l//lloo... 1) f. It follows that
II fullo = IIf al 1lf foll 211 folloo
Then, we directly establish the existence of constant Ca3 such that
Vonnen 1f(x) fo(x)))x(dx))
}
1/2
s css WWo--bo W.""
+ Wy"')
Choosing C17 = C33, C18 = C29, C19 = C30 and C2n = C31 completes the proof.
21

---

B.2 Proof of Lemma 4
In this subsection, we first formally define the architectures of modified recurrent neural networks (MRNNs) and their sparse variants (SMRNNs). First, we define a new activation function that operates on specific components instead of all components. Given an index set C N, the modified activation function øy is defined as
(s(s)
I, otherwise.
0r(s);
st
Using this modified activation function, for any E [L] define the modified recurrent operation
7/ (V-i); (V-)) + Bi Cc); Here, A.. Bi E RwxW & E ERWW and 71 -
(V-1)
={ a(AiNt (V-1) (V-1))
Biv + Ei); if I,
(AIT)- (V--1) Biv + +ū)) otherwise.
We denote by MRNN and,, (W, L, s)) the class of neural network functions defined with the structure of RNNs but composed of modified recurrent layers instead of original recur- rent layers. It is clear that RAN dg.d, (W, L,, K) S MRNN dn,dy (W, L, K, Similarly, SMRNN dz.dy (W, L, K, s) is defined the class following the structure of SRNNs, also utilizing modified recurrent layers.
The proof builds on Lemma 10 and Lemma 11, which construct two intermediate representations. Specifically, Lemma 10 shows that an SFNN can equivalently expressed as an SMRNN prediction function, and Lemma 11 further shows that SMRNN can be represented by an SRNN. Proofs of these lemmas are deferred to Appendix B.5 and Appendix B.6, respectively.
Lemma 10. For any SFNN E SFNN 4g x (W, L, s), there exists SMRNN prediction function f E SMRNN d. ((dx + 1)W, N + L, + 3(dx + 1)WN) such that
f(x) f(X), X [0, 10,1,
Lemma 11. For any SMRNN prediction function f SMRNN de (W, L, there exists SRNN prediction function f SRNNN de (W + 1,2L, + 6LW) such that
f(x) = f(X), X E [0. 1
Proof of Lemma 4. From Lemma 10, we know there exists SMRNN prediction function f EMMRNN,.1. ((dx + +1J,, N + L, 3(dx + 1)WN N) such that
j(x) f(X), X [0, 1
Then, applying Lemma 11 to the function 7, we obtain SRNN prediction function f
SRNN dz ((1x 1)W + 1,2L 2N, + 3(dx + I)WN + 6(dx + I) N)W), such that
f(x) f(X), X E [0, 1|
This completes the proof.
B.3 Proof of Lemma 6
To establish the lower bound, we utilize the Lipschitz continuity Pe Pelxx....XX (-), which implies that PeIXi.....X. (u) 2 p/2 when lu| < p2(2)) By applying Lemma S6 from the supplementary materials of Padilla et al. [2022], we derive the following result:
P p²
R₂ R- 2 min
E Emin [[ ...,.... XN) ..,...,,))....
4' 166
22

---

(f(X1..... .x)..xx) fo XN)) (B.1)
Furthermore, as both f and fu are bounded by K, the squared difference satisfies the inequality (J(Xi...., -Na Xx))² (X Xx)- Xx)..
By incorporating this result into (B.1) and assuming K 2 1, we arrive at the desired lower bound for the excess quantile risk.
For the upper bound, we have the following decomposition:
p-(u - v) p+ v(T iuu - 0)) / ( 8) 1(u S 0))
Recall that (f) E6pp (Yn f(X Xx))], Choosing and I.
,NN,, we have
R2( - (fo) =Elp, (YN -
-rr(e]
1 Pel X (s)dsd
=F
1f
P
2
where the last inequality follows from Assumption 1.
B.4 Proof of Lemma?
To prove this lemma, we first introduce the definition of the uniform covering number. Definition 5 (Uniform covering number). Let d E N and F= {f [ + R} be class of vellev
functions on X. For e >> 0, the uniform covering number of F under the supremum norm is defined as
Na(e; F)
sup (x).,...,)) 2d
N(e, 1lo,,, Flxi...xah)
where F Fixi,...x. {(f(x))).....)) f(xx)) E S Rd and N(e, || bloo, W) the e-covering number subset W S Fd under the supremum norm
lloe.
Given the above definition, we next present two lemmas characterizing the covering numbers of RNNs and SFNNs.
Lemma (RNN covering number bound). Let Fx a L, K).
Then we have
Kn
log (c. lloe, Fx) S 2 L log max L} log
Lemma 13 (SFNN covering number bound). Let FN = SFNN de, (W, L, K, s). Then we have
Kn
log N (s, Fx) S sa log(W log
The following lemma shows that any SRNN can be represented by an SFNN, which will be used to characterize the covering number of SRNNs.
Lemma 14. For any SRNN prediction function f SRNN de LW, s), there exists an SFNN SFNN X.Ndd ((2N - 1)W, (N 1)L 2, 2Ns 2N"WL) such that
f(x) f(X), X [0. 1|
Finally, we provide two technical lemmas that serve as essential components in the proof.
Lemma 15 (Theorem 7.3 in Bousquet [2003]). Let Zi...... be i.i.d. random vectors drawn from the distribution of Z, and J measurable class of functions such that EE = 0 for any f F. Assume sup || < A and let a be positive constant such that o² > supser Ef2 (Z). Then, for any x 0,
a
2x
4Ar
P
sup fEF
f(2))
2
2E
sup fEF
Σ
f(2i)
+ +
(B.2)
n1
12
3n
23

---

Lemma 16 (Corollary 5.1 in Cheroozhukov et al. [2014]). Denote S = ([0, 1 x R)" and let E S be i.dd random vectors drawn from the distribution of Z E S. Let F be measurable class of functions S + R. to which a measurable envelope F is attached. Assume that || < oc, any positive constant such that su (Z) S x² > F1 Furthermore, we assume that there exist constants A 2 e and v 2 such that sup& [c ll F)).: (A/c)" for any 0 < I 1, where the supremum is taken over all n-discrete probability measures Q on S and
Zn Z......
and let a >
(e, F) is the e-covering number of under the (2( norm. Then,
AllFFl2
2 vn
A||
E
sup fEF Vn
Σ
log
0
0
where F = mixisisn F(z:)-
Proof of Lemma 7. Recall that F(6) = (s E EF foll2 5.}.
Since Pr (-) is Lipschitz function, we have
sup SEF(8)
1g(f, 2:)1 < sup }EF(6)
1f(x x:) -t xx)] S 2KK.
Therefore, supjer(6) l0(f.z.) - B2g(f. f, 2) < 4K =: A. Moreover,
B(19(f.2:) < sup FEF(8)
sup
Bllf (41...... N..S-,N.....,.............. 88.
JEF(8)
which further implies
Ells(s - Ezg(J, 2)195
ElggS 2:1 s 00 =2
sup SEF(6)
sup SEF(6)
Denoting E(6) = ESupsef(M) In EE= ((1,z))
- Bzg(J, 2) Lemma 15 gives
P
sup (EF(A)
Σ
((f - -Ez(ff 2 2E(0) +
2r +
4Ar
(B.3)
3n
==1
for any 2 0.
Now, establish an upper bound of the expectation E(8). We denote Mn (o(ff 21), F())).
. Case 1: F = RNN (W, K)
Combining the Lipschitz continuity of p+ (-) and Lemma 12 gives that for any e (0,
log N, (s Il,, & log N₂ (c, F) S W2I2 log max (W, L) log Kn
Also, the Lipschitz property of p+() implies that F 2K is envelope function of M,,(5).. Thus, for any discrete probability measure Q supported on points, we have
log (e ll2, (2) S 2 log max(W, L) log
where the La norm is taken with respect to Q. Applying Lemma 16, we have
L2 log max W, L}
nK
2 log max W,, L}
nK
E(0) $ V
log
+ 2K
log
S 66, + 53
for any 8 2 1/n. Thus, when 8 2 5g. we have E(5) S 56g. By combining this and (B.3), there exists universal positive constant C34 such that
sup (EF(8)
Σ (ff2,)) Exig(S, Z)] T 1=1
2 034 -8 66 + V
SE
24

---

holds for any 0 S x I no² and 8 2 56.
Case 2: F = SRVN d2 (W, I., K). Combining the Lipschitz continuity of pr(+), Lemma 13, and Lemma 14 yields that for any e. E (0, K) and x WL,
Kn
log Nn. (e, -lo,, M,, log N₂ (e, (
(W log
Also, the Lipschitz property of p+(+) implies that F = 2K is an envelope function of (o). Thus, for any discrete probability measure Q supported on points, we have
log (c12f122 I1,, Mx.()) x sL log(W log
n1
2e
where the La is taken with respect to Q. Applying Lemma 16, we have
ST.
nn (nk) 0
ST. log(W
nk 0
E(6) S 01
log
+ 2K
log
n
n
$ 86, + o
for any 8 2 1/n. Thus, when 8 2 66. have E(8) S dōg. By combining this and (B.3), there exists
a universal positive constant C35 > 0 such that
P sup (ET(8)
g 2.) Ex 2)] 2 c35 - 8n + Sem
holds for any 0 S x S no² and 8 2 5s. Choosing Czs maxfcss, cas} completes the proof.
B.5 Proof of Lemma 10
-(N) For any given SFNN j. we first construct an SMRNN Pax+s such that T0N+1 is equivalent to the first layer of J. Subsequently, building upon Te sequentially construct Tex by directly leveraging the remaining
layers of f. ensuring that Tex+z is equivalent to Proof of Lemma 10. Our construction begins with a key result for building the initial SMRNN 78841
For any given Au, AR,N E Blxd= for E [W] and bias vector e E R.W, we construct an SMRNN TeN+1 = TN+1 o ONN o ......o1 T op: Rd XN + 1IW++ of width (dx + 1)W, depth N ++ 1. The layers {r
N+1 and are defined explicitly as follows for each t E [N] and I E [N].
A
B1,1
v-())
19 (V) =
r (t-1)
(V)
+
E (dW++))
A
Bwa
v(())
I
x(()
Or,
(0) (V) = +0), (x)
(())
E R (dx++11)
:
Id
(()
Or
0
Ode.dx Odn.1 01,dx
Idx b
Odw.!! ) b = En A,.AA E Rixd,
where A =
Bbd
[Auu 1S1,15N" = ( 1-1
{
i i [ (dx +
c= (osema) and Diss
1SKJSN"
otherwise
25

---

Combining Lemma 6 in Song et al. [2023] and Step 2 of the proof of Lemma 8 in Jiao et al. [2024]
gives that T0x+1 constructed above satisfies the following property:
j== A15x(0 +c)
N
j= x2,x +02)
(N) dd++
(X) =
(B.4)
E Aw,x(() + tow)
OdgW;1
In this construction, we have T(p) = dx. T(rN=1) S W + dx. and for / [N]:
T(n) = W(1 + dx +1+ |lbyllo)
- 2W (dx + 1).
We bound the sparsity of rox++ by
N-1
T(rex
T(m)+(()) s 2d +1)+ W + dx.
(B.5)
By applying this result given SFNN f. parameterized (vec( B1) Jc T ...... vee B1) [.,.) we proceed to define the corresponding parameters required for initializing the SMRNN F8 Here, B1 W B₁ WW for I 2...., C1 and for each i E [W] and j [N], we define Ar, (B 1:jd,- The notation (B1)a, (j- (111
1, Bg E RIXW C. E RW for I E [L 1],
and C9 Bd Specifically,
the i-th row of B, corresponding to entries from column (j 1)dx + 1 to jdx. By applying these definitions to the construction described (B.4), we obtain 10x It then follows from (B.4) that
refers to the sub-
0 (
B vec vec(())+ (
(N) 0044 (X) =
OdeW.1
We next construct the remaining L - layers using the 2-nd to I--th layers of the SFNN f. Define
(V)
(00)
{
o (Biv(t + +0) (iv()) +c, 1= L,
1= 2, 2....,,. L 1, a(V) (Qv ..... Qv())),
where
B OWLw.W
Ow.dww O&,Widww
B
for 1 - L
BL
dey OddxtIjW-d.M -dW
Br.
Odd+1IW- d
C₁
for 1 = 2,...,
Oaww.1
(O(dFtlw-ds)))
Q Tay Odd )
Then SMRNN ToNtL = =NONN+L . FN+2 o Tox+1 satisfies
(N)
(X) = f(x).
Finally, we compute the sparsity of PaN+z.
T(Fox=n)) =T(ronan) +
T(#++)+ +T())
1=2
26

---

22 + + +
[IB:llo + (c)
<d +(f),
where the first inequality follows from (B.5) and dy
B.6 Proof of Lemma 11
Proof of Lemma 11. This proof is adapted from the proof of Lemma 3 in Song et al. [2023]. We first show that any modified recurrent layer can be represented by the composition of sequence of recurrent layers and linear maps. Specifically, for any modified recurrent layer RdxN - RdxN
satisfying p(0 (X) o,( (Ar(((1)) + Bx(1) +c), there exist linear maps p RaxN Rd++1)xN
and IR 1)+1) - R along with recurrent layers r2. R(d+1xxN - (d++)xx such that
f(X) qor2 (XX), for X e X
(B.6)
where C IR is compact set.
We present the construction of rı, 12, P: directly as follows. For any t E [M], we choose
la
(11))
(x) =
x(1) =
O1dd
0
10 (V) 0 B
v(t)
+ 20 do1
1 (V)
(V) = of
[Arl (Arl--)(()) (V) + vil)
- k
02+1
(V) = la Oa,1 ( 11
lo-k
-lo-
v(0)
0
where V [(d+1)xN v(t) is the t-th column of V, A (
Ā
)
IN
-lo-k 10+))
Io-k
0
and 2g max supxexN (AF(--)(() + (x +))l
By the choice of zo, for any E [N]. we have
(ri p) (1) (X)
B
(1)
0
+
0
+ 2n 1d+1
20l4 (1x(1) + 2)
20
We next prove that our construction satisfies the following property by induction
7(0(())s)) 1:sk
(r2on (X) = 20 10- + F(t) (X) x)++d
for t E [N].
(B.7)
20
For / = 1, we have
(120r (2)
(zols
(n)p) (x)
(0=t11+)
(Bxx))
+ + 20ld
2011
+
04+1-K/
27

---

of (Bx(2) + 1.)1
20ld-k + (Bx(1) + c) k+idd
20
F()) (x) 1:k
20 1a- + 7(2) X 1idd l:d
20
where the second-to-last equality follows from the choice of 20- Assuming the induction hypothesis holds for - 1, we have
(Bx(1) + +C)
(11 (x) + (220L (0411-)) =
20 (Bx(1) + c)&+Idd
With the above align, we further obtain
(r2 OTI pp)
(X)
p(1 1) X 1:k6
20 2old-k + (Bx(1) + o)s+id ))
(Bx(") + +C)1..
20 20ld-k 7((-1)(X8
X
1++1dd
20
(Ar 1) (xx)kz (AAp(t-)) (x) x))++1dd
20ld-k + (Bx + C) ))
(Bx(1) + 2)
20
7(t) (1))) 1:ke
20 20-a-k + 7(t)) (X)))))) x))+1:d
20
where the last equality follows from the choice of 20- Thus, we have demonstrated that our construction satisfies (B.7). By the definition of for any / E [N] we have
(gor2 O (X)
Jk =( Ta O4 ) (
Ta-k -la-k ) (
7(0) x))x 1:k
20ld- + 7(1) (X) (x)(()) x)++1:
20
(X).
Thus, construction satisfies (B.6).
Applying (B.6) to any SMRNN ř. = -ořL o op with depth L and width W and setting W. we know that there exist 2 recurrent layers ..,..., F2L and 2L linear maps P PL, Q1 ........ QL such that for any X E [0, 1 XNN
(x = Fr.
x
= (qr.21 1PL) (g1r2r1p1) 7 X
=(ğqL) r22 (T2L- 1PLL-1 o (r3p2q1) o 02 or (X)
By setting q = Tgr, p = PIδ. u2 = 1P411 for I [L - 1}. and re = q o r2r. f21-1
F2L-2 ... or4ofs Ti p. we obtain
ro (x (X), for X E [0, 1|
It is clear that ro has depth 2L and width W + 1. Finally, we analyze the sparsity, for E 1).
((22+1) ((22+114419) ((r241) Billo ق 20w 1w +
(B.8)
28

---

T(22) = ||All + (W + + k),
(B.9)
| ||Arllo W.
(B.10)
Combining (B.8), (B.9), and (B.10), obtain
L-1
T(ro) (To)) +
((2)) +T(r)) Σ ((244)
2L
=T(ri)) + T((n)) +
T(n)
1=1
ST(P1)) + T(9g))
(I|Arllo + || + ( + 4LW
1=1
=T(p) + T + + + Σ T(h) + 4LW +
ST 6LW,
where the first inequality follows from +
+1 S W and C + zowwllo 1 |allo + W.
B.7
Proof of Lemma 12
Our covering number bound for RNNs is consistent with Lemma 12 in Jiao al. [2024], which itself based on Theorem 7 in Bartlett et al. [2019]. For completeness, we provide full proof of the covering number bound for RNNs.
Proof of Lemma 12. By Lemma 9 in Jiao et al. [2024], we have
FN RNN dr (W, K)
1)W, (N FNN d.NNd..N.. (1N
K).
Then it follows that
log N,, FN) S log Nn
1,W, (N + 1)L,
(B.11)
By Theorem 7 in Bartlett et al. [2019], for any W, L, K 2 0 have
Kn
log N,, 6, IINN NN (W, WW, k)) x W log max max(W, T.) log
(B.12)
Combining (B.11) with (B.12) yields the desired bound and completes the proof.
B.8 Proof of Lemma 13
We first present the following lemma and definitions of Vapnik-Cherronenkis dimension and Pseudo- dimension of real-valued function class
Lemma 17 (Lemma 17 in Bartlett et al. [2019]). Suppose W M and let Pi...., PM be polyno- mials of degree at most D in W variables. Then we have
| ( (a))) a R 2(2cMD/W)
Definition 6 (VC-dimension). Let F Rdx + be class of real-valued functions. The VC- dimension of F, denoted as VCdim(F), is the largest integer E N for which there exist points (Xi Xm) E Rd xm such that
sgn(f(x,n) ((s (f = 2
29

---

Here, sgn(-) is the sign function defined as
1, if x > 0,
sgn(1)
00 if x S O.
Furthermore, the quantity
[1)(m) = (s sgn(f((m))) E F}
referred to as the growth function of the funetion class F.
Definition 7
Given real-valued function class F Rd - R, the pseudo-
dimension, denoted as Pdim((F) is the largest m N for which there exist (XI,..., X,n E IR xm
and (y1,... (mm) E R such that for any (bi,...,bun bm) E {0,1)* mm, there exists / E J such that for any E [m]:
f(x:) 144
b; 1.
Proof of Lemma 13. We prove this bound in three steps. First, we reformulate the problem as an inequality depending only on the VC-dimension, using the relationships between the pseudo- dimension, the VC-dimension, and the covering number. Second, we adapt the proof of Theorem 7 Bartlett al. [2019] to establish a bound for the growth function of the prediction function class derived from sparse neural network with a fixed zero- poo-paraeeter configuration. Finally, we account for all possible sparse structures derive the overall VC-dimension bound for the prediction function class.
Step 1: From Covering Number to VC-dimension.
By Theorem 12.2 in Anthony and Bartlett [2009], for any E (0, K), have
Pdiml FN)
K
No (1 Fl),
Σ
If n < Pdim(Fx), then
Pdim(
"enK"
(a 11o, FN)
1
On the other hand, if 2 PdimmFN). then
enK (Pdim(FN)
Pdim(FNr
enk
Pdim(FFx)
(e, Fx)) <
Thus, we have
enk
log No. (e: |1,o FN) < Pdim(FN) log
(B.13)
By Theorem 14.1 in Anthony and Bartlett [2009], there exists function class FN+1 generated from the SFNN with width W, length L + and sparsity + 1 such that
Pdim(FN) < VCdim(Fxx11)
(B.14)
In the remainder of the proof, we restrict our analysis to the case of N for simplicity. The results obtained can be directly extended the case of N +
Step Cardinality Bound for SFNNs with fixed sparsity structure.
First, we formally define the function class generated from the SENN with fixed sparsity structure as follows:
F (fo E E 0 = 0 if = 0, Vi e [M].
30

---

Here, M EE=' Wt El=o W1 is the total number of parameters in the SFNN and & E {0.1)A a binary vector satisfying llEllo Each entry of & corresponds to an entry of 0 in the SFNN: f then 0, treated as free parameter; otherwise, ç = 0, 0, is fixed to zero. Thus, FV represents the subclass FN constrained by the sparsity pattern specified by &.
In this step, we consider a function class F% with a fixed parameter &. Since FV is a parameterized function class whose parameter vector has M zero entries, any function in FV can be represented as g(x, 0¹) where 0' E R Given m 2 s and ; our goal is to derive an upper bound for the cardinality of the following set:
(sgn(g(xx, (gmp. f))) 10 R-}
Recall that the ReLU activation function (ax) max 0} is piecewise-defined function. Con- sequently, (x, 0*) is also a piecewise-defined function with respect to 0'. To apply Lemma 17, we partition IR into disjoint regions S (/1..... ]p). where 1 I 0 for + and i,j (D]. Within each region 4, for all E [m], the function g 8¹) reduces to a polynomial function. We denote the number of free parameters in the first layers by s' We construct this partition recursively with the goal for S₁ to satisfy the following properties for any / E [L]:
1. Ifl = then |S|-1| otherwise, |8,-1| - |S:221 . 2 /2eWmll/s
2. For any I St-1. w E [Wi]. and j E [m]] the input to the w-th neuron at level / is polynomial function fı,w., of degree most which depends on s variables.
We first select S. = IR. The input to any neuron at the first layer is an affine function degree at most 1, depending at most o variables. Consequently, Sn satisfies both required properties.
Given Sn- that satisfy the two properties, we proceed to construct Sn- For any region E Sn-1- Since fi is polynomial function of degree most n for each [W.] and [m]. we derive the following inequality by applying Lemma 17.
Hsgn(ft..,,(P)) E [W.),j [m], 0' E Ty 2(2eW,,mn/s2
This result implies that can partition I into most 2(2eW (s})) disjoint regions. Within each region, all these polynomials maintain the same sign. Consequently, for any j E [m] the output of any neuron at the ni-hh level remains a polynomial of degree at most n.
We define Sn be the set of all disjoint regions generated from all / E Sh-1- Since each layer connecting the n-th and (n + 1)-th layers can increase the degree by at most one, for any E [m]. the input any neuron the (n + 1)-th layer is a polynomial of degree most + 1, depending 2+1 variables. Thus, Sn satisfies both of the required properties. By induction, we have at most s
L-1
S S II 2(2eW
II
2(2eW' 222e''lss) 1=1
1=1
where W = max(W,dy) As implied by Lemma 17, the cardinality of any partition in St-1 is bounded by ((2eW"ml/ss)) By AM-GM inequality, we obtain
(sgn(g(xi, 0')) "..).....gn(((.,.......))))))).. 0')) O' R*H
2(2eW" ml/ /s/
1=1
(2eW'm m ΣE 1 El
<24
(B.15)
El=
Step 3: Obtain bound on VCdim(Fx).
31

---

In this step, we consider all possible fixed sparsity structures &, of which there are at most (w) - M
Summing (B.15) over all such structures yields, for any (x1, Xm) E IR xm
(sgn([(x1)),
(xx())))..,s(((((mxm)) 2eW'm m Ei=1
() |/ FNNH
Elsr
Σ
2
Σ $/$
l|Ell==
Ii-1
2eWW'm m 11=1 1
max IEllo=s
Ei=1 $
(B.16)
Let Emax be the sparsity structure that attains the maximum in (B.16):
2eW' m
Et
Eis
Emax
arg max
Σi=1 S''
I llo==
Lx Emax
22W'm Dtsr Ei cmax
The upper bound in (B.16) can thus be expressed as M
Choosing mn VCdim(Fx), have
2 VCdim(Tx) S M 2eW'VCdim(Fx) Et= 1)
Eh 2 t
Σt= 86
By Lemma 18 in Bartlett al. [2019], obtain
VCdim(FNi s + s log M + ② S Emax logz JoseeeW )
< + s logz (W2L) L) + Ls log, 2oW"
$ SL log( ( 1.2).
(B.17)
Finally, ensure the completeness of this proof, we need to show VCdim(FN) 2 8. This result follows directly from Theorem 3 in Bartlett et al. [2019].
Combining (B.13), (B.14), and (B.17), we obtain the desired bound.
B.9 Proof of Lemma 14
This proof is adapted from the proof of Lemma 9 in Jiao et al. [2024]. As first step, we demonstrate that for any recurrent layer r: RWxX x + RWxx satisfying (()) (X) (Ar(t-2) Bx (x +c), there exists an SFNN with width (2.N 1)W and length N + 1, such that
,()) (X)
(x)
f(x)
(B.18)
x(x)(
Here, X E WxN, We construct f as f = fN+1 fN o... fi with the parameters defined as follows
B 0 0
0
0
C
O -1
0
Á
b
0 0
0 0
-I
32

---

0oo.!
I
-1
0
0
b
j = 2. N,
A
-B
1
1
AN+1=
b = 0.
-
In these block matrices, each block of size W W. Å₁ has 2N rows and N columns of blocks.
For = 2,..., N, A, has 2j - 4 identity matrices in the upper left comner and 2N 2j identity matrices in the lower right corner, and by has (2j 2)W zeros above AN+1 has N rows and 2.V columns of blocks.
We verify that this construction satisfies (8.18) through direct calculation.
x(1)
x(2)
f(x) fN+1 fN fxo..off
x(()
(()((()))
o
x(-r((())
a(p(i-1) (x))
a(-r(i-1) (x))
fN+1 o fN o... fit1
(x)
((x(+t1)
ol a(-x((t1)
a(x(N))
a(-x((i))
((rl)(x))
o a(-r("(())
-NN+
(x(N-1)(x))
of a(-r(N-1)(2) (x))
(())(x)
x()
(x)
(()(
This derivation relies on the property (xx) - (-x)= =x x, which applied sequentially across the
last three equalities. Thus, our construction for satisfies (B.18).
33

---

Next, we extend this construction to complete SRNN та = q ..... 01 p with width W and depth L. For an input sequence [0,1]***N, we can find SFNNS fi. fr. with width (2N - 1 and length N + 1, as well as layers p and q, such that
aod fea xx)
fo........ ho - fr ... ... hoh
(XX SOAA
(X) Caiiii Collad Colope ())
O (x)a (ollddo tu)
ofo ..... o Jr. q =
- Bolandd (do 44
(x ())(dd (i To w.................... . 14
) (x) lel lo....... 12) 14)
do =
(T),
= 1qopz 1.. o (x).
Here, the layer p constructed map X to the stacked outputs of p applied to each input x(() and the final layer ạ applies to the N-th segment of its input vector. Choosing f = qo oo. o ff op completes the construction.
Finally, we analyze the sparsity
T()) -
+ +2 2N +
1=1
S 2.NT(rr) (re) + WL.
In conclusion, we obtain
TESENN
- (1W, (N 1)L+ 2,2Ns +2N
C GDP Growth Analysis
Analyzing GDP growth is cornerstone of economic research, providing crucial insights into economy's overall health and trajectory. Given its capacity to effectively capture and forecast upside and downside economic risks, QR has become indispensable tool in GDP growth analysis [Adrian et al., 2019].
To further evaluate the practical advantages of RNN-based QR estimators under potentially nonstation- ary settings, we conduct an empirical analysis using real GDP data We compare the out-of-sample prediction performance of RNN- SRNN- FNN-, and QRF-based estimators for one-quarter-ahead GDP forecasting. For all models, the input sequence length is fixed at N = 4, corresponding to four consecutive quarterly GDP observations x,-3 to x,, representing target variable y, is the GDP level at time + All models were implemented according to the specifications detailed in Section 4.
year of lagged values. The
The dataset partitioned chronologically, with the most recent 30% reserved for testing. For all NN-based estimators, the preceding 70% is used for training and validation. During training, 20% of this subset is held out for validation, and early stopping is applied based on the validation check loss.
*We utilize U.S. GDP growth data from April 1947 to December 2024, accessible at https://fred. stlouisfed. stwouisfed.oorg/eeeeeeeeeeeeee wtt.:///erreeee1/11112226656555.
34

---

All models are trained by minimizing the empirical check loss and evaluated on the test set using the same metric.
Table presents the mean empirical check loss for each model at quantile levels 7 = 0.1, 0.25, 0.5, 0.75, and 0.9. The results indicate that the RNN- and SRNN-based estimators consistently outperform both the FNN and QRF methods across most quantile levels. Furthermore, the SRNN achieves predictive aceuracy comparable that of the standard RNN while inducing sparsity. This suggests that the sparse architecture does not compromise performance.
Model
T = 0.1
7 = 0.25
T == 0.5
T == 0.75
T 0.9
QRF
0.849
1.225
1.410
1.246
0.911
FNN
0.867
1.180
1.773
2.505
2.657
RNN 0.835
1.113
1.349
1.154
0.904
SRNN
0.837
1.700
1.211
1.200
0.898
Table 2: Out-of-sample prediction errors at different quantiles for GDP growth analysis.
To assess the robustness of RNN-based estimators in QR in the nonstationary environment, we extend our analysis include in-sample estimators, enabling a direct comparison of in-sample and out-of-sample performance for both RNN and SRNN methods. Adrian et al. [2019] showed that GDP growth volatility primarily driven by lower quantiles, while upper quantiles remain relatively stable over time. Motivated by this, focus our comparison quantile levels 0.25, where tail behavior is most prominent. For the in-sample predictions, the full dataset used for training and validation. The comparative results, shown in Figure 3, reveal striking alignment between sample and out-of- sam-lle quantile estimates. In these figures, the r-axis represents time, while the y-axis denotes the GDP growth rate. This consistency is particularly noteworthy, given that major financial crises (c.g., 2007-2009), which represent substantial tail events, are absent from the data used for out-of-sample estimation. These findings highlight the reliability and generalizability of RNN-based QR methods, even under adverse and previously unseen market conditions.
0.05, 0.1, and
0.25
-1
-2
-2
-3
-3
-4-
-4
-5
Oul-o/-sample no-sample
-5
In- Ra-sample
1950
1960
1970
1980
1990
2000
2010
1950
1960
1970
1980
1990
2000
2010
Figure 3: In-sample and samplee comparison for RNN and SRNN models.
D Hyperparameter Sensitivity Analysis
In this section, we conduct additional experiments varying network depth L, width W, and pruning ratio on Model I in Section 4 under the f2.25 distribution assess hyperparameter sensitivity.
We first vary E 4, and W {16, 32, 64, 128. 256, 512, 1024], evaluating each (L. W) configuration under the same data and training budget. As shown in Table 3, performance generally improves as both the network depth L and width W increase. When either L or W small, the model tends to underfit, leading to worse MSE. In contrast, when both are sufficiently large, the model achieves best performance, suggesting that adequate depth and width are essential for capturing the heavy-tailed characteristics of the 25 distribution.
35

---

L W
16
64
128
256
512
1024
0.056
0.054
0.054
0.053
0.060
0.065
0.097
0.055
0.051
0.051
0,050
0.054
0.058
0.073
0.055
0.050
0.052
0.050
0.051
0.052
0.070
0.050
0.050
0.050
0.052
0.050
0.050
0.057
0.060
0.050
0.053
0.050
0.049
0.049
0.056
0.050
0.050
0.050
0.050
0.049
0.049
0.056
Table 3: SRNN MSE across different and W.
We next study sparsity while fixing other settings as Section MSE exhibits non-monotonic trend with respect to the pruning ratio: moderate pruning improves performance, suggesting an effective regularization effect, whereas aggressive pruning degrades performance- likely optimum that balances compactness and representational capacity.
reported in Table 4, the
likely due to the removal of critical parameters. These results indicate an interior
Pruning ratio
0.2
0.3
0.4
0.5
0.6
0.7
0.8
MSE
0.046
0.044
0.040
0.039
0.038
0.036
0.038
Table Effect of pruning ratio on SRNN MSE.
36

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: The abstract and introduction accurately stated the contributions and scope of the paper.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation long is clear that these goals are not attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: Our work focuses on fixed quantile level. We discuss this limitation of our results and future direction in Section
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (c.g.. independence assumptions, noiseless settings, model sell-speciiiatiinn asymptotic approximations only holding locally). ). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-tex, system might not be used reliably to provide closed captions for online lectures because it fails handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimesss.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community, Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
37

---

Answer: [Yes]
Justification: All the theoretical results the paper are accompanied by a full set of assumptions, with complete and correct proofs included in the Appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proof' in the paper should be numbered and cross referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors encouraged provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: We clearly disclose all information necessary for reproducing the presented experimental results, thereby ensuring their reproducibility.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, other means that appropriate to the research performed.
also be provided via detailed
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) If the contribution is a new model (e.g., a large language model), then there should either way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reprodueibility may be tricky in some cases, in which authors are welcome to describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path reproducing or verifying the results.
5. Open access to data and code
38

---

Question: Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [No]
Justification: Our experiments easily reproducible, and we provide all necessary informa- tion for their re-implementation in the paper. As this paper focuses on theoretical analysis, the code is not provided.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips.c cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for a new open-source benchmark).
The instructions should contain exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https: / /nips. cwtppbliccguiiee..e wwalcc/guiee///ee //co//euumiisiionnPPP.
for more details.
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental setting/details
Question:
Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary understand the results?
Answer: [Yes]
Justification: This paper provides all details essential for understanding the results in Section 4, Section 5, Appendix C, and Appendix D.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: This paper reports appropriate information about the statistical significance.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emoor bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
39

---

The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [No]
Justification: Our experiments utilized relatively small datasets with minimal computational requirements that can be easily reproduced on standard personal computers, making detailed specifications of compute resources unnecessary.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary failed experiments that didn't make into the paper).
9. Code ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https: //nerr//////rs..../
Answer: [Yes]
Justification: We have reviewed the NeurIPS Code of Ethics and ensured that our paper conforms to it.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification:
Due its theoretical nature, this paper does not discuss societal impacts.
Guidelines:
40

---

The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No. they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negatíve societal impacts include potential malicious or unintended uses (e.g.. disinformation, generating fake profiles, surveillance), fairness considerations (c.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of MI.).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: do not make publicly available any data or models that pose a high risk of misuse.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: The datasets used in this study are the DJIA and GDP, both of which are publicly available.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
41

---

For scraped data from a particular source (e.g., website). the copyright and terms of service of that souree should be provided.
assets are released, the license, copyright information, and terms use in the package should be provided. For popular datasets, paperswi thcode. com comddataees has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
We don't introduce any new assets in this paper.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/mmodel as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification:
The paper does not involve crowdsoureing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsoureing nor research with human subjects.
42

---

Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if' applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Does the paper describe the usage of LLMs if is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Question:
Answer: [NA]
Justification: Our research does not involve LLMs as any important, original, or non-standard components.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non-standard components.
Please refer to our LLM policy (https /tt////eurpppppepe mweuripp.cc/Coofeeeeeeeee ec/eonnnrrencee/22255/ for what should or should not be described.
43


Reviews:
[{"review_id": "wfQ9ZIDNcg", "reviewer": "Reviewer_aMFs", "summary": "This paper studies nonparametric quantile regression using recurrent neural networks (RNNs) and sparse recurrent neural networks (SRNNs) to approximate the true conditional quantile function, under compositional  model.", "questions": "Why consider sparse RNNs? What are the benefits of sparse RNN over RNN ?\n For me, a neural network class with sparsity constraint is hard to train in practice.", "limitations": "This is a theoretical paper; however, the primary approximation results are merely a straightforward combination of existing ones. Consequently, the contribution lacks sufficient novelty.", "rating": 3, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 3, "originality": 2, "strengths_and_weaknesses": "Strengths:\nThe paper derives nearly optimal rates for ERM estimators based on ReLU RNNs and SRNNs under certain mixing conditions on the data.\n\nWeaknesses:\nThe main result hinges on bounding both the approximation error and the generalization error. For the generalization error under mixing sequences, the authors rely on existing results. As for the approximation error, the contribution does not meet the publication standards of NeurIPS, as it appears to be a straightforward combination of existing results from Kohler and Langer (2021), Schmidt-Hieber (2020), and Jiao et al. (2024).", "was_revised": true, "final_justification": "Thanks for the author’s reply. I am not convinced by the evidence presented regarding sparse constraint neural network training, as it is theoretically difficult to provide guarantees for this type of training. I will maintain my score.", "rebuttal": "Thank you for the insightful feedback. Below, we clarify our key contributions and highlight how our approach differs from existing work, in the hope of addressing your concerns.\n\n---\n**W1.** The main result hinges on bounding both the approximation error and the generalization error. For the generalization error under mixing sequences, the authors rely on existing results. As for the approximation error, the contribution does not meet the publication standards of NeurIPS, as it appears to be a straightforward combination of existing results from Kohler and Langer (2021), Schmidt-Hieber (2020), and Jiao et al. (2024).\n\n**A1.** We thank the reviewer for their careful reading and understanding of the theoretical results. We acknowledge that our approximation bound for RNNs builds on known approximation results in the literature (e.g.,  Schmidt-Hieber (2020); Kohler and Langer (2021); Jiao et al. (2024)). However, our results for SRNNs are not a direct combination of existing work. In particular, Lemmas 4 and 14 in the supplementary material provide a nuanced characterization of the equivalence between SRNNs and SFNNs, grounded in a detailed analysis of neural network sparsity. To the best of our knowledge, this equivalence has not been established in prior work. Moreover, it is more general and nuanced than Proposition 2 in Jiao et al. (2024).\n\nFurthermore, our oracle inequalities (which include generalization error bounds) for quantile regression under dependence are entirely new and technically nontrivial. To better situate these contributions, we begin with a concise overview of our analytical approach, highlighting the key innovations as they arise at each step.\n\n- **Novel decomposition.** We begin by introducing a donut-shaped decomposition which has not appeared in the context of nonparametric regression with neural networks. Specifically, we introduce donut-shaped sets (Line 514), which allow us to decompose the probability bound $\\mathbb{P}(\\|\\hat{f} - f_0\\|_2>\\delta\\_\\star)$ into manageable components by bounding each term  $\\mathbb{P}(\\hat{f}\\in\\mathcal{D}_i)$ separately. This decomposition is different from the direct argument used in Eq. (16) of Jiao et al. (2024), whose analysis crucially relies on the nonnegativity of the squared loss—namely, $$\\mathbf{E}[(f(x_t)-y_t)^2-(f_0(x_t)-y_t)^2]=\\mathbf{E}[(f(x_t)-f_0(x_t))^2]\\ge0,$$ where $f_0$ denotes the true function and $y_t=f_0(x_t)+\\epsilon_t$. The check loss function, however, does not satisfy this property. To address this, we develop the novel decomposition described above.\n\n- **Refined blocking technique.**  Next, we introduce a refined blocking technique to relate the mixing sequence to its i.i.d. counterparts, which differs from the approach in Jiao et al. (2024). In Step 2 of their proof of Theorem 13, part of the data within each partition is discarded to facilitate analysis under the squared loss, leveraging the nonnegativity property noted above. Consequently, we design a new partitioning procedure (Figure 2) that retains all observations, enabling a more comprehensive analysis under dependent data. Moreover, while Lemma 16 in Jiao et al. (2024) plays a central role in their argument, it cannot be directly applied to our setting. Instead, we employ a probabilistic counterpart (Lemma 5) to carry out our analysis.\n\n- **Sharper inequality.** Finally, we develop a novel and sharp empirical process inequality (Lemma 7) that underpins the tightness of our oracle inequality. One reason the result of [1] lacks tightness is that their analysis relies on a non‑sharp application of the Bernstein inequality. To address this limitation, we derive Lemma 7 by building on Theorem 7.3 of [2] and Corollary 5.1 of [3].\n\nIn summary, while we build upon foundational results, we respectfully highlight that our results for SRNNs and oracle inequalities (including generalization error) involve nontrivial theoretical innovation. We will emphasize these contributions more clearly in the revised version.\n\n---\n**Q1.** Why consider sparse RNNs? \n\n**A2.** Thank you for this question. Sparse neural networks—both in theory and practice—have received increasing attention due to their benefits in generalization and computational efficiency. \n- **Theoretical side.**\n    - **Reduction of hypothesis space.** In the neural network context, sparsity significantly reduces the effective hypothesis space, yielding improved generalization bounds ([4-5]; Schmidt-Hieber (2020); [6]). For example, sparsity constraints shrink covering numbers and Rademacher complexities, directly tightening oracle inequalities and sample‑complexity results ([7-8]). \n    - **Lottery ticket hypothesis.** Theoretical work on pruning ([9-10]) demonstrates that large dense networks often contain much smaller sparse subnetworks—the so‑called “lottery tickets”—that achieve comparable accuracy.\n- **Practical side.** \n    - **Reduction of computational cost.** Sparsity drastically reduces the number of active parameters, which cuts FLOPs and shrinks model size without harming accuracy. Early pruning and compression studies like [11-12] demonstrated that pruning dense weights into sparse patterns can significantly cut model size while maintaining accuracy.\n    - **Deployment on resource-constrained hardware.** Sparse neural networks enable deployment on edge devices and low-power environments by reducing memory bandwidth and compute load ([13-14]).\n    \nIn summary, the combined theoretical and practical advantages provide strong motivation for our focus on SRNNs.\n\n---\n**Q2.** What are the benefits of sparse RNN over RNN? \n\n**A3.**\nThank you for this thoughtful question. In recent times, there has been an increasing preference for more complex neural network architectures, which has led to a growing consideration of sparse neural network structures. Our work demonstrates the advantages of sparse RNNs from both theoretical and practical perspectives.\n\n- **Theoretical side.** A comparison of Theorems 4 and 5 reveals that while RNNs and SRNNs require a comparable number of neurons to achieve the optimal rate, SRNNs admit a more complex architecture that better reflects practical network designs.\n\n- **Practical side.** In Section 4 and Appendix C, we complement the theory with empirical evidence. Numerical experiments (Figure 1) and a real‑world GDP forecasting study (Table 1, Figure 3) show that SRNNs match the predictive performance of standard RNNs while pruning the smallest 40% of parameters. This parameter reduction directly translates into lighter models, lower storage costs, and improved efficiency.  \n\nIn summary, SRNNs retain the modeling power of dense RNNs but deliver tangible benefits.\n\n---\n**Q3.** For me, a neural network class with sparsity constraint is hard to train in practice.\n\n\n**A4.**  Thank you for sharing this perspective. Sparse training has become much more accessible in recent years, with the development of user-friendly tools and methods.\n\n- **Built‑in support in common frameworks.** Popular libraries like PyTorch already include straightforward pruning utilities (e.g., `torch.nn.utils.prune`), which allow users to take a standard trained network and prune weights or even whole neurons with just a few lines of code. Similarly, TensorFlow’s Model Optimization Toolkit provides off‑the‑shelf functions for both unstructured and structured sparsity. These tools have made it easy for practitioners to experiment with sparsity without having to implement complex algorithms from scratch.\n\n- **Gradual and structured methods.** Beyond these basic utilities, there is also a growing body of research on more sophisticated but still practical techniques ([15-16]). For example, [15] proposed structured sparsity methods that remove entire neurons or channels.\n\nIn summary, modern frameworks make it easy to experiment with sparsity, and established methods from pruning and gradual sparsification enable SRNNs to be trained with minimal additional complexity.\n\n---\n\n\nWe will revise the paper to more clearly highlight our key contributions. If our responses have satisfactorily resolved your questions and addressed your concerns, we kindly ask you to consider updating your score. Thank you again for your time and thoughtful feedback!\n\n\n**Reference**\n\n[1] Deep quantile regression: Mitigating the curse of dimensionality through composition, ArXiv 2021\n\n[2] Concentration inequalities for sub-additive functions using the entropy method, In Stochastic Inequalities and Applications 2003\n\n[3] Gaussian approximation of suprema of empirical processes, The Annals of Statistics 2014\n\n[4] Approximation and estimation for high-dimensional deep learning networks, ArXiv 2019\n\n[5] Adaptivity of deep ReLU network for learning in Besov and mixed smooth Besov spaces: optimal rate and curse of dimensionality, ICML 2019\n\n[6] Quantile regression with ReLU networks: Estimators and minimax rates, JMLR 2022\n\n[7] Spectrally-normalized margin bounds for neural networks, NeurIPS 2017\n\n[8] A PAC-Bayesian approach to spectrally-normalized margin bounds for neural networks, ICLR 2018\n\n[9] The lottery ticket hypothesis: Finding sparse, trainable neural networks, ICLR 2019\n\n[10] Proving the lottery ticket hypothesis: Pruning is all you need, ICML 2020\n\n[11] Learning both weights and connections for efficient neural networks, NIPS 2015\n\n[12] Deep Compression: Compressing Deep Neural Networks with Pruning, Trained Quantization and Huffman Coding, ICLR 2016\n\n[13] GPU Kernels for Block-Sparse Weights, ArXiv 2017\n\n[14] Eager pruning: Algorithm and architecture support for fast training of deep neural networks, ISCA 2019\n\n[15] Learning structured sparsity in deep neural networks, NeurIPS 2016\n\n[16] Exploring sparsity in recurrent neural networks, ICLR 2017"}, {"review_id": "BpfgsUq2Fr", "reviewer": "Reviewer_iRDm", "summary": "This paper analyzes the statistical behavior of quantile regression (QR) estimators based on recurrent neural networks (RNNs) and their sparse variants (SRNNs) for time series prediction with dependent data. In contrast to earlier work that often assumes independent and identically distributed (i.i.d.) data or short memory processes, this study considers the more general case of stationary $\\beta$-mixing sequences, including both exponential and algebraic decay. The authors provide non-asymptotic error bounds and prove oracle inequalities for empirical risk minimization (ERM) using these neural architectures. They also derive minimax-optimal convergence rates (up to logarithmic factors), showing that statistical dependence slows down learning rates in a way that depends on the mixing rate parameter $r$. These results demonstrate that neural QR methods remain effective in dependent settings and extend existing theory on neural network regression.", "questions": "### Questions and Suggestions for the Authors\n\n- **Origin and intuition of the $(1 - 1/r)$ factor**:  \n  The convergence rate under algebraic $\\beta$-mixing is slower by a factor of $(1 - 1/r)$ compared to the i.i.d. case. Could you explain where this term comes from in your analysis? Is this slowdown an expected outcome when using standard blocking arguments to handle dependence, or is it specific to your proof technique? A clearer intuition for how the mixing decay parameter $r$ enters the rate would help readers interpret the result.\n\n\n- **Theoretical role of the blocking parameter $\\ell$**:  \n  The error bounds rely on a blocking technique with parameter $\\ell$ to control dependence under $\\beta$-mixing. Could the authors provide more theoretical insight into how $\\ell$ should scale with $n$ and the mixing rate $r$ to optimize the trade-off between approximation and dependence errors? Additionally, is the use of blocking sharp or conservative in this context — for example, could a different dependence-handling argument lead to faster rates?\n\n\n- **Estimating the mixing rate $r$ in practice**:  \n  The theoretical results depend on whether the process is exponentially or algebraically mixing, and on the value of $r$. How could $r$ be estimated or checked from real data? If the wrong type of mixing is assumed, how much would that affect model selection or prediction accuracy?\n\n- **Effect of nonstationarity**:  \n  The results assume that the process is strictly stationary. However, in real applications, some time series may show slow changes in their distribution. How much do your results rely on the stationarity assumption? Could the method be extended to models that are only locally or piecewise stationary?\n\n- **Effect of stronger dependence (smaller $r$)**:  \n  Since the theory suggests that performance depends on $r$, it would be helpful to show a small simulation or example that changes the strength of dependence (e.g., by using different values of $r$). This would make the theoretical results easier to understand for applied researchers.\n\n- **Failure cases**:  \n  All examples shown in the paper favor your method. Could you include some experiments where the method might not do well, such as very short time series, weak signals, nonstationary data, or very extreme quantile levels?\n\n- **Computation time and scalability**:  \n  It would help to compare the time it takes to train and evaluate each method (QRF, FNN, RNN, SRNN). \n\n- **Hyperparameter sensitivity**:  \n  Neural network methods can be sensitive to their architecture. A small study that shows how performance changes when:\n  - the number of layers $L$ is changed,\n  - the width $W$ is increased or decreased,\n  - the pruning ratio (40%) is changed,  \n  would help readers understand how stable the method is and how it should be tuned.", "limitations": "Yes", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "### Strengths and Weaknesses\n\n**Strengths**\n\n- **Theoretical novelty under dependence**: The paper extends quantile regression theory to stationary $\\beta$-mixing sequences and provides convergence rates that adapt to both exponential and algebraic mixing decay. This represents a meaningful and nontrivial advance beyond existing results for i.i.d. or weakly dependent data.\n\n- **Precise characterization of mixing effects**: The paper clearly shows how the decay exponent $r$ of the $\\beta$-mixing coefficients influences learning rates. In particular, algebraic mixing slows the rate by a factor of $(1 - 1/r)$, and the rates gradually approach those for the i.i.d. case as $r \\to \\infty$.\n\n- **Transparent connection to neural architecture**: The convergence results are explicitly linked to the width and depth of the RNN/SRNN architectures, as well as to sparsity in the SRNN case. These modeling choices are well motivated and consistent with prior work on neural function approximation.\n\n- **Clear comparison with related work**: The paper positions itself clearly with respect to existing literature, emphasizing the novelty of its dependent setting and the use of recurrent network structures.\n\n**Weaknesses**\n\n- **Limited intuition for role of $r$**: Although the dependence of the convergence rate on $r$ is mathematically clear, the paper does not offer much explanation for why this dependence arises. A more detailed discussion of the blocking technique and the trade-off between bias, variance, and dependence would improve understanding.\n\n- **Strong theoretical assumptions**: \n  The paper's theoretical guarantees rely on structural assumptions that may be restrictive or difficult to verify in practice. Specifically:  \n  - The target quantile function is assumed to lie in a hierarchical interaction model with known smoothness and sparsity structure ($\\mathcal{H}_{d_x \\times N}^l(\\mathcal{P}, K)$; Assumption 2), which may not hold in real data.  \n  - The conditional distribution of the error term must be Lipschitz and uniformly bounded away from zero (Assumption 1), excluding many heavy-tailed or discontinuous noise distributions.  \n  - The data-generating process is assumed to be strictly stationary and $\\beta$-mixing with known decay, which is difficult to assess empirically and may be violated in real-world time series with regime changes or long memory.\n\n  While such assumptions are common in the theoretical literature, their practical relevance and verifiability should be more thoroughly discussed in the paper.\n\n- **No empirical illustration of theoretical effects**: The paper does not include any experiments showing how convergence or generalization degrades as a function of the mixing rate $r$, or comparing performance between i.i.d. and dependent settings.\n\n- **Restrictive stationarity assumption**:  \n  The theoretical results require the sequence $\\{\\mathbf{z}_t\\}$ to be strictly stationary and $\\beta$-mixing. While this is a more general dependence structure than assuming i.i.d. samples, the stationarity assumption may still be too strong for many practical time series settings, including financial and macroeconomic data, where trends, regime shifts, or heteroscedasticity are common. The paper does not examine the robustness of the method to nonstationary violations or provide guidance on verifying this assumption empirically.", "was_revised": true, "final_justification": "### Final Justification\n\nAfter reading the rebuttal and engaging with the authors’ detailed responses, I am satisfied that the main concerns raised in my initial review have been thoughtfully addressed. I maintain my **score of Accept**, based on the following points:\n\n- **Clarification of theoretical assumptions**:  \n  The authors acknowledged that assumptions such as strict stationarity and known $\\beta$-mixing rates are strong, but also standard in related work. They provided additional references and clarified how these assumptions affect the rates and theoretical guarantees.\n\n- **Improved intuition and technical explanation**:  \n  The revised explanation of the $(1 - 1/r)$ rate slowdown and the theoretical role of the blocking parameter $\\ell$ help make the analysis more transparent. The clarification of how $\\ell$ balances the stochastic and dependence errors was particularly helpful.\n\n- **New empirical evidence**:  \n  The authors added experiments illustrating the impact of temporal dependence (via decreasing $r$) on convergence and performance. They also examined failure modes at extreme quantile levels and provided helpful runtime comparisons.\n\n- **Scalability and sensitivity**:  \n  Additional experiments on hyperparameter sensitivity and computational cost across dataset sizes improve the practical understanding of the method.\n\n- **Commitment to revision**:  \n  The authors indicated they will revise the manuscript to incorporate these clarifications and results. These revisions will significantly improve the clarity and completeness of the paper.\n\nWhile some idealized assumptions remain, I find the overall contribution both theoretically solid and empirically validated. The paper represents a valuable advance in distributional regression under temporal dependence and is suitable for publication.", "rebuttal": "Thank you for your valuable feedback and for the careful review of our paper. Below, we address your concerns regarding the weaknesses and provide detailed answers to your questions.\n\n---\n**W1.** **Limited intuition for role of $r$**.\n\n**A1.** Thank you for your insightful observations. The parameter $r$ governs the strength of dependence between observations: larger values of $r$ correspond to weaker dependence, while smaller values indicate stronger temporal dependence. This directly affects the convergence rate of our estimators — sequences with weaker dependence converge more quickly. We will include a detailed discussion in the revision.\n\n---\n**W2.** **Strong theoretical assumptions**.\n\n**A2.** Thank you for this insightful observation. We agree the assumptions are strong, though they are standard in recent theoretical works [1] on quantile regression.\n- **Hierarchical Interaction Model.** We acknowledge that the intrinsic structure is typically unknown in practice. Theoretically, if the target function only satisfies a Hölder condition, the rate depends on the ambient dimension $d$; but if it exhibits intrinsic low-dimensional structure, the rate depends on the intrinsic dimension $\\gamma^\\star\\ll d$.\n- **Noise Density.** This assumption is common in theoretical analyses of quantile regression—spanning linear, kernel-based, and neural models. It can often be relaxed, especially for excess risk bounds, and our framework covers many heavy-tailed distributions (e.g., Cauchy, t, Pareto).\n- **Data-Generating Process.** We agree these assumptions may not always hold in real-world time series. We will clarify this limitation in the revision and discuss implications and potential extensions to locally or piecewise stationary processes.\n---\n**W3.** **No empirical illustration of theoretical effects**.\n\n**A3.** We appreciate this helpful suggestion. The mixing rate $r$ indeed plays a key theoretical role in convergence and generalization. To illustrate its impact, we conducted a new experiment with four nonlinear AR models and four i.i.d. models at quantile level $\\tau=0.5$:  \n- **$AR_i$:**  \n$$X_t = \\alpha_i (X_{t-1} + X_{t-2})+\\epsilon_t,\\ \\epsilon_t \\sim \\mathcal{N}(0,1),\\ X_{-1}, X_{-2}\\sim\\mathcal{N}(0,\\beta_i).$$\n- **$IID_i$:**  \n$$Y_t\\sim\\mathcal{N}(0,\\beta_i).$$  \n\nWe set $\\alpha_i$ and $\\beta_i$ as the $i$‑th elements of $\\{0.1,0.2,0.3,0.4\\}$ and $\\{45/44,10/9,35/26,15/7\\}$. For each $i$, $IID_i$ is the stationary distribution of $AR_i$. By Theorem 1 of [2], each $AR_i$ is exponentially $\\beta$‑mixing with $r=1$, while the i.i.d. cases correspond to $r=\\infty$.  \n\nIn each experiment, we use an SRNN with $W=100$ and $L=3$ as the predictor. As shown in Table 1, stronger dependence clearly slows convergence, reflected in the higher MSE values of AR models versus their i.i.d. counterparts.  \n\n**Table 1: MSE comparison for different $\\beta$‑mixing rates**\n|Models|$\\alpha_i=0.1$|$\\alpha_i=0.2$|$\\alpha_i=0.3$|$\\alpha_i=0.4$|\n|-|-|-|-|-|\n|IID|0.0006|0.0006|0.0008|0.001|\n|AR|0.0015|0.0022|0.0045|0.0046|\n\n---\n**W4.** **Restrictive stationarity assumption**.\n\n**A4.** We thank the reviewer for this thoughtful comment. We agree that strict stationarity and $\\beta$‑mixing, though standard in theory, are strong in practice. Due to space limits, we placed our real‑data evaluation—a GDP forecasting task—in Appendix C. Despite the dataset’s nonstationary nature, our method performs robustly, providing empirical evidence that the estimators remain effective beyond these idealized assumptions. \n\n---\n**Q1.** **Origin and intuition of the $(1-1/r)$ factor slowdown**.\n\n**A5.** Thank you for the thoughtful question. This slowdown naturally arises from the standard blocking technique used to handle dependence. In our analysis, the block length $\\ell$ depends on the mixing rate $r$ and sample size $n$, balancing the error between dependent and independent sequences. When optimizing $W_0L_0$ to trade off approximation and generalization errors, the convergence rate inevitably introduces the $(1-1/r)$ term. We will clarify this derivation in the revision.\n\n---\n**Q2.** **Theoretical role of the blocking parameter $\\ell$**.\n\n**A6.** Thank you for this insightful question. In our analysis, we set $\\ell$ to balance the stochastic error $\\delta_b$ and dependence error $\\beta((\\ell-1)N+1)$, achieving the smallest tail bound in Eqs. (3.3)–(3.4). We take $\\ell \\asymp (\\log n)^{1/r}$ for the exponential case and $\\ell \\asymp n^{1/r}\\log n$ for the algebraic case. A similar dependence on $\\ell$ and $r$ appears in [3], though focused on squared loss. Our choice is sharp in the exponential case since the rate is optimal, but may be conservative for the algebraic case, as the lower bound there remains unknown.\n\n---\n**Q3.** **Estimating the mixing rate $r$ in practice**.\n\n**A7.** Thank you for this thoughtful question. Under exponential mixing, $r$ appears only inside logarithmic terms, so the convergence rate is largely insensitive to its exact value. As a result, not knowing the value of $r$ precisely has minimal impact on performance. In contrast, under algebraic mixing, the rate includes a multiplicative factor of $n^{-(1-1/r)\\gamma^\\star/(2\\gamma^\\star+1)}$, which makes it directly and significantly affected by $r$. Fortunately, several methods for estimating $r$ have been proposed [4–6], and we refer readers to those works for details.\n\n---\n**Q4.** **Effect of nonstationarity**.\n\n**A8.**  Thank you for this valuable question. Our current theoretical analysis is based on the assumption of strict stationarity, which serves as an initial attempt to derive convergence guarantees under temporal dependence. However, we also conducted a nonstationary GDP forecasting experiment in Appendix C, where our method performs robustly. This provides empirical evidence that the estimators remain effective even beyond the idealized assumptions.\n\n---\n**Q5.** **Effect of stronger dependence (smaller r)**.\n\n\n**A9.** We appreciate this helpful suggestion and have added the new experiment as detailed in **A3**.\n\n---\n**Q6.** **Failure cases**.\n\n**A10.** We thank the reviewer for this important comment. We agree that both RNNs and SRNNs perform poorly at extremely low quantile levels. To illustrate this, we conducted additional experiments on Model 2 at five quantile levels, where the noise is generated from a $t_{2.25}$ distribution. As reported in Table 2, although RNN and SRNN outperform FNN and QRF, they still experience a notable performance drop at these extreme quantiles.\n\n**Table 2: MSE comparison under very extreme quantile levels**\n\n|Methods|$\\tau$=0.005|$\\tau$=0.01|$\\tau$=0.02|$\\tau$=0.03|\n|:-:|:-:|:-:|:-:|:-:|\n|QRF|78.4|39.8|29.4|20.8|\n|FNN|94.7|49.7|25.7|16.8|\n|RNN|67.3|34.7|19.9|13.8|\n|SRNN|67.4|34.7|19.8|13.7|\n\n---\n\n**Q7.** **Computation time and scalability.** \n\n**A11.** We thank the reviewer for this suggestion. To compare computational cost, we ran experiments on datasets of different sizes. Table 3 summarizes the results: each column shows the number of training/testing samples, and each cell reports training/testing time (seconds).  \n\nKey observations:\n- **Training time.** FNN is fastest across all sizes. QRF trains faster than RNN and SRNN on small datasets, but slows sharply as data grows. SRNN consistently takes longer than RNN due to pruning overhead.\n- **Testing time.** All NN models evaluate far faster than QRF.\n- **Scalability.** NN models scale efficiently for inference, while QRF’s evaluation time grows steeply with dataset size.\n\n**Table 3: Computation time comparison under different dataset sizes**\n\n|Methods|500/5000|1000/10000|2000/20000|4000/40000|8000/80000|\n|:-:|:-:|:-:|:-:|:-:|:-:|\n|QRF|0.652/1.035|1.344/3.337|3.286/11.226|8.532/51.767|20.869/182.0|\n|FNN|0.284/0.004|0.125/0.003|0.124/0.007|0.271/0.017|0.373/0.036|\n|RNN|1.759/0.007|2.881/0.009|2.909/0.021|1.179/0.029|3.062/0.064|\n|SRNN|2.469/0.013|5.090/0.017|4.691/0.050|1.414/0.048|6.074/0.147|\n\n---\n**Q8.** **Hyperparameter sensitivity**.\n\n**A12.** Thank you for this suggestion. To assess hyperparameter sensitivity, we ran additional experiments varying network depth ($L$), width ($W$), and pruning ratio on Model 1 under the $t_{2.25}$ distrubution. \n\n- **Depth and width (Table 4):**  Prediction accuracy declines slightly as model size grows, then stabilizes once the network is sufficiently wide and deep.\n- **Pruning ratio (Table 5):** Moderate pruning (up to 70%) improves accuracy, with performance peaking at a pruning ratio of 0.7 before slightly degrading.\n\n**Table 4: Accuracy across different $L$ and $W$**\n||$W=16$|$W=32$|$W=64$|$W=128$|$W=256$|$W=512$|$W=1024$|\n|:---:|:---:|:---:|:---:|:---:|:---:|:---:|:---:|\n|$L=2$|0.056|0.054 |0.054 |0.053|0.06|0.065| 0.097|\n|$L=3$|0.055|0.051|0.051 |0.05| 0.054|0.058| 0.073|\n|$L=4$|0.055|0.05|0.052 |0.05| 0.051|0.052| 0.07|\n|$L=5$|0.05|0.05|0.05| 0.052|0.05|0.05|0.057|\n|$L=6$|0.06|0.05|0.053| 0.05 |0.049|0.049|0.056|\n|$L=7$|0.05|0.049|0.049 |0.05|0.049|0.049|0.056|\n\n**Table 5: Effect of pruning ratio on SRNN accuracy**\n|Pruning ratio|MSE|\n|:-:|:-:|\n|0.2 | 0.046|\n|0.3| 0.044|\n|0.4| 0.040|\n|0.5| 0.039|\n|0.6| 0.038|\n|0.7| 0.036|\n|0.8| 0.038|\n\n\n---\n\nWe will revise the paper according to your comments to help readers understand this work better. Thanks again for your insightful review and appreciation of our work!\n\n**References:**\n\n[1] Quantile regression with ReLU networks: Estimators and minimax rates, JMLR 2022\n\n[2] Mixing properties of ARMA processes, Stochastic processes and their applications 1988\n\n[3] Approximation Bounds for Recurrent Neural Networks with Application to Regression, ArXiv 2024\n\n[4] Estimating beta-mixing coefficients, PMLR 2011\n\n[5] Inferring the mixing properties of an ergodic process, ArXiv 2021\n\n[6] Estimating beta-mixing coefficients via histograms, Electron. J. Statist 2015"}, {"review_id": "SDTRv09H1b", "reviewer": "Reviewer_nz3r", "summary": "This paper studies nonparametric quantile regression using recurrent neural networks (RNNs) and sparse recurrent neural networks (SRNNs) to approximate conditional quantile functions with hierarchical interaction structure. The authors demonstrate that RNN and SRNN estimators with ReLU activation achieve optimal nonparametric convergence rates under stationary observations. They derive approximation error bounds for the hierarchical interaction model.", "questions": "- Question 1. How do the convergence rate results translate into practical guidance for model selection, hyperparameter tuning, or confidence interval construction in quantile regression applications?\n\n- Question 2. What methods do the authors propose for hypothesis testing, confidence bands, or uncertainty quantification given their theoretical framework?\n\n- Question  3. Beyond the brief mention in simulations, what detailed empirical evidence demonstrates the practical advantages of RNNs and SRNNs over existing nonparametric quantile regression methods?\n\n- Question  4. How do the authors justify applying their method to GDP forecasting when GDP data typically exhibits non-stationary behavior, violating the stationarity assumptions required by their theory?\n\n- Question  5. Given the disconnect between the abstract's emphasis on GDP forecasting and the main text's purely theoretical focus, how should the authors restructure the paper to better align the promised applications with the delivered content?", "limitations": "Yes", "rating": 2, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 2, "significance": 1, "originality": 2, "strengths_and_weaknesses": "### Strengths \n- Theoretical results are provided. \n\n### Weaknesses\n- The usefulness of theoretical results are unclear. Theoretical results do not provide any statistical inference methods. \n- The practical applications of the proposed method are not examined in the main text. The paper appeals only to a limited audience and would be more suitable for a specialized journal.\n- Please review the literature properly. For instance, one of the important reference is \n   - White, H. (1992). Nonparametric estimation of conditional quantiles using neural networks. In Computing Science and Statistics: Statistics of Many Parameters: Curves, Images, Spatial Models (pp. 190-199). New York, NY: Springer New York.\n - The abstract of the paper mentions an application of GDP forecasting but the main text has no information about it at all. The organization of the paper should be reconsidered. Also, GDP usually do not satisfy the stationary conditions made in the paper and thus the authors might need another application.", "was_revised": true, "final_justification": "The authors kept shifting arguments to avoid answering my questions, which leads to unsatisfactory conclusion. Also, the authors admitted that the paper have little implication for real world applications.", "rebuttal": "We sincerely appreciate your comments. Our theoretical study is motivated by the practical success of RNNs in handling sequence and time series data. Below, we address the weaknesses you identified and provide detailed responses to your specific questions.\n\n---\n**W1.** The usefulness of theoretical results are unclear. Theoretical results do not provide any statistical inference methods.\n\n**A1.** We appreciate the important question. As detailed in **A2.**, our primary focus in this paper has been to establish the fundamental theoretical properties of our proposed methods. Besides, we will include a dedicated discussion on hyperparameter tuning and conformal prediction in the revision, as detailed in **A5**.\n\n---\n**W2.** The practical applications of the proposed method are not examined in the main text. The paper appeals only to a limited audience and would be more suitable for a specialized journal.\n\n**A2.** We thank the reviewer for this comment. The primary focus of our paper is to establish the fundamental theoretical properties of RNNs and SRNNs. We clarify the contributions of our work in terms of techniques and results separately below.\n- **Result.** To the best of our knowledge, this is the first work to establish a minimax‑optimal convergence rate for nonparametric quantile regression using RNNs and SRNNs under stationary observations. \n- **Techniques.** (i) **Novel decomposition.** Our analysis is based on a new donut-shaped decomposition. This decomposition is different from the direct approach used in Eq. (16) of [1], whose analysis crucially relies on the nonnegativity of the squared loss. (ii) **Sharper inequality.**  Compared to the existing results [2], we develop a novel and sharper empirical process inequality (Lemma 7) that underpins the tightness of our oracle inequality. \n\nWe acknowledge that practical guidance is a highly important aspect for real-world applications. However, a detailed treatment of these topics was not the focus of this foundational theoretical work. We note that many works focusing primarily on theory—sometimes with minimal or no experimental components—have been published in top-tier machine learning conferences and journals, including NeurIPS[3], ICLR[4], and JMLR[5]. These examples demonstrate that theoretical contributions of independent significance are valued by the NeurIPS community. We believe that NeurIPS values fundamental theoretical advancements, and our work provides significant contributions to the understanding of nonparametric quantile regression with neural networks, which is a topic of broad interest in the machine learning community.\n\n---\n**W3.** Please review the literature properly. \n\n**A3.** We appreciate the reference and will cite and discuss White, H. (1992) in the revision. \n\n---\n**W4.** The abstract of the paper mentions an application of GDP forecasting but the main text has no information about it at all. The organization of the paper should be reconsidered. Also, GDP usually does not satisfy the stationary conditions made in the paper and thus the authors might need another application.\n\n**A4.** Thank you for this insightful comment. As noted in **A8**, we will reorganize the paper accordingly and, as discussed in **A7**, include an additional experiment that satisfies the stationarity conditions.\n\n---\n**Q1&2** How do the convergence rate results translate into practical guidance for model selection, hyperparameter tuning, or confidence interval construction in quantile regression applications? What methods do the authors propose for hypothesis testing, confidence bands, or uncertainty quantification given their theoretical framework?\n\n**A5.** We appreciate these insightful questions. The primary goal of our paper is to establish the core theoretical properties of RNNs and SRNNs, as discussed in **A2**.\n\nWe also recognize the practical importance of these considerations. In practice, once a solid theoretical foundation is set, researchers naturally move on to issues like tuning and uncertainty quantification. While a full exploration of these topics might extend beyond the typical scope of a conference paper, we are committed to addressing this valuable feedback. \n\n**Hyperparameter tuning.**  As formalized in Theorems 4 and 5, our theoretical results provide principled guidance for selecting the width and depth of the neural networks—as well as the sparsity level for SRNNs—such that the resulting estimators achieve minimax-optimal convergence rates.\n\n**Conformal prediction.** We introduce a novel split conformal RNN quantile regression algorithm leveraging our proposed RNN-based methodology. Building upon the principles outlined in [6], Algorithm 1 constructs a valid prediction interval for $y_{n+1}$ given new input data $x_{n+1}$. Additionally, we plan to conduct a new experiment to evaluate the potential impact of $\\beta$-mixing on our procedure.\n\n**Algorithm 1: Split Conformal RNN Quantile Regression.**\n\n**Input:**\n* Sequentially stationary observations $\\\\{(x_t, y_t)\\\\}_{t=1}^n$, where any consecutive $N$ observations follow the same distribution.\n* Miscoverage level $\\alpha \\in (0,1)$, $\\alpha_{lo}=\\alpha/2$ and $\\alpha_{hi}=1-\\alpha/2$.\n* Algorithm $\\mathcal{A}$: Quantile regression with RNN.\n\n**Process:**\n1.  Form new dataset $\\\\{(X\\_t, y\\_{t+N-1})\\\\}\\_{t=1}^{n-N+1}$, where $X\\_t=(x\\_t,\\ldots,x\\_{t+N-1})\\in \\mathbb{R}^{d\\_x\\times N}$, for $t\\in[n-N+1]$.\n2.  Choose $n_0\\in[n-N+1]$ and construct two disjoint sets $\\mathcal{I}_1=[n_0]$ and $\\mathcal{I}_2=[n-N+1]/[n_0]$.\n3.  Fit two conditional quantile functions: $\\\\{\\hat{q}\\_{\\alpha\\_{lo}}, \\hat{q}\\_{\\alpha\\_{hi}}\\\\} \\leftarrow \\mathcal{A}(\\\\{(X\\_i, y\\_i) : i \\in \\mathcal{I}\\_1\\\\})$.\n4.  Compute $E\\_i:=\\max \\\\{\\hat{q}\\_{\\alpha\\_{lo}}(X\\_i)-y\\_i, y\\_i-\\hat{q}\\_{\\alpha\\_{hi}}(X\\_i)\\\\}$ for $i \\in \\mathcal{I}\\_2$.\n5.  Compute $Q_{1-\\alpha}(E, \\mathcal{I}_2)$, the $(1-\\alpha)(1+1/|\\mathcal{I}_2|)$-th empirical quantile of $\\\\{E_i : i \\in \\mathcal{I}_2\\\\}$.\n\n**Output:**\n- Prediction interval $C(x\\_{n+1}) = [\\hat{q}\\_{\\alpha\\_{lo}}(X\\_{n-N+2}) - Q\\_{1-\\alpha}(E, \\mathcal{I}\\_2), \\hat{q}\\_{\\alpha\\_{hi}}(X\\_{n-N+2}) + Q\\_{1-\\alpha}(E, \\mathcal{I}\\_2)]$ for unseen input $x\\_{n+1}$, where $X\\_{n-N+2}=(x\\_{n-N+2},\\ldots,x\\_{n+1})$.\n\n---\n**Q3.** Beyond the brief mention in simulations, what detailed empirical evidence demonstrates the practical advantages of RNNs and SRNNs over existing nonparametric quantile regression methods?\n\n**A6.** Thank you for this insightful question. For the stationary setting, we conducted an additional experiment satisfying the stationarity conditions (see **A7**) and found that both RNN‑ and SRNN‑based estimators consistently outperform FNN and QRF methods. For the non‑stationary setting, we further verify the robustness of RNN‑ and SRNN‑based estimators in Appendix C.\n\n---\n**Q4.** How do the authors justify applying their method to GDP forecasting when GDP data typically exhibits non-stationary behavior, violating the stationarity assumptions required by their theory?\n\n**A7.** Thanks for the important question. As noted, GDP data exhibits non-stationary behavior, yet our empirical results show that our method remains robust, providing strong evidence that the proposed estimators are effective even beyond the idealized stationarity assumptions of our theoretical framework.\n\nTo further support our theoretical results under stationarity, we conducted an additional experiment using the Dow Jones Industrial Average dataset from January 1, 2000, to December 31, 2020. Since the log returns of stock prices are widely regarded as approximately stationary [7-8], this dataset is well-suited for evaluation.\n\nWe compare the out-of-sample predictive performance of four estimators for 30-business-day-ahead forecasting, using data from January 1, 2000, to December 31, 2019, for training, and the remaining data for testing. All models follow the specifications in Section 4.\n\nTable 1 presents the mean empirical check loss at five quantile levels. The results show that the RNN and SRNN estimators consistently outperform FNN and QRF methods, with the SRNN achieving predictive accuracy similar to the RNN while also inducing sparsity.\n\n\n\n**Table 1: Out-of-sample prediction errors at different quantiles.**\n|Models|$\\tau=0.1$|$\\tau=0.25$|$\\tau=0.5$|$\\tau=0.75$|$\\tau=0.9$|\n|:-:|:-:|:-:|:-:|:-:|:-:|\n|QRF|0.456|0.698|0.817|0.616|0.365|\n|FNN|0.538|0.735|0.810|0.662|0.404|\n|RNN|0.410|**0.640**|0.760|0.562|0.306|\n|SRNN|**0.406**|0.647|**0.759**|**0.561**|**0.305**|\n\n---\n**Q5.** Given the disconnect between the abstract's emphasis on GDP forecasting and the main text's purely theoretical focus, how should the authors restructure the paper to better align the promised applications with the delivered content?\n\n**A8.** We appreciate this important point. We will revise the paper's structure by adding the new experiment on stationary observations as discussed in **A7** into the main text.\n\n---\nIf our responses have satisfactorily addressed your concerns, we would be grateful if you would consider updating your score. Thank you again for your time and thoughtful review.\n\n\n**References:**\n\n[1] Approximation Bounds for Recurrent Neural Networks with Application to Regression, ArXiv 2024\n\n[2] Deep quantile regression: Mitigating the curse of dimensionality through composition, ArXiv 2021\n\n[3] Improved particle approximation error for mean field neural networks, NeurIPS 2024\n\n[4] Generalization bounds and model complexity for kolmogorov–arnold networks, ICLR 2025\n\n[5] Quantile regression with ReLU networks: Estimators and minimax rates, JMLR 2022\n\n[6] Conformalized Quantile Regression, NeurIPS 2019\n\n[7] Value at risk forecasting with the ARMA-GARCH family of models in times of increased volatility, IES Working Paper 2011\n\n[8] Does the Baltic Dry Index predict economic activity in South Africa? A review from 1985 to 2016, South African Journal of Economic and Management Sciences 2018"}, {"review_id": "XQRN13LbcO", "reviewer": "Reviewer_USe5", "summary": "This (mainly theoretical) contribution studies quantile regression using recurrent neural networks (RNNs). Both fully-connected and sparse RNNs are covered. Generalisation bounds for the empirical risk minimisers are derived under a compositional assumption on the target function. The paper shows advantages of RNNs (sparse and non-sparse) in a simulation study.", "questions": "- provide a more complete literature review (for instance ArXiv:2009.08848 is not cited and discussed)\n- work out a theoretical explanation why RNNs are better than FNNs. \n- compare the methods on real data", "limitations": "Yes", "rating": 5, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strength:\n- The paper combines several areas in the literature, including quantile regression, recent approximation results for RNNs, and hierarchical function classes. \n- The mathematics is solid, the paper is well-written, and although this is mainly a theoretical contribution, the paper is accessible to a wider audience. \n\n\nWeaknesses:\n- There is a close connection between RNNs and FNNs (Proposition 2 in Jiao et al. 24, mentioned on l.78 of the article). This means that approximation theory just carries over. I do not see in the theory why RNNs are really necessary here and why not simply taking an FNN (in the simulations a big difference can be observed)\n- Although sparse networks are mentioned as a main contribution, it seems rather straightforward\n  \n\n- typo in l.107: $\\mathbf{v}_L^{(1)},\\ldots, \\mathbf{v}_L^{(N)}$ not defined\n- can you change $\\log n_{\\max}/n_{\\min}$ to $\\log (n_{\\max})/n_{\\min}$ everywhere?", "was_revised": true, "final_justification": "The provided answers to the referees are very honest and clear. My understanding is that they mainly confirm my earlier assessment. Therefore, I will stick to the original scores. For score 6, one would need a compelling theoretical reason showing the superiority of RNNs compared to FNNs as observed in the simulations.", "rebuttal": "We sincerely appreciate your thoughtful feedback. Below, we respond to the weaknesses you identified and provide detailed answers to your specific questions.\n\n\n---\n**W1.** There is a close connection between RNNs and FNNs (Proposition 2 in Jiao et al. 24, mentioned on l.78 of the article). This means that approximation theory just carries over. I do not see in the theory why RNNs are really necessary here and why not simply taking an FNN (in the simulations a big difference can be observed).\n\n**A1.** We thank the reviewer for this insightful comment. Theoretically, if constant factors are disregarded, we agree that the approximation theory for RNNs and FNNs can appear similar. However, our simulations clearly demonstrate the practical advantages of RNNs. While theory may not capture all nuances, RNNs' recurrent structure is inherently better suited for sequential data and capturing temporal dependencies than FNNs. Our empirical results show that this architectural difference leads to substantial performance gains, particularly in time series forecasting. Thus, despite theoretical similarities in approximation bounds, the practical benefits of RNNs, as evidenced by our simulations, justify their necessity for superior performance in sequential data applications.\n\n---\n**W2.** Although sparse networks are mentioned as a main contribution, it seems rather straightforward.\n\n**A2.** We thank the reviewer for this comment. While we acknowledge that sparse network theories build on non-sparse ones, our contribution to sparse networks in the context of nonparametric quantile is both novel and nontrivial. In particular, Lemmas 4 and 14 in the supplementary material provide a nuanced characterization of the equivalence between SRNNs and SFNNs, grounded in a detailed analysis of neural network sparsity. To the best of our knowledge, this equivalence has not been established in prior work. Moreover, it is more general and nuanced than Proposition 2 in [1]. \n\nFurthermore, our central technical innovation lies in establishing oracle inequalities (Theorem 3), including generalization error bounds, for quantile regression under β-mixing dependence—a setting not previously addressed in the literature. To better situate these contributions, we present a concise overview of our analytical approach, highlighting the key innovations as they arise at each step.\n\n- **Novel decomposition.** We begin by introducing a donut-shaped decomposition which has not appeared in the context of nonparametric regression with neural networks. Specifically, we introduce donut-shaped sets (Line 514), which allow us to decompose the probability bound $\\mathbb{P}(\\|\\hat{f}-f_0\\|_2>\\delta\\_\\star)$ into manageable components by bounding each term  $\\mathbb{P}(\\hat{f}\\in\\mathcal{D}_i)$ separately. This decomposition is different from the direct approach used in Eq. (16) of [1], whose analysis crucially relies on the nonnegativity of the squared loss—namely, $$\\mathbf{E}[(f(x_t)-y_t)^2-(f_0(x_t)-y_t)^2]=\\mathbf{E}[(f(x_t)-f_0(x_t))^2]\\ge0,$$ where $f_0$ denotes the true function and $y_t=f_0(x_t)+\\epsilon_t$. The check loss function, however, does not satisfy this property. Consequently, we develop the novel decomposition introduced here.\n\n- **Refined blocking technique.**  Next, we introduce a refined blocking technique to relate the mixing sequence to its i.i.d. counterparts, which differs from the approach in [1]. In Step 2 of their proof of Theorem 13, part of the data within each partition is discarded to facilitate analysis under the squared loss, leveraging the nonnegativity property noted above. However, the check loss function used in quantile regression does not possess this property. Consequently, we design a new partitioning procedure (Figure 2) that retains all observations, enabling a sharper and more comprehensive analysis under dependent data. Moreover, while Lemma 16 in [1] plays a central role in their argument, it cannot be directly applied to our setting. Instead, we establish and employ a probabilistic counterpart (Lemma 5) to carry out our analysis.\n\n- **Sharper inequality.** Finally, we develop a novel and sharp empirical process inequality (Lemma 7) that underpins the tightness of our oracle inequality. One reason the result of [2] lacks tightness is that their analysis relies on a non‑sharp application of the Bernstein inequality. To address this limitation, we derive Lemma 7 by building on Theorem 7.3 of [3] and Corollary 5.1 of [4].\n\nWe will clarify and highlight these novel aspects more explicitly in the revision.\n\n---\n**W3.** typo in l.107: $\\mathbf{v}_L^{(1)}, \\dots,\\mathbf{v}_L^{(N)}$  not defined\n\n**A3.** Thank you for pointing it out. We will make the corrections in the revision.\n\n---\n**W4.** can you change $\\log n_{\\max}/n_{\\min}$ to $\\log (n_{\\max})/n_{\\min}$  everywhere?\n\n**A4.** Thank you for the advice. We will make the suggested changes in the revision.\n\n---\n**Q1.** provide a more complete literature review (for instance ArXiv:2009.08848 is not cited and discussed).\n\n**A5.** Thank you for pointing this out. We will include and discuss [arXiv:2009.08848]. We will also expand the literature review accordingly.\n\n---\n**Q2.** work out a theoretical explanation why RNNs are better than FNNs.\n\n**A6.** We thank the reviewer for this excellent and insightful point regarding the theoretical distinctions between RNNs and FNNs. We will provide a more detailed discussion on this in the revised manuscript.\nAs the reviewer correctly points out, our theoretical framework, particularly concerning approximation capabilities, is based on the equivalence between RNNs(SRNNs) and FNNs(SFNNs) in certain contexts (Lemma 4 and Lemma 14 of the supplementary material). Consequently, from a purely theoretical perspective, the performance of both architectures might appear similar. However, as we indeed observed in our simulations and real-data experiments, RNNs consistently demonstrate superior performance in practice. This practical advantage stems from the specialized recurrent structure of RNNs, which allows them to effectively capture and model temporal dependencies and sequential patterns inherent in time series data. \n\nA rigorous theoretical proof remains an open problem, and we plan to continue exploring it in future work. Meanwhile, we offer a simple example illustrating why RNNs can be more expressive than FNNs. Consider an RNN with $d_x=d_y=L=1$ and $W=N=2$. In this configuration, the RNN can exactly compute the function  \n$$f(x_1,x_2) =\\max(x_2+\\max(x_1,0),0)+\\max(x_2-1+\\max(x_1-1,0),0).$$By contrast, any FNN with $d_x=W=2$ and $L=d_y=1$ cannot represent $f$.  \n\nThe obstruction arises from the number of non-differentiable points in the induced function. Fix $x_2=-9$; then the univariate function $f(x_1,-9)$ exhibits three non-differentiable points. However, any such FNN can induce, for each fixed $x_2$, at most two breakpoints in $x_1$.\n\n---\n**Q3.** compare the methods on real data.\n\n**A7.** We thank the reviewer for this valuable suggestion. We agree that demonstrating the practical advantages of our methods on real-world data is crucial. In our current submission, we have already included an empirical evaluation on real data: a GDP forecasting task. Due to space constraints, this experiment was placed in Appendix C. In this evaluation, we compare our proposed methods against existing nonparametric quantile regression approaches, specifically QRF and FNN. The results from this experiment already demonstrate the superior performance of our approach compared to these baselines.\nTo further strengthen our empirical evidence, we conducted an additional experiment using the Dow Jones Industrial Average dataset covering the period from January 1, 2000 to December 31, 2020. Consistent with prior studies, the log returns of stock prices are widely regarded as approximately stationary, making this dataset well-suited for our evaluation [5–6]. \n\n\nWe compare the out‑of‑sample predictive performance of RNN‑, SRNN‑, FNN‑, and QRF‑based estimators for 30‑business‑day‑ahead forecasting. The training set covers the period from January 1, 2000 to December 31, 2019, with the remaining data reserved for testing. All model implementations follow the specifications provided in Section 4.\n\n\nTable 1 presents the mean empirical check loss for each model at five different quantile levels. The results indicate that the RNN- and SRNN-based estimators consistently outperform both the FNN and QRF methods across most quantile levels. Furthermore, the SRNN achieves predictive accuracy comparable to that of the standard RNN while inducing sparsity.\n\n**Table 1: Out-of-sample prediction errors at different quantiles.**\n|Models| $\\tau=0.1$ | $\\tau=0.25$|$\\tau=0.5$|$\\tau=0.75$|$\\tau=0.9$|\n|:-:|:-:|:-:|:-: |:-:|:-:|\n|QRF|0.456|0.698|0.817|0.616|0.365|\n|FNN|0.538|0.735|0.810|0.662|0.404|\n|RNN|0.410|**0.640**|0.760|0.562|0.306|\n|SRNN|**0.406**|0.647|**0.759**|**0.561**|**0.305**|\n\n---\nWe are grateful for your meticulous review and will carefully revise the identified issues and typos to improve the clarity and readability of the paper.\n\n\n**References:**\n\n[1] Approximation Bounds for Recurrent Neural Networks with Application to Regression, ArXiv 2024\n\n[2] Deep quantile regression: Mitigating the curse of dimensionality through composition, ArXiv 2021\n\n[3] Concentration inequalities for sub-additive functions using the entropy method, In Stochastic Inequalities and Applications 2003\n\n[4] Gaussian approximation of suprema of empirical processes, The Annals of Statistics 2014\n\n[5] Value at risk forecasting with the ARMA-GARCH family of models in times of increased volatility, IES Working Paper 2011\n\n[6] Does the Baltic Dry Index predict economic activity in South Africa? A review from 1985 to 2016, South African Journal of Economic and Management Sciences 2018"}]

Author rebuttal:

Meta-review: This paper studies nonparametric quantile regression using recurrent neural networks (RNNs) and sparse RNNs (SRNNs) to approximate conditional quantile functions with hierarchical compositional structure. The authors establish convergence rates under mixing conditions and derive approximation bounds using connections between sparse feedforward and recurrent architectures. The theoretical results are supported by simulations and an empirical application to GDP forecasting.

While the reviewer team did not reach full consensus, I find the overall contribution to be sound and well-motivated. The theoretical development appears technically solid, and the approach adds value to the literature by bridging quantile regression and recurrent neural architectures under a principled framework.

I recommend acceptance. It would be beneficial for the authors to incorporate feedback from the reviewers and discussion period into the camera-ready version to further improve clarity and positioning.
Decision: Accept (poster)