Optimal estimation of variance in nonparametric regression with random design
Yandi Shen Chao Gao Daniela Witten Fang Han
englishbabel


=4

























    










lemmaLemma

propositionProposition
thmTheorem
definitionDefinition
corollaryCorollary
assumptionAssumption
remarkRemark




thmbis[1]
  
   thm
  
  
propbis[1]
  
   proposition
  

propbisbis[1]
  
   proposition
  


  
lemmabis[1]
  
   lemma
  
  
corollarybis[1]
  
   corollary
  
  
assumptionbis[1]
  
   assumption
  
  
assumptionbisbis[1]
  
   assumption
  




eq(<ref>)
chapChapter <ref>
secSection <ref>
algoAlgorithm <ref>
figFig. <ref>
tabTable <ref>
rmkRemark <ref>
clmClaim <ref>
defDefinition <ref>
corCorollary <ref>
lmmLemma <ref>
lemmaLemma <ref>
propProposition <ref>
appAppendix <ref>
exExample <ref>
exerExercise <ref>
solnSolution <ref>
condCondition <ref>





ℙ
1
 
 
 
 
#1#1
#1#1
#1 
#1 

#1
#1#2#1#2 



∑
∫
/


















































ℚ



























 







||
{}
()
[]
⌈⌉




%̂ŝ%̅s̅%̃s̃












apalike


















 


1]Yandi Shen
2]Chao Gao
3]Daniela Witten
4]Fang Han
[1]
University of Washington
 

ydshen@uw.edu

[2]
University of Chicago
 

chaogao@galton.uchicago.edu

[3]
University of Washington
 

dwitten@uw.edu

[4]
University of Washington
 

fanghan@uw.edu





Optimal estimation of variance in nonparametric regression with random design
    [
    
=============================================================================




Consider the heteroscedastic nonparametric regression model with random design

    Y_i = f(X_i) + V^1/2(X_i)ε_i,    i=1,2,…,n,

with f(·) and V(·) α- and β-Hölder smooth, respectively. 
We show that the minimax rate of estimating V(·) under both local and global squared risks is of the order

    n^-8αβ/4αβ + 2α + β∨ n^-2β/2β+1,

where a∨ bmax{a,b} for any two real numbers a,b. This result extends the fixed design rate n^-4α∨ n^-2β/(2β+1) derived in <cit.> in a non-trivial manner, as indicated by the appearances of both α and β in the first term. 

In the special case of constant variance, we show that the minimax rate is n^-8α/(4α+1)∨ n^-1 for variance estimation, which further implies the same rate for quadratic functional estimation and thus unifies the minimax rate under the nonparametric regression model with those under the density model and the white noise model. To achieve the minimax rate, we develop a U-statistic-based local polynomial estimator and a lower bound that is constructed over a specified distribution family of randomness designed for both ε_i and X_i.



Keywords: variance estimation, nonparametric regression, random design, minimax rate, U-statistics. 




§ INTRODUCTION




Consider the model 

    Y_i = f(X_i) + V^1/2(X_i)ε_i,   i=1,2,…,n,

where {X_i}_i=1^n are independent and identically distributed (i.i.d.) univariate random design points, and {ε_i}_i=1^n are i.i.d. with zero mean, unit variance, and are independent of {X_i}_i=1^n. In this paper, we study the optimal estimation of V(·) under both local and global squared risks. 






 Variance estimation is a fundamental statistical problem <cit.> with wide applications. It is useful in, for example, construction of confidence bands for the mean function, estimation of the signal-to-noise ratio <cit.>, and selection of the optimal kernel bandwidth <cit.>.
 
When {X_i}_i=1^n are fixed, estimation of V(·) in (<ref>) has been studied extensively in the literature via residual-based methods <cit.> and difference-based methods <cit.>. One important heuristic from previous studies is that, compared to residual-based methods, difference-based methods are able to achiever a smaller bias and subsequently a smaller mean squared error by avoiding direct estimation of the mean function. More precisely, when X_i = i/n, i=1,…,n and f(·) and V(·) in (<ref>) are α- and β-Hölder smooth, respectively, <cit.> proposed a difference estimator which achieved the optimal rate of the order n^-4α∨ n^-2β/2β + 1 under both local and global squared risks. 

In contrast, our study focuses on the case where {X_i}_i=1^n are i.i.d. random design points on the real line. For this, we show that when f(·) and V(·) in (<ref>) are α- and β-Hölder smooth, respectively, the minimax rate of estimating V(·) is of the order n^-8αβ/4αβ + 2α + β∨ n^-2β/2β+1 under both local and global squared risks. This result has several noteworthy implications:


  * The minimax rates in random and fixed design settings share a common component, n^-2β/2β+1, as well as the same transition boundary α = β/(4β + 2).

  * For α < β/(4β + 2), a faster rate is achievable with a random design.

  * Unlike the fixed design setting, for α < β/(4β + 2), α and β are now both present in the first term of the minimax rate in the random design case.


































We now discuss in more detail this minimax rate. 
The upper bound of the minimax rate is achieved by smoothing pairwise differences via local polynomial regression, the former of which is formulated via U-statistics. Our analysis of this estimator hence relies on the four-term Bernstein inequality in <cit.>, and unlike classic kernel methods, requires no smoothness assumption on the design density. 

For the lower bound, due to the appearances of both α and β in the non-trivial n^-8αβ/4αβ + 2α + β part of the minimax rate and the additional randomness of {X_i}_i=1^n, the derivation is much more involved than its counterpart in the fixed design setting. We tackle the first difficulty of entangled α and β via a proper localization technique in the construction of the mean function f(·), depicted in Figure <ref> in Section <ref>. The second difficulty caused by the randomness of {X_i}_i=1^n is resolved with a new trapezoid-shaped construction of the mean f(·), aided by a result due to <cit.> on the sparse multinomial distribution. This result helps characterize the asymptotic behavior of the locations of {X_i}_i=1^n and plays a key role in our proof, but to our knowledge has not been well used in the nonparametric statistics literature. 

In the special case of constant variance, (<ref>) is reduced to


    Y_i = f(X_i) + σε_i,    i=1,2,…,n,


and the goal becomes estimation of σ^2. 





In this case, the problem is linked to estimation of a quadratic functional, which has been studied in depth in the other two benchmark nonparametric models, the density model <cit.> and the white noise model <cit.>. In the density model, one observes an i.i.d. univariate sequence {X_i}_i=1^n from some unknown density f(·), and the goal is to estimate ∫ f^2(x)dx. In the white noise model, one observes a continuous-time process from dY_t = f(t)dt + n^-1/2 dW_t for t∈[0,1] with W_t  a standard Wiener process. The goal is to estimate ∫_0^1 f^2(t)dt. Under an α-smoothness condition on f(·), the minimax rate in both of the aforementioned two cases is n^-8α/(4α+1)∨ n^-1 (cf. Theorem 1(ii) and 2(ii) in <cit.>, Theorem 4 in <cit.>). 

Following <cit.>, a quadratic functional of interest under (<ref>) with random design is

    Q ∫ f^2(x)p_X(x)w(x)dx,

where p_X(·) is the unknown design density and w(·)≥ 0 is some known weight function. Assuming in (<ref>) that f is α-Hölder smooth, we show that 












the minimax rate of estimating σ^2 and Q (when σ^2 is unknown) is n^-8α/(4α+1)∨ n^-1, thereby unifying the minimax rate of quadratic functional estimation in all three benchmark nonparametric models. 









In this paper, we also provide extensions of (<ref>) to multivariate cases, with a focus on the multivariate nonparametric regression model




    Y_i    = f(X_i) + σε_i,    i=1,2,…, n,

and the nonparametric additive model

    Y_i    = ∑_k=1^d f_k(X_i,k) + σε_i,    i=1,2,…,n,

in both fixed and random designs. Here, X_i  (X_i,1,…,X_i,d)^⊤, i=1,…,n, for some fixed positive integer d. Regarding the fixed design, we consider two types, namely, the grid design (GD) and the diagonal design (DD). With a total of n design points, the former places them on a regular grid in the d-dimensional cube [0,1]^d while the latter only places design points on the diagonal. Details are given in Sections <ref> and <ref>. 


























































We summarize the minimax rates in all of the aforementioned models in Table <ref>. 

The rest of the paper is organized as follows. Section <ref> presents the simple model (<ref>) with constant variance. Section <ref> discusses its heteroscedastic extension   (<ref>). Section <ref> discusses the multivariate nonparametric regression model (<ref>), the additive model (<ref>), and several other extensions of our main results. The essential lower bound proof of the minimax rate n^-8α/(4α+1)∨ n^-1 under model (<ref>) is presented in Section <ref>, with the rest of the proofs given in a supplement.

The notation used throughout the paper is as follows. For any positive integer n, [n] denotes the set {1,2,…,n}. For any real number a, we use a to denote the smallest integer greater than or equal to a, and a the largest integer strictly smaller than a. For any positive integer d, 0_d denotes the zero vector of dimension d and 𝐈_d denotes the identity matrix of dimension d. For a real vector x, x and x_∞ denote its Euclidean and infinity norms, respectively. For a real matrix 𝐀, we use 𝐀, 𝐀_F, and |𝐀| to denote its spectral norm, Frobenius norm, and determinant, respectively. For an m-times differentiable function f:→ with some positive integer m, we use f^(k) to denote its kth derivative for k=1,2,…,m. For identically distributed random variables X_i and X_j, we use _X_i(·) and p_X_i(·) to denote the distribution and density of X_i,  to denote X_i - X_j, and p_(·) to denote the density of X_i - X_j. Similar notation _X_i(·), p_X_i(·), , p_(·) applies to identically distributed random vectors X_i and X_j. For a positive integer d and μ∈^d,Σ∈^d× d, N_d(μ,Σ) stands for the d-dimensional normal distribution with mean μ and covariance Σ. We will drop the subscript d for simplicity when d=1. Φ(·) and φ(·) represent the standard normal distribution and density. More generally, we will write φ_μ,σ^2(·) as the density for the normal distribution with mean μ and variance σ^2. For two probability measures , defined on a common space (Ω,A), (,) denotes their total variation distance, that is, (,) sup_A∈A*(A)-(A). For two real sequences {a_n} and {b_n}, a_n≲ b_n if |a_n|≤ C|b_n| for some positive absolute constant C. We say a_n≍ b_n if a_n≲ b_n and b_n≲ a_n. 
 





































§ HOMOSCEDASTIC CASE



To illustrate some of the main ideas developed in this paper, we begin with a discussion of the elementary univariate homoscedastic nonparametric regression model (<ref>):

    Y_i = f(X_i) + σε_i,    i=1,2,…,n.

Here, {X_i}_i=1^n are i.i.d. copies of a univariate random variable X, f(·) belongs to an α-Hölder class that will be specified soon, and {ε_i}_i=1^n are i.i.d. copies of a variable ε with zero mean and unit variance and are independent of {X_i}_i=1^n. Both the mean function f(·) and the distribution of {X_i}_i=1^n are assumed unknown. 

Model (<ref>) has been extensively studied using residual-based and difference-based methods; see, among many others, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>. A related functional estimation problem has also been studied in semiparametric models <cit.>. Most of the previous studies focus on the case of fixed design, especially the equidistant design with X_i=i/n, i∈[n], for which the minimax rate of estimating σ^2 under an α-Hölder smoothness constraint on f(·) is known to be n^-4α∨ n^-1 (cf. Theorems 1 and 2 in <cit.>).





In detail, let I be a fixed (possibly infinite) interval on the real line. Define the Hölder class Λ_α,I(C_F) on I as follows:

    Λ_α,I(C_F) {f:     for all  x, y∈ I and  k = 0,…, α,
       *f^(k)(x)≤ C_F and *f^(α)(x) - f^(α)(y)≤ C_F|x-y|^α^'},

where α^'α - α. Denote the support of X as (X).


Define 




the joint distribution class P_cv,(X,ε) (where “cv" stands for “constant variance") with the following conditions:


  (a) X satisfies (X)⊂ I.

  (b) X has density p_X(·) and there exists a fixed positive constant C_0 such that 

    sup_x∈ p_X(x) ≤ C_0.






  (c)There exist two fixed constants δ_0 > 0 and c_0 > 0 such that for any  0 < δ < δ_0, there exists a set U_δ⊂[-1,1] such that

    λ(U_δ) ≥ c_0    and   inf_u∈U_δp_(uδ) ≥ c_0,

where λ(·) represents the Lebesgue measure on the real line, and = X_i - X_j.

  (d)ε^4 ≤ C_ε for some fixed positive constant C_ε.



Note that no smoothness condition is placed on the density of X. Condition (c) essentially requires the density p_ to be “dense" around 0, and is strictly weaker than a uniform lower bound of p_ over a fixed neighborhood of 0. It also follows from the following sufficient condition on the marginal density p_X(·) (see Lemma <ref> in the supplement for the justification):


  (c^') X is compactly supported (taken to be [0,1] without loss of generality). There exists some positive constant c_0 and subset S⊂[-1,1] with Lebesgue measure λ(S)≥ 3/4 such that p_X(t)≥ c_0 uniformly over t∈ S.

In particular, (c^') covers the uniform distribution on [0,1] and the distribution of X in the lower bound construction in the proof of Theorem <ref>.


The rest of the section is devoted to proving, for any fixed positive constants C_F and C_σ, the following minimax rate: 

    inf_σ^2sup_f∈Λ_α,I(C_F)sup_σ^2 ≤ C_σsup__(X,ε)∈P_cv,(X,ε)*σ^2 - σ^2^2 ≍ n^-8α/(4α+1)∨ n^-1,

where _(X,ε) denotes the joint distribution of (X,ε), and σ^2 ranges over all estimators of σ^2.













 §.§ Upper bound


The upper bound is achieved by a difference estimator based on U-statistics (with convention 0/ 0 =0):

    σ̂^2 n 2^-1∑_i<jK_h(X_i - X_j)(Y_i-Y_j)^2/2/n 2^-1∑_i<jK_h(X_i-X_j).

Here, K_h(·)  K(·/h)/h, where h = h_n is a bandwidth parameter satisfying h_n↓ 0 as n→∞, and K(·) is a symmetric density kernel supported on [-1,1] that satisfies 

    M_K ≤inf_|u|≤ 1K(u) ≤sup_|u|≤ 1 K(u) ≤M_K

for two fixed constants M_K and M_K; one example is the box kernel K(u) = 1*|u|≤ 1/2 which satisfies (<ref>) with M_K = M_K = 1/2. 

The following error bound is derived via the exponential inequality for degenerate U-statistics due to <cit.>.



























































Suppose the kernel K(·) in σ^2 is chosen such that (<ref>) is satisfied with constants M_K and M_K, and the bandwidth h_n is chosen as 

    h_n≍
    n^-2/(4α+1),     0 < α < 1/4,
    
    n^-1,    α≥ 1/4.

Then, under (<ref>) with random design, it holds that

    sup_f∈Λ_α,I(C_F)sup_σ^2 ≤ C_σsup__(X,ε)∈P_cv,(X,ε)*σ̂^2 - σ^2^2≤ C*n^-8α/(4α + 1)∨ n^-1,

where C is some fixed positive constant that only depends on M_K, M_K, α,C_F,C_σ and C_0,c_0,C_ε in P_cv,(X,ε).

















The error rate in Theorem <ref> is achieved by choosing the optimal bandwidth h_n to balance the “bias-variance" decomposition:

    **σ̂^2 - σ^2^2^1/2≲ h_n^2(α∧ 1) + 1/nh_n^1/2,
 
where a∧ bmin{a,b} for any two real numbers a,b. The bias term h_n^2(α∧ 1) reflects the second-order effect of the unknown mean on variance estimation, which has been noted by <cit.> and <cit.>. The variance part follows from the fact that there is an average number of n^2h_n pairs of (i,j) such that |X_i - X_j|≤ h_n. We note that the same “bias-variance" decomposition has appeared in quadratic functional estimation in the density model and Gaussian sequence model <cit.>. See Section <ref> for a more detailed discussion.     




While most of the previous works are in the context of fixed design, <cit.> considered constant variance estimation with random design, and their estimator (formula (1.4) therein) is almost identical to our σ^2. Under certain assumptions (Assumptions 1 and 2 and (2.4) - (2.7) therein), they show that their estimator is root-n consistent and asymptotically normal. However, as commented in the first paragraph on p. 184 of their paper, their condition (2.7) is only satisfied when the mean function smoothness α is strictly larger than 1/4, and no analysis is provided below this threshold. Our minimax rate n^-8α/(4α+1)∨ n^-1 therefore confirms that α≥ 1/4 is indeed the minimal requirement for any variance estimator to be root-n consistent and we also demonstrate the optimality of σ^2 for 0 < α < 1/4.






Finally, in (<ref>), we have assumed that the smoothness index α is known. If it is unknown, then the variance can be estimated adaptively via Lepski-type methods <cit.>. This is discussed in more detail in Section <ref>.





 §.§ Lower bound


The derivation of the lower bound in (<ref>) is much more involved. In particular, the construction in the fixed design setting (cf. Theorem 2 in <cit.>) cannot be extended to the random design case, since the spike-type construction of f(·) located at each deterministic design point leads to a sub-optimal rate in the random design setting. To achieve a sharp rate, we have to exploit the randomness of {X_i}_i=1^n; this requires us to handle a highly convoluted alternative hypothesis  that no longer leads to a product measure of {Y_i}_i=1^n given each realization of {X_i}_i=1^n in LeCam's two-point method. This calls for a careful analysis of the locations of {X_i}_i=1^n. 



We now sketch a proof of the n^-8α/(4α+1) component in (<ref>) for 0 < α < 1/4, with a particular emphasis on where the difference arises with the fixed design setting. 





The proof can be roughly divided into two steps. In the first step, we construct a two-point testing problem with the null being a Gaussian (H_0) and the alternative a Gaussian location mixture (H_1). In the second step, we approximate the Gaussian location mixture (H_1) by a location mixture with compact support (H_1), which, unlike the alternative in the first step, belongs to the considered model class. 

We start by introducing the construction of f(·), σ^2, ε, and X under the null H_0 and the alternative H_1 in the first step. For each n, let

    h_n ≍ n^-2/(4α + 1),   θ_n^2 ≍ h_n^2α,   and    N 1/(6h_n),

and divide the unit interval [0,1] into N intervals of length 6h_n, with n large enough and h_n chosen such that N is a positive integer. 


  Choice of f(·): Under H_0, let f≡ 0. Under H_1, let f(·) be a piecewise trapezoidal function on the N intervals. That is, for each i∈[N], f takes on a value of h_n^αr_i on the intervals [(6i-5)h_n, (6i-1)h_n] and then linearly decreases to zero on the two endpoints 6(i-1)h_n and 6ih_n, with {r_i}_i=1^N i.i.d. standard normal variables.

  Choice of σ^2: Under H_0, let σ^2 = 1+θ_n^2. Under H_1, let σ^2 = 1.

  Choice of ε: Under both H_0 and H_1, let ε∼N(0,1).

  Choice of X: Under both H_0 and H_1, let {X_i}_i=1^n be uniformly distributed over the union of the upper bases of the trapezoids, that is, over ⋃_i=1^N[(6i-5)h_n, (6i-1)h_n].

See Figure <ref> for an illustration of the construction.




 In contrast to the spike-type construction of f(·) in the fixed design setting, our construction is trapezoid-shaped, 
which guarantees a maximal variation in the mean to compensate for the difference in the variance under the null and alternative. 
This is unnecessary in the fixed design setting since the point of maximal variation in the mean (center of each spike) can be directly placed at each fixed X_i = i/n, resulting in n evenly spaced spikes in f(·). 


Denote the joint distribution of {(X_i, Y_i)}_i=1^n under H_0 and H_1 by _0 and _1 with respective density p_0 and p_1.  Under the above construction, conditional on {X_i}_i=1^n, {Y_i}_i=1^n are distributed as 

    H_0: p_0({Y_i}_i=1^n|{X_i}_i=1^n) = ∏_i=1^nφ_0,1+θ_n^2(Y_i)

and 

    H_1: p_1({Y_i}_i=1^n|{X_i}_i=1^n) = ∏_j=1^N ∫*∏_{i:b_i = j}φ_h_n^α v, 1(Y_i)φ(v)dv,

where {b_i}_i=1^n is the location index sequence of {X_i}_i=1^n defined as 

    b_i := j    if  X_i∈[(6j-5)h_n, (6j-1)h_n],

which characterizes which trapezoid each X_i falls into. Using Lemma <ref> that will be stated in Section <ref>, one can then upper bound 

    (_0, _1)  = (_0({Y_i}_i=1^n|{X_i}_i=1^n), _1({Y_i}_i=1^n|{X_i}_i=1^n)) ≲θ_n^2nh_n^1/2,

which can be made smaller than a sufficiently small constant c by choosing h_n sufficiently small.


The second step of the proof aims to find a sequence of bounded random variables {r_i}_i=1^N to replace the standard normal sequence {r_i}_i=1^N in _1, so that for each realization of {r_i}_i=1^N, the corresponding f(·) in the alternative is α-Hölder smooth with a fixed constant. Then, denoting the distribution of {r_i}_i=1^N as , one wishes to approximate the conditional distribution _1({Y_i}_i=1^n|{X_i}_i=1^n) in H_1 by _1({Y_i}_i=1^n|{X_i}_i=1^n) with density 

    p_1({Y_i}_i=1^n|{X_i}_i=1^n) = ∏_j=1^N ∫*∏_{i:b_i = j}φ_h_n^α v, 1(Y_i)(dv)

in H_1. 
Even with the aid of moment matching techniques already established in the literature, upper bounding (_1,_1) is still nontrivial. 
Specifically, unlike in the fixed design setting, now with high probability the conditional distribution of {Y_i}_i=1^n given {X_i}_i=1^n is no longer a product measure. This is because multiple X_i's could fall into the same trapezoid in the construction of f(·). This can be handled relatively easily in the first step since there we only have to analyze the pairwise correlation of Y_i| X_i and Y_j | X_j depending on whether X_i and X_j fall into the same trapezoid, but it is much less tractable in the second step. More specifically, in order to match moments, we now have to divide the X_i's into groups based on their memberships among the trapezoids, which naturally requires us to

monitor the locations of {X_i}_i=1^n, and in particular the number of X_i's that fall into the same trapezoid. This is possible by observing that the memberships of {X_i}_i=1^n now follow a sparse multinomial distribution (n^2/(4α + 1) bins, n balls) so that a result in <cit.> can be applied.

This allows us to show that with high probability the maximum number of X_i's in each trapezoid is bounded by a fixed constant,

which, along with Lemma <ref> in Section <ref>, allows us to calculate 

    (_1,_1)≲ nθ_n^2p

for p  1 + 1/4α. This indicates that (_1,_1) is smaller than some sufficiently small constant c. Then, by the triangle inequality,

    (_0, _1) ≤(_0, _1) + (_1, _1) ≤ 2c.


Details of the above derivation will be given in Section <ref>. The resulting lower bound is as follows.



Under (<ref>) with random design, it holds that

    inf_σ^2sup_f∈Λ_α,I(C_F)sup_σ^2 ≤ C_σsup__(X,ε)∈P_cv,(X,ε)*σ^2 - σ^2^2 ≥ c*n^-8α/(4α+1)∨ n^-1,

where c is some fixed positive constant that only depends on α,C_F, C_σ and C_0, c_0,C_ε in P_cv,(X,ε), and σ^2 ranges over all estimators of σ^2.




It remains an open problem to prove a lower bound rate that is strictly slower than n^-1 over the sub-class of P_cv,(X,ε) with more regular designs, which includes in particular the uniform design on [0,1]. We conjecture that in this case, n^-8α/(4α+1)∨ n^-1 is still the minimax rate in view of analogous results in quadratic functional estimation <cit.>.






§ HETEROSCEDASTIC CASE



We now study the heteroscedastic model (<ref>),

    Y_i = f(X_i)  + V^1/2(X_i)ε_i ,   i=1,2,…,n,

where {X_i}_i=1^n are i.i.d. copies of X on the real line, f(·) and V(·) are α- and β-Hölder smooth on the fixed (possibly infinite) interval I, respectively, and {ε_i}_i=1^n are i.i.d. copies of ε with zero mean and unit variance and are independent of {X_i}_i=1^n. As in Section <ref>, smoothness indices α and β are assumed known, while f(·),V(·), and the distribution of X are unknown. 
For any estimator V(·), the estimation accuracy is measured both locally via 

    R_1(V, V;x^*) *V(x^*) - V(x^*)^2

at a point x^* in the support of X, (X), and globally via 

    R_2(V, V) ∫*V(x) - V(x)^2_X(dx)

with _X the distribution of X. 


Model (<ref>) has been studied in, for example, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, with a focus mainly on the fixed design case. An exception is <cit.>, with which we draw a detailed comparison in Remark <ref> below. Theorems 1 and 2 in <cit.> established a minimax rate of the order n^-4α∨ n^-2β/(2β+1) under equidistance design X_i = i/n, i∈[n] when f(·) and V(·) are α- and β-Hölder smooth on [0,1]. 





Define P_vf,(X,ε) (where “vf" stands for “variance function") as follows:


  (a) X satisfies (X)⊂ I.

  (b) X has density p_X(·), and there exists a fixed positive constant C_0 such that 

    sup_x∈ p_X(x) ≤ C_0.


  (c) There exist fixed positive constants c_0 and δ_0 such that

    inf_x^*∈(X)p_X(x^*) ≥ c_0    and 
       inf_0<δ<δ_0inf_x^*∈(X)λ**u∈[-1,1]: x^* + δ u∈(X)≥ c_0,

where λ(·) is the Lebesgue measure on the real line.

  (d)ε^4 ≤ C_ε for some fixed positive constant C_ε.

One can readily verify that P_vf,(X,ε)⊂P_cv,(X,ε), with the latter defined in the beginning of Section <ref>. Compared to P_cv,(X,ε), Condition (c) in P_vf,(X,ε) is posed on the marginal density and support of X, 
since in the variance function case we require a sufficient number of close pairs (X_i, X_j) around each target x^*. We also note that, as in P_cv,(X,ε), no smoothness assumption is posed on the design density in P_vf,(X,ε). 

The rest of the section is devoted to proving, for any fixed positive constants C_F and C_V, the following minimax rates



    inf_Vsup_f∈Λ_α,I(C_F)sup_V∈Λ_β,I(C_V)sup__(X,ε)∈P_vf,(X,ε)sup_x^*∈(X) R_1(V,V;x^*)    ≍ n^-8αβ/4αβ + 2α + β∨n^-2β/2β+1,
    
    inf_Vsup_f∈Λ_α,I(C_F)sup_V∈Λ_β,I(C_V)sup__(X,ε)∈P_vf,(X,ε) R_2(V,V)    ≍ n^-8αβ/4αβ + 2α + β∨ n^-2β/2β+1,

where _(X,ε) denotes the joint distribution of (X,ε), and V(·) ranges over all estimators of V(·).





































 §.§ Upper bound


We now propose an estimator of V(x^*) for some fixed x^*∈(X) by combining pairwise differences with local polynomial regression. We first introduce some notation. Let ℓ be the largest integer strictly smaller than β and 

    q(u)  (1,u,u^2/2!,…, u^ℓ/ℓ!)^⊤.
 




For any 1 ≤ i < j≤ n, define 

    D_ij (Y_i-Y_j)^2/2,  X_ij (X_i +X_j)/2,  and  K_ij K_h_1(X_i - X_j)K_h_2(X_ij - x^*),

where h_1,h_2 are two bandwidths. Define an (ℓ+1)×(ℓ+1) matrix



















    𝐁_n n 2^-1∑_i<jq*X_ij - x^*/h_2q^⊤*X_ij - x^*/h_2K_ij

and 𝐁_n^* as its adjugate such that 𝐁_n𝐁_n^* = 𝐁_n^*𝐁_n = |𝐁_n|𝐈_ℓ+1. For example, when ℓ = 1, we have

    𝐁_n = 
    [ s_0 s_1; s_1 s_2 ]
    ,  𝐁_n^* =
    [  s_2 -s_1; -s_1  s_0 ],    and    |𝐁_n| = s_0s_2 - s_1^2,

where

    s_k n 2^-1∑_i<j*X_ij - x^*/h_2^kK_ij,      k=0,1,2.




Following <cit.>, we propose a robust local polynomial estimator:

    V_(x^*) n 2^-1∑_i<j D_ij(|𝐁_n| + τ_n)^-1q^⊤(0)𝐁_n^*q*X_ij - x^*/h_2K_ij,

where τ_n is some sufficiently small positive constant that decays to 0 polynomially with n. 
Let 

    w_ijn 2^-1q^⊤(0)𝐁_n^*q*X_ij - x^*/h_2K_ij   and  w_ij w_ij/(|𝐁_n| + τ_n).

Then, it holds that V_(x^*) = ∑_i<jw_ijD_ij, ∑_i<jw_ij = |𝐁_n|, and

    ∑_i<jw_ij(X_ij - x^*)^k = ∑_i<jw_ij(X_ij - x^*)^k = 0,      k=1,2,…,ℓ.

The last property (<ref>) is referred to as the reproducing property of local polynomial estimators (cf. Proposition 1.12 in <cit.>).

















































Suppose the kernel K(·) in V_ is chosen such that (<ref>) holds with constants M_K and M_K, τ_n ≍ n^-κ for some fixed constant κ≥ 1, and the bandwidths h_1,h_2 are chosen as

    (h_1,h_2) ≍*n^-2β/4αβ + β + 2α, n^-4α/4αβ + β + 2α,     0< α < β/4β +2,
    *n^-1, n^-1/2β + 1,    α≥β/4β +2.

Then, under (<ref>) with random design, it holds that

    sup_f∈Λ_α,I(C_F)sup_V∈Λ_β,I(C_V)sup__(X,ε)∈P_vf,(X,ε)sup_x^*∈(X) R_1(V_, V;x^*) ≤C *n^-8αβ/4αβ + β + 2α∨ n^-2β/2β +1

and

    sup_f∈Λ_α,I(C_F)sup_V∈Λ_β,I(C_V)sup__(X,ε)∈P_vf,(X,ε) R_2(V_, V) ≤ C*n^-8αβ/4αβ + β + 2α∨ n^-2β/2β +1,

where C is some fixed positive constant that only depends on M_K, M_K,α,β,C_F,C_V and C_0,c_0,C_ε in P_vf,(X,ε).









Variance function estimation in (<ref>) with fixed design X_i = i/n, i∈[n], has been studied in <cit.>. There the minimax rate is 





    inf_Vsup_f∈Λ_α,[0,1](C_F)sup_V∈Λ_β,[0,1](C_V)sup_ε^4 ≤ C_εsup_x^*∈[0,1] R_1(V, V;x^*) ≍ n^-4α∨ n^-2β/(2β + 1),
    
    inf_Vsup_f∈Λ_α,[0,1](C_F)sup_V∈Λ_β,[0,1](C_V)sup_ε^4 ≤ C_ε R_2(V, V) ≍ n^-4α∨ n^-2β/(2β + 1),

with the integral in R_2 under the Lebesgue measure on [0,1]. Comparing the above result with the error rate in Theorem <ref>, we see that the transition boundary in both the fixed and random design settings is α = β/(4β+2). When α≥β/(4β + 2), V(·) under both R_1 and R_2 can be estimated at the classic nonparametric rate n^-2β/(2β + 1) as if the mean function f(·) were known. When α < β/(4β + 2), a faster rate can be achieved in the random design case. This can be intuitively understood by the fact that, by constrast to the fixed design case, a significant portion of pairs have distance smaller than 1/n in the random design setting.





As has been noted in <cit.>, in the fixed design setting, estimating the variance (function) by smoothing the squared residuals obtained from pre-estimation of the mean function f(·) is sub-optimal. The same conclusion also applies to the random design setting. Since the design being fixed or random has no first-order effect on the estimation of the mean, the above method only achieves the rates n^-4α/(2α+1)∨ n^-1 in variance estimation and n^-4α/(2α+1)∨ n^-2β/(2β+1) in variance function estimation, neither of which is minimax optimal. 



Unlike in the fixed design case, once below the threshold α = β/(4β + 2), α and β are now both present in the minimax rate in the random design case, suggesting that the smoothness of V(·) always has an effect on its estimation. This is because variance function estimation in the random design setting is essentially a “two-dimensional" problem, where we have to jointly choose two optimal neighborhood sizes to characterize the closeness between (i) each X_i and X_j;  and (ii) every pair (X_i,X_j) and each target point x^*. By contrast, in the fixed design setting, the distance between X_i and X_j is constrained to be no smaller than 1/n, and thus cannot be jointly optimized with the distance between (X_i,X_j) and x^*.



One might wonder whether the following Nadaraya-Watson type estimator can be used to establish the upper bound in Theorem <ref>:

    V_NW(x^*) ∑_i<j K_h_1(X_i - X_j)K_h_2(X_ij - x^*)D_ij/∑_i<j K_h_1(X_i - X_j)K_h_2(X_ij - x^*),

where K(·) is now chosen to be a higher-order kernel to further reduce bias when β > 1. It turns out that the analysis of V_NW requires an extra assumption on the smoothness of the density p_X(·) which can be completely avoided with V_. Moreover, it is well-known that local polynomial estimators have good finite sample properties and boundary performances when X is compactly supported <cit.>.















<cit.> considered minimax estimation of the variance function (and more generally, its derivatives) in the context of nonparametric regression with random design. We focus on the comparison of their results on variance function estimation with ours. Their lower bound (Theorem 1 therein) is proved independent of the smoothness level of the mean function and upper bound (Theorem 4 therein) is proved under sufficient smoothness on the mean function. Therefore their minimax rate is only comparable to the n^-2β/(2β + 1) component in ours.  In this case, their lower bound of the order n^-(2β - 1)/(2β) is proved over the following class of variance function:

    S_β*1 + ∑_k=1^∞δ_ke_k: |δ_k|≲ k^-β

for any β > 1, where {e_k}_k=1^∞ is an arbitrary basis on L^2([-π,π]). Moreover, continuous differentiability of the error density is required in their paper. In contrast, we pose no smoothness conditions on the error density, and neither S_β nor S_β+1/2 can be embedded in the β-Hölder class Λ_β considered in our setting (e.g., f(x)=|x| with domain [-π,π] belongs to S_2 but is not 1.5- or 2-Hölder smooth since it is not differentiable at the origin). In summary, the results in <cit.> neither imply nor contradict the n^-2β/(2β+1) part in our minimax rate, and our results are more refined since they characterize the exact elbow α = β/(4β+2) and also the minimax rate below this threshold.

















 §.§ Lower bound


The following are matching lower bounds to Theorem <ref>.



Under (<ref>) with random design, for any x^*∈(X), 

    inf_Vsup_f∈Λ_α,I(C_F)sup_V∈Λ_β,I(C_V)sup__(X,ε)∈P_vf,(X,ε) R_1(V, V;x^*) ≥ c *n^-8αβ/4αβ + β + 2α∨ n^-2β/2β +1,

where c is some fixed positive constant that only depends on α,β, C_F, C_V and C_0,c_0,C_ε in P_vf,(X,ε), and V ranges over all estimators of V.




Under (<ref>) with random design, 

    inf_Vsup_f∈Λ_α,I(C_F)sup_V∈Λ_β,I(C_V)sup__(X,ε)∈P_vf,(X,ε) R_2(V, V) ≥ c*n^-8αβ/4αβ + β + 2α∨ n^-2β/2β +1,

where c is some fixed positive constant that only depends on α,β, C_F, C_V and C_0,c_0,C_ε in P_vf,(X,ε), and V ranges over all estimators of V.









Due to the appearances of both α and β in the nontrivial n^-8αβ/4αβ + β + 2α part of the minimax rate, proving the above two results is more involved than proving Theorem <ref>. In particular, it takes an extra step of localization in the construction of the mean function f(·) as well as V(·). More precisely, for the lower bound at a target point x^* in Theorem <ref>, our construction of both f(·) and V(·) only has variation within a small neighborhood of x^*. Such localized construction is not necessary in the fixed design setting, since when proving the n^-4α component therein (see Remark <ref>), the variance function can simply be taken as a constant. 

In what follows, we give a proof sketch of the nontrivial n^-8αβ/(4αβ + β + 2α) component of the lower bound in Theorem <ref> for α < β/(4β + 2); the proof of Theorem <ref> can be seen as an extension of Theorem <ref> via a standard construction of multiple hypotheses. We assume the support of X is contained in I = [0,1], and for clarity of illustration, here we present the construction for an interior point x^*∈ (0,1)⋂(X). The proof works for boundary points as well.


We continue to adopt the two-step approach introduced in the proof sketch of Theorem <ref> in Section <ref>. The second step is very similar with the help of Lemmas <ref> and <ref>, so we will focus on the construction under the null H_0 and alternative H_1 in the first step. Choose the parameters

    h_1 ≍ n^-2β/4αβ + β + 2α,    h_2 ≍ n^-4α/4αβ + β + 2α,    and   θ_n^2 =  h_1^2α = h_2^β

so that h_2/h_1→∞ as n→∞.


  Choice of V(·): Under H_0 let V≡ 1. Under H_1, let V(·) be one minus a smooth bump function around x^* with width h_2 and height h_2^β so that V(x^*) = 1-θ_n^2.

  Choice of f(·): Under H_0 let f≡ 0. Under H_1, let f(·) be a “local" version of the design in Theorem <ref>. That is, f takes on a value of 0 outside of [x^* - h_2, x^* + h_2], and inside that h_2-neighborhood of x^*, f is piecewise trapezoidal with upper base length 2h_1, lower base length 4h_1 and height {h_1^αr_i}_i=1^N for a standard normal sequence {r_i}_i=1^N with N h_2/(2h_1) a positive integer. 

  Choice of ε: Under both H_0 and H_1, let ε∼N(0,1). 

  Choice of X: Under both H_0 and H_1, let X be uniformly distributed on the union of [0,1]\ [x^*-h_2,x^* + h_2] and the upper bases of all the trapezoids inside [x^* - h_2, x^* + h_2]. 


See Figure <ref> for an illustration of H_1. 

Under the above construction, the squared distance between the null and alternative hypotheses (1 - (1 - θ_n^2))^2 = θ_n^4 ≍ n^-8αβ/4αβ + β + 2α is the desired minimax rate. Using Lemma <ref>, we can show that

    (_0, _1) ≲θ_n^2nh_1^1/2h_2^1/2≤ c

for some sufficiently small c, where _0 and _1 represent the joint distribution of {(X_i, Y_i)}_i=1^n under H_0 and H_1, respectively. The detailed proof is presented in the supplement.





§ DISCUSSION


The two univariate models (<ref>) and (<ref>) discussed in the previous two sections raise natural questions about possible extensions to the multivariate setting. In what follows, we first present some partial results in this direction in the sense of (<ref>) and (<ref>). We then establish some connections between our study and quadratic functional estimation and variance estimation in the linear model. Lastly, we discuss two more extensions of (<ref>) in the direction of adaptive estimation and mean function with inhomogeneous smoothness. Throughout, consider C_F, C_σ, C_0,c_0,C_ε to be fixed positive constants.



 §.§ Multivariate nonparametric regression


Consider the following multivariate version of (<ref>):

    Y_i = f(X_i) + σε_i,    i=1,2,…, n,

where {X_i}_i=1^n = {(X_i,1,…,X_i,d)^⊤}_i=1^n are i.i.d. copies of X = (X_1,…,X_d)^⊤ in ^d for some fixed positive integer d, {ε_i}_i=1^n are i.i.d. copies of ε with zero mean and unit variance and are independent of {X_i}_i=1^n, and f:^d→ belongs to a d-dimensional anisotropic Hölder class with smoothness index α = (α_1,…,α_d)^⊤ defined below. The goal is to estimate σ^2 with f(·) and the distribution of X as nuisance parameters. This problem has been studied in <cit.>, <cit.>, <cit.>, to name a few, again with a focus on the fixed design setting. 


Let I_1,…,I_d be d fixed (possibly infinite) intervals on  and let I be their Cartesian product I_1×…× I_d⊂^d. Following <cit.> and <cit.>,  we define an anisotropic Hölder class Λ_α,I(C_F) on I as follows. For any x∈I and k∈[d], let f_k(·|x_-k) denote the univariate function y↦ f(x_1,…,x_k-1,y,x_k+1,…,x_d), with x_-k defined as x without the kth component. Then, Λ_α,I(C_F) is defined as all f:I↦ such that

    max_1≤ k≤ dmax_0≤ j≤α_ksup_x∈I*f_k^(j)(·|x_-k)_∞≤ C_F

and

    max_1≤ k≤ dsup_x∈Isup_y_1,y_2∈ I_k*f_k^(α_k)(y_1|x_-k) - f_k^(α_k)(y_2|x_-k)/|y_1 - y_2|^α_k^'≤ C_F,

where again α_k is the largest integer strictly smaller than α_k and α_k^'α_k - α_k. Let (X) be the support of X.









Define P_mcv,(X,ε) (where “mcv" stands for “multivariate constant variance") as the multivariate counterpart of P_cv,(X,ε):


  (a) X satisfies (X)⊂I.

  (b) X has density p_X(·) and there exists a fixed positive constant C_0 such that  

    sup_u∈^d p_X(u) ≤ C_0 .


  (c)There exist two fixed constants δ_0 > 0 and c_0>0 such that for any δ∈^d that satisfies δ_∞ < δ_0, there exists a set UU_δ⊂[-1,1]^d such that

    λ(U_δ)≥ c_0    and   inf_u∈U_δp_(u_1δ_1,…,u_dδ_d)≥ c_0,

where λ(·) represents the Lebesgue measure on ^d. 

  (d)ε^4 ≤ C_ε for some fixed positive constant C_ε.



For an upper bound on the minimax risk, we propose the following multivariate extension of (<ref>) via a product kernel (again with convention 0/0 = 0):

    σ_d^2 n 2^-1∑_i<j*∏_k=1^d K_h_k(X_i,k - X_j,k)(Y_i - Y_j)^2/2/n 2^-1∑_i<j*∏_k=1^d K_h_k(X_i,k-X_j,k),

where K(·) is a kernel chosen to satisfy (<ref>), and {h_k}_k=1^d is a kernel bandwidth sequence. 








































In the following results, we will use α to denote the harmonic mean of the d-dimensional smoothness index α, i.e. α d/(∑_k=1^d 1/α_k). This quantity is known as the effective smoothness in classical problems such as anisotropic density estimation <cit.> and anisotropic function estimation <cit.>.  



Suppose 0< α_k ≤ 1, k∈[d]. Suppose the kernel K(·) in σ_d^2 is chosen such that (<ref>) is satisfied with constants M_K and M_K, and the bandwidth sequence is chosen as h_k ≍ n^-2α/(α_k(4α+d)) for all k∈[d]. Then, under (<ref>) with random design, it holds that

    sup_f∈Λ_α,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P_mcv,(X,ε)*σ_d^2 - σ^2^2 ≤ C*n^-8α/(4α + d)∨ n^-1,
 
where C is some fixed positive constant that only depends on M_K, M_K,α,C_F,C_σ and C_0,c_0,C_ε in P_mcv,(X,ε). 








Under (<ref>) with random design, it holds that

    inf_σ^2sup_f∈Λ_α,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P_mcv,(X,ε)*σ^2 - σ^2^2 ≥ c*n^-8α/(4α + d)∨ n^-1,
 
where c is some fixed positive constant that only depends on α,C_F, C_σ and C_0,c_0,C_ε in P_mcv,(X,ε), and σ^2 ranges over all estimators of σ^2. 






We note that Proposition <ref> is only proved for α_k∈(0,1], k∈[d]. The general case when α_k is possibly larger than 1 is much more involved due to the difficulty in the random design analysis. Propositions <ref> and <ref>, combined, imply that the minimax rate is n^-8α/(4α + d)∨ n^-1 for α_k∈(0,1], k∈[d]. In particular, when f is in an isotropic α-Hölder class (0<α≤ 1), this rate becomes n^-8α/(4α + d)∨ n^-1. We also remark that a different estimator achieving the rate n^-8α/(4α + d)∨ n^-1 over an isotropic α-Hölder class has been briefly sketched in <cit.>.

For completeness, we also state without proof some results for model (<ref>) in the fixed design setting. In particular, we consider the following two types of fixed designs in the d-dimensional unit cube [0,1]^d, namely, the grid design (GD):

    (X_(i_1,…,i_d),1,…, X_(i_1,…,i_d),d) = (i_1/n^1/d, …, i_d/n^1/d), 
     
    (i_1,…,i_d) ∈ [n^1/d]×…×[n^1/d]

assuming n^1/d is an integer, and the diagonal design (DD):

    (X_i,1,…,X_i,d) = (i/n, …, i/n),   i∈[n].

Here for any positive integer n, [n] denotes the set {1,2,…,n}. Let α_maxmax_k∈[d]α_k and α_minmin_k∈[d]α_k. The first result for (GD) is a simple modification of the isotropic result in <cit.> by taking differences along the smoothest direction with index α_max. The second result can be readily deduced from the fact that Y_i = f(i/n) + σε_i, i∈[n], where f(x) f(x,…,x) is α_min-Hölder smooth. 




Under (<ref>) with fixed design (GD), it holds that

    inf_σ^2sup_f∈Λ_α,[0,1]^d(C_F)sup_σ^2≤ C_σsup_ε^4≤ C_ε*σ^2 - σ^2^2 ≍ n^-4α_max/d∨ n^-1

up to some fixed positive constant that only depends on α, C_F,C_σ,C_ε, where σ^2  ranges over all estimators of σ^2.




Under (<ref>) with fixed design (DD), it holds that

    inf_σ^2sup_f∈Λ_α,[0,1]^d(C_F)sup_σ^2≤ C_σsup_ε^4≤ C_ε*σ^2 - σ^2^2 ≍ n^-4α_min∨ n^-1

up to some fixed positive constant that only depends on α, C_F,C_σ,C_ε, where σ^2 ranges over all estimators of σ^2.


When f(·) belongs to an isotropic α-Hölder class, Proposition <ref> implies the minimax rate n^-4α/d∨ n^-1 derived in <cit.>. Comparison with the random design rate n^-8α/(4α + d)∨ n^-1 thus shows that, for 0 <α≤ 1, a faster rate is again achievable in the random design setting for α < d/4. 





 §.§ Nonparametric additive model



Consider variance estimation in the additive model (<ref>):

    Y_i = ∑_k=1^d f_k(X_i,k) + σε_i,    i=1,2,…,n,

for some fixed integer d ≥ 2, where {ε_i}_i=1^n are i.i.d. with zero mean and unit variance and are independent from {X_i}_i=1^n = {(X_i,1,…,X_i,d)^⊤}_i=1^n in the random design setting. Unlike Section <ref>, we specify d≥ 2, since the minimax rate in the fixed design (GD) has completely different behavior for d=1 and d≥ 2 (see Proposition <ref> below).




  §.§.§ Fixed design

We first consider the two fixed designs (GD) and (DD) defined in (<ref>) and (<ref>). For both designs, we consider an error distribution class with only a finite fourth moment condition. We start with (GD), where by iteratively taking pairwise differences, one is able to estimate the variance at the parametric rate n^-1 without any smoothness assumption on the additive components {f_k}_k=1^d. For simplicity, we illustrate this idea with d = 2 with two additive components f(·) and g(·), and assume that √(n) is an even number. In this case,

    Y_i,j = f*i/√(n) + g*j/√(n) + σε_i,j,    (i,j)∈[√(n)]× [√(n)],
  
where {ε_i,j}_i,j∈[√(n)] are i.i.d. with zero mean and unit variance. By taking the pairwise difference in the first dimension, we have

    Y_(i_1,i_2), j Y_i_1,j - Y_i_2,j = f*i_1/√(n) - f*i_2/√(n) + σ(ε_i_1,j - ε_i_2,j)

for all j∈[√(n)] and (i_1,i_2)∈[√(n)]× [√(n)] such that i_1≠ i_2. Taking again the pairwise difference in the second dimension, we have

    Y_(i_1,i_2),(j_1,j_2) Y_(i_1,i_2),j_1 - Y_(i_1,i_2),j_2 = σ(ε_i_1,j_1 - ε_i_2,j_1 - ε_i_1,j_2 + ε_i_2,j_2)
 
for all (i_1,i_2,j_1,j_2)∈[√(n)]× [√(n)]×[√(n)]×[√(n)] such that i_1≠ i_2 and j_1≠ j_2. Clearly, we have Y_(i_1,i_2),(j_1,j_2) = 0 and (Y_(i_1,i_2),(j_1,j_2)) = 4σ^2. Let m √(n)/2 and define I{(1,2), (3,4), …, (2m-1, 2m)} with cardinality m. Then, for the set of data points {Y_(i_1,i_2),(j_1,j_2)}_(i_1,i_2),(j_1,j_2)∈I with cardinality m^2 = n/4, it can be readily verified that they are i.i.d. with mean 0 and variance 4σ^2. Therefore, with Y̅ defined as the sample average of {Y_(i_1,i_2), (j_1,j_2)}_(i_1,i_2),(j_1,j_2)∈I, the sample variance estimator,

    σ_add, GD^2 1/n∑_(i_1,i_2), (j_1,j_2)∈I*Y_(i_1,i_2),(j_1,j_2) - Y̅^2,

achieves the parametric rate n^-1. A similar derivation holds for general d.

















Suppose d≥ 2. Under (<ref>) with fixed design (GD), it holds that

    inf_σ^2sup_f_k,k∈[d]sup_σ^2≤ C_σsup_ε^4≤ C_ε*σ^2 - σ^2^2 ≍ n^-1

up to some fixed positive constant that only depends on C_σ and C_ε, where σ^2  ranges over all estimators of σ^2, and the first supremum is taken over all functions defined on [0,1] for each k∈[d].


Now we move on to the design (DD), where we assume each additive component f_k in (<ref>) is α_k-Hölder smooth on [0,1] with some fixed constant C_F. In this case, the model can equivalently be written as 

    Y_i = f(i/n) + σε_i,    i=1,2,…,n,

where f∑_k=1^d f_k is α_min-Hölder smooth. Therefore, the univariate estimator and lower bound in <cit.> can be directly applied. 











Under (<ref>) with fixed design (DD), it holds that

    inf_σ^2sup_f_k∈Λ_α_k,[0,1](C_F),k∈[d]sup_σ^2≤ C_σsup_ε^4≤ C_ε*σ^2 - σ^2^2 ≍ n^-4α_min∨ n^-1

up to some fixed positive constant that only depends on C_F,C_σ,C_ε, where σ^2  ranges over all estimators of σ^2.

Comparison of Propositions <ref> and <ref> shows that, in contrast to grid design (GD) and random design below, there is no gain from an additive structure in the mean function for the diagonal design (DD). 
 


  §.§.§ Random design

We now discuss (<ref>) with a random design for {X_i}_i=1^n when f_k is α_k-Hölder smooth on some fixed set I_k for each k∈[d]. Since a shift in the mean does not affect the estimation of variance, we assume f_k(X_1,k) = 0 for each k∈[d] for simplicity. 




Recall the definition of P_cv,(X,ε) in the beginning of Section <ref>. Define the joint distribution class P_add,(X,ε) (where “add" stands for “additive") as:


   For each k∈[d], the joint distribution of (X_k,ε) belongs to P_cv,(X,ε) and the components of X are mutually independent.










In view of Theorem <ref>, the following lower bound is immediate.


Under (<ref>) with random design, it holds that

    inf_σ^2sup_f_k∈Λ_α_k,I_k(C_F),k∈[d]sup_σ^2≤ C_σsup__(X,ε)∈P_add,(X,ε)*σ^2 - σ^2^2 ≥ c*n^-8α_min/4α_min+1∨ n^-1,

where c is a fixed positive constant that only depends on α, C_F,C_σ and C_0,c_0,C_ε in P_add,(X,ε), and σ^2  ranges over all estimators of σ^2.
 



We now describe a procedure that matches the lower bound in Proposition <ref>, but depends crucially on mutual independence. For illustrative purposes, we again consider the case of only two additive components f(·) and g(·), which are α- and β-Hölder smooth, respectively. Let X and W denote the two covariates. For each i∈[n], define

    ε^X_i f(X_i) + σε_i    and   ε^W_i  g(W_i) + σε_i,

and their corresponding variances

    σ_X^2  f^2(X) + σ^2   and   σ_W^2  g^2(W) + σ^2.

Clearly, we have ε^X_i = 0 and ε^W_i = 0, and ε^X_i and ε^W_i are independent of g(W_i) and f(X_i), respectively. Now, notice that the additive model in (<ref>) can be equivalently viewed as Y_i = f(X_i) + ε^W_i. Thus by applying the univariate kernel estimator defined in (<ref>) to {(Y_i,X_i)}_i=1^n, which we denote as σ^2_W, one obtains

    *σ^2_W - σ_W^2^2 ≤ C(n^-8α/(4α + 1)∨ n^-1)

for some fixed positive constant C. Similarly, defining σ^2_X as σ^2_W, one has

    *σ^2_X - σ_X^2^2  ≤ C(n^-8β/(4β + 1)∨ n^-1).

Lastly, under a finite fourth moment assumption on ε, a sample variance estimator of {Y_i}_i=1^n, denoted as σ^2_Y, achieves the parametric rate n^-1 in estimating the total variance (Y), which can be decomposed as f^2(X) +  g^2(W) + σ^2.



Consequently, we have shown that the method-of-moments estimator 

    σ_moment,2^2σ_X^2 + σ_W^2 - σ^2_Y

achieves the optimal rate in Proposition <ref>.




We summarize the above derivation for the natural extension σ_moment,d^2 to general d.



Under (<ref>) with random design, it holds that 

    sup_f_k∈Λ_α_k,I_k(C_F),k∈[d]sup_σ^2≤ C_σsup__(X,ε)∈P_add,(X,ε)*σ_moment,d^2 - σ^2^2 ≤ C*n^-8α_min/4α_min+1∨ n^-1,

where C is some fixed positive constant that only depends on α,C_F,C_σ and C_0,c_0,C_ε in P_add,(X,ε).


Propositions <ref> and <ref> together imply the minimax rate over P_add,(X,ε), which further illustrates the fact that an additive structure in the mean function could possibly avoid the “curse of dimensionality" in variance estimation. 
However, we note that our results crucially rely on the mutual independence condition. It is still largely unclear if the same minimax rate could apply to the general case without this condition, though a discussion of an interesting connection to variance estimation under linear models shall be made in Section <ref>. 



 §.§ Connection to quadratic functional estimation











We now formally state the connection between quadratic functional estimation and variance estimation in (<ref>), the first of which has been studied in, for example, <cit.>, <cit.>, <cit.>, and <cit.>.

Recall the definition of Q in (<ref>) with some non-negative weight function w(·). Squaring both sides of   (<ref>), multiplying by w(X_i), and then taking the expectation, one has

    *Y_i^2w(X_i) = *f^2(X_i)w(X_i) + σ^2(w(X_i)ε_i^2) = Q + σ^2 w(X_i).
 
Under a finite fourth moment assumption on ε, both *Y_i^2w(X_i) and w(X_i) can be estimated at the parametric rate via the sample mean estimator, and σ^2 can be estimated via σ^2 in (<ref>) with rate n^-8α/(4α + 1)∨ n^-1 under the quadratic risk. Therefore, the estimator

    Q1/n∑_i=1^n Y_i^2w(X_i) - *1/n∑_i=1^n w(X_i)·σ^2

achieves the same rate n^-8α/(4α + 1)∨ n^-1. In fact, it is not possible to improve upon this rate since if there exists an estimator Q with a faster convergence rate, then the “conjugate" estimator of σ^2 defined as

    σ^2 max*1/n∑_i=1^n Y_i^2w(X_i) - Q/1/n∑_i=1^n w(X_i), 0·1*1/n∑_i=1^n w(X_i) > 0

will also converge to σ^2 at a faster rate, violating the lower bound in Theorem <ref>. 

The following result summarizes the derivation. Recall the definition of P_cv,(X,ε) in the beginning of Section <ref>.



Suppose the weight function w(·) in the definition of Q is uniformly bounded on . Then, it holds that

    inf_Qsup_f∈Λ_α,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P_cv,(X,ε)*Q - Q^2 ≍ n^-8α/(4α + 1)∨ n^-1

up to some fixed positive constant that only depends on w(·), α, C_F, C_σ and C_0,c_0,C_ε in P_cv,(X,ε), where Q  ranges over all estimators of Q.





 §.§ Connection to the linear model


Throughout this paper, we have treated the distribution of X as a nuisance parameter. Interestingly, when we do know the distribution of X, variance estimation in nonparametric regression with random design becomes substantially easier with the aid of parallel work in the high-dimensional linear model <cit.>. We first elaborate on this point using the simple model (<ref>), and then formulate corresponding results for (<ref>) and (<ref>).

By applying the inverse of the distribution function F of X, (<ref>) can be equivalently written as

    Y_i = f(U_i) + σε_i,    i=1,2,…,n,

where {U_i}_i=1^n = {F(X_i)}_i=1^n are i.i.d. uniform on [0,1], and f(·)  f∘ F^-1(·) is still α-Hölder smooth under Lipschitz continuity on F^-1. Then, using a wavelet expansion for Hölder classes (cf. Proposition 2.5 in <cit.>), one has

    Y_i = f_1(U_i) + ∑_j=1^2^Jψ_j(U_i) + σε_i,    i=1,2,…,n,

where {ψ_j}_j=1^∞ is an L_2-orthonormal wavelet basis under the Lebesgue measure on [0,1], and f_1(·) is the remainder term after truncation at resolution J = J_n which satisfies f_1_∞ = O(2^-α J_n). Let ψ (ψ_1,…, ψ_2^J) and assume without loss of generality that ψ = 0_2^J, since a mean shift does not affect the estimation of variance. Moreover, due to the orthonormality of {ψ_j}_j=1^∞, we have (ψ) = (ψψ^⊤) = 𝐈_2^J. Following <cit.> and <cit.>,  the estimator

    σ^2_proj1/n-1∑_i=1^n (Y_i - Y̅)^2 - n 2^-1∑_i<j Y_iY_j ψ^⊤(U_i)ψ(U_j)

has a variance term of the order (2^J_n + n)/n^2 and a bias term of the order 2^-2α J_n. Therefore, by choosing the optimal truncation level 2^J_n≍ n^2/(4α+1), σ^2_proj recovers the optimal rate n^-8α/(4α+1)∨ n^-1 in Theorem <ref>.













Define σ^2_proj, d (with tensor wavelet basis) and σ^2_proj, add as the natural extensions of σ^2_proj under (<ref>) and (<ref>), respectively (see the proofs of Propositions <ref> and <ref> in the supplement for exact definitions). In the wavelet expansion, we will use J_k to denote the truncation level for the kth component of f(·) in (<ref>) and f_k in (<ref>), and we use F_k to denote the marginal distribution of X_1,k. Recall that α= d/(∑_k=1^d 1/α_k) for α = (α_1,…,α_d)^⊤.



Suppose the distribution of X is known with (X)⊂I for some fixed set I⊂^d, and F_k^-1(·) is Lipschitz continuous for all k∈[d] with some fixed positive constant. Then, when 2^J_k is chosen to be of the order n^2α/(α_k(4α+d)) for k∈[d] in σ^2_proj,d, it holds that 




    sup_f∈Λ_α,I(C_F)sup_σ^2≤ C_σsup_ε^4≤ C_ε*σ^2_proj,d - σ^2^2 ≤ C*n^-8α/(4α + d)∨ n^-1,

where C is some fixed positive constant that only depends on α,C_F,C_σ,C_ε, and the distribution of X.




Suppose the distribution of X is known with (X)⊂ I_1×…× I_d for some fixed intervals I_1,…,I_d on the real line, and F_k^-1(·) is Lipschitz continuous for all k∈[d] with some fixed positive constant. Then, when 2^J_k is chosen to be of the order n^2α_k/(4α_k+1) for k∈[d] in σ^2_proj,add, it holds that 




    sup_f_k∈Λ_α_k,I_k(C_F),k∈[d]sup_σ^2≤ C_σsup_ε^4≤ C_ε*σ^2_proj,add - σ^2^2 ≤ C*n^-8α_min/4α_min+1∨ n^-1,

where C is some fixed positive constant that only depends on α,C_F,C_σ,C_ε, and the distribution of X.


As in the classical setting of mean function estimation via orthogonal series, the difference of the rates in Propositions <ref> and <ref> is clearly explained by the number of wavelet bases used to approximate f in (<ref>) and {f_k}_k=1^d in (<ref>). We also note that, quite interestingly,  Proposition <ref> gives results beyond the case 0<α_1,…,α_d≤ 1 considered in Proposition <ref>, and Proposition <ref> does not rely on the mutual independence of the components of X.



 §.§ Adaptive estimation of constant variance


In this subsection, we consider adaptive estimation of the variance σ^2 in model (<ref>). This is achieved by a Lepski-type procedure <cit.>. Let σ^2(h) be the estimator in (<ref>) with an explicit dependence on the bandwidth parameter h. For any given sample size n and fixed positive constant δ, define two positive integers m_1 and m_2 such that 2^-m_1≤ n^-1≤ 2^-m_1+1 and 2^-m_2-1≤ n^-(2-δ)≤ 2^-m_2, and define the following dyadic grid

    H_δ*2^-j: m_1≤ j≤ m_2, j∈ℤ.

Then, define the estimator σ_adapt^2 σ^2*h_δ with

    h_δmax*h∈H_δ: *σ^2(h) - σ^2(h^')≤τ (log n)^1/2n^-1(h^')^-1/2, ∀ h^'∈H_δ, h^' <h

for some sufficiently large positive constant τ. If the set being maximized is empty, we will take h_δ = n^-(2-δ).

To state the error bound of σ_adapt^2, we need the following variant P^adapt_cv, (X,ε) of the distribution class P_cv, (X,ε) considered in Theorem <ref>, where we replace the finite fourth-moment assumption (d) therein by the stronger sub-Gaussian tail condition:


  (d^') There exist some fixed positive constants C_1,ε and C_2,ε such that exp(tε)≤ C_1,εexp(C_2,εt^2) for any t∈ℝ.

A similar exponential moment assumption has been made in the context of adaptive estimation under fixed design (cf. Theorems 1 and 2 in <cit.>).




For any given sufficiently small fixed α_* > 0, fix some δ_*∈(0, 8α_*/(4α_*+1)). Suppose the kernel K(·) in σ_adapt^2 = σ^2*h_δ_* is chosen such that (<ref>) is satisfied with constants M_K and M_K, and τ in h_δ_* is chosen to be sufficiently large (only depending on δ_*, C_1,ε, C_2,ε). Then, under (<ref>) with random design, it holds uniformly over all α≥α_* that

    sup_f∈Λ_α,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P^adapt_cv, (X,ε)*σ_adapt^2 - σ^2^2 ≤ C**log n/n^2^4α/(4α+1)∨ n^-1,

where C is some fixed positive constant that only depends on δ_*, M_K,M_K, C_F, C_σ, and C_0,c_0,C_1,ε,C_2,ε in P^adapt_cv, (X,ε).




The following proposition shows that the extra poly-logarithmic term cannot be removed.



Let ϕ_n,α (log n/n^2)^2α/(4α+1)∨ n^-1/2 for any α > 0 and positive integer n. Consider any fixed positive α_* and α_*≤α_1<α_2<∞. Then, for any sufficiently large n and sufficiently small fixed positive constant c, any estimator σ^2 will satisfy that, if

    sup_f∈Λ_α_2,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P^adapt_cv, (X,ε)*(σ^2 - σ^2)/ϕ_n,α_2^2 ≤ c,

then

    sup_f∈Λ_α_1,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P^adapt_cv, (X,ε)*(σ^2 - σ^2)/ϕ_n,α_1^2 ≥ c.





The above two results combined are in line with analogous adaptation results in quadratic functional estimation <cit.>.




































































§ PROOF OF THEOREM <REF>




We will only prove the lower bound n^-8α/(4α + 1) in the regime 0<α<1/4 since for α≥ 1/4, the rate reduces to the parametric rate n^-1 and the proof is straightforward. Throughout the proof, C represents a generic sufficiently large positive constant and c represents a generic sufficiently small positive constant always taken to be smaller than 1/4. Both C and c only depend on α,C_F,C_σ,C_ε, C_0,c_0 and might have different values for each occurrence. By appropriately rescaling the parameters in the lower bound construction, without loss of generality, we assume that the sample size n and the constants C_F,C_σ,C_ε, C_0 are sufficiently large, c_0 is sufficiently small, and [0,1]⊂ I. 

We will make use of Le Cam's two point method. Introduce the following constants:

    θ_n^2  h_n^2α cn^-4α/(4α+1)   and    N  N_n  1/(6h_n),

where we tune the constant c in h_n so that N is a positive integer. We now specify f(·), distribution of X and distribution of ε in the null and alternative hypotheses, H_0 and H_1, respectively.



   Choice of σ^2: Under H_0, let σ^2 = 1+θ_n^2. Under H_1, let σ^2 = 1.

   Choice of ε: Under both H_0 and H_1, let ε∼N(0, 1).

   Choice of X: Under both H_0 and H_1, let X be uniformly distributed on the union of the intervals [(6i-5)h_n, (6i-1)h_n] for i∈[N].

   Choice of f(·): Under H_0, let f ≡ 0. Under H_1, let f take the value h_n^αr_i on [(6i-5)h_n, (6i-1)h_n], where {r_i}_i=1^N are N i.i.d. symmetric and bounded random variables with distribution  satisfying

    ∫_-∞^∞ x^j(dx) = ∫_-∞^∞ x^jφ(x)dx,    j=1,…, q,

where q is some fixed odd integer strictly larger than 1+1/(2α). Let f be 0 at points 6(i-1)h_n for i∈[N], and then linearly interpolate f for the rest of the unspecified points on [0, 1].

See Figure <ref> for an illustration. In the definition of f(·) under H_1, the existence of the distribution  is guaranteed by Lemma <ref>, and the range of {r_i}_i=1^N, which we denote as B, only depends on α. 

Clearly, σ^2 ≤ C_σ under both H_0 and H_1. Moreover, f(·) under both H_0 and H_1 belongs to Λ_α,[0,1](C_F) due to the boundedness of {r_i}_i=1^N in H_1. Next, we show that the joint distribution of (X,ε) belongs to P_cv,(X,ε). Condition (d) clearly holds and Condition (a) holds with I = [0,1]. Condition (b) holds as well by the fact that p_X(u) = 3/2 for u∈[(6i-5)h_n, (6i-1)h_n] for i∈[N] and p_X(u) = 0 otherwise. Lastly, for Condition (c), it holds by the convolution formula that for any 0 < u < 1/2

    p_(u)    = ∫_u^1 p_X(t)p_X(t-u)dt ≥∑_i=u/(6h_n)+1^N∫_(6i-5)h_n^(6i-1)h_n p_X(t)p_X(t-u)dt
       ≥∑_i=u/(6h_n)+1^N  3/2·3/2· 2h_n ≥3/8 - 9h_n ≥1/4

for sufficiently large n. Here, the second inequality follows from the fact that for any fixed t∈[(6i-5)h_n, (6i-1)h_n], p_X(t) = 3/2 and p_X(t-u)=0 on a subset with Lebesgue measure at most 2h_n. By symmetry of , Condition (c) also holds with δ_0 = 1/2 and U_δ≡ [-1,1].



Denote by σ^2_i, f_i, _i,(X,ε), i=0,1, the choice of σ^2, f, and _(X,ε) under H_0 and H_1, respectively. Let π be the distribution on Λ_α,I(C_F)  such that f_1∼π. Moreover, let _σ^2,f,_(X,ε) represent the expectation with respect to the model (<ref>) with parameters σ^2, f, _(X,ε). Then, we have

    inf_σ^2sup_f∈Λ_α,I(C_F)sup_σ^2≤ C_σsup__(X,ε)∈P_cv,(X,ε)*σ^2 - σ^2^2
       ≥inf_σ^2*1/2_σ^2_0, f_0, _0,(X,ε)*σ^2 - σ^2^2 + 1/2∫_σ_1^2,f,_1,(X,ε)*σ^2 - σ^2^2dπ(f)
       ≥inf_σ^2*1/2_σ^2_0, f_0, _0,(X,ε)*σ^2 - σ^2^2 + 1/2_σ_1^2,f_1,_1,(X,ε)*σ^2 - σ^2^2,
 

where the first inequality follows by lower bounding the maximum risk with Bayes risk with prior π. In what follows, we will use _0 and _1 to denote the joint distribution of {Y_i,X_i}_i=1^n under H_0 and H_1, respectively. Note that the choice of θ_n^2 in (<ref>) leads to the desired lower bound under the quadratic loss. Therefore, adopting the standard reduction scheme with Le Cam's two point method (cf. Theorem 2.2 in <cit.>), it suffices to show that (_0,_1) ≤ c < 1. To show this, let {r_i}_i=1^N be N i.i.d. standard normal random variables, and _1 be the joint distributions of {X_i,Y_i}_i=1^n under H_1 with {r_i}_i=1^N replaced by {r_i}_i=1^N. Then, by triangle inequality, we have

    (_0, _1) ≤(_0, _1) + (_1, _1).
 
We will show (_0, _1) ≤ c and (_1, _1) ≤ c seperately. 

For the first inequality, define x (x_1,…,x_n), dx dx_1… dx_n and similarly for y and dy. Denote p_0, p_1, and p_1 as the densities of _0, _1, and _1 with respect to the Lebesgue measure. Then, we have

    (_0, _1)    = 1/2∫∫*p_0(x, y)-p_1(x,y) dxdy
       = ∫ p(x)dx*1/2∫*p_0(y|x) - p_1(y|x)dy
       =  ∫ p(x)dx(_0(y|x), _1(y|x)),

where p(x)∏_i=1^np_X(X_i) stands for the common density of {X_i}_i=1^n under _0 and _1. Note that under _0, y|x∼N_n(0, Σ_0), with Σ_0 = (1+θ_n^2)𝐈_n. Define {b_i}_i=1^n to be the location index sequence of {X_i}_i=1^n taking values in [N], that is, 

    b_i = j   if   X_i ∈[(6j-5)h_n, (6j-1)h_n].

 Then, due to the symmetry of {r_i}_i=1^N and design of the nonparametric component f, it holds that under _1, y|x∼N_n(0, Σ_1), with (Σ_1)_ii = 1 + h_n^2α = 1 + θ_n^2 and (Σ_1)_ij = h_n^2α1{b_i = b_j} for i≠ j. Define N_0∑_i≠ j1{b_i = b_j}. Since Σ_1 is positive definite (see Lemma <ref> in the supplement), we have by Lemma <ref> that

    (_0(y|x), _1(y|x)) ≤ Cθ_n^2/1+θ_n^2N_0^1/2≤ Cθ_n^2N_0^1/2.

Note that N_0 is a random variable that depends on {X_i}_i=1^n, and by (<ref>) and Jensen's inequality we have 

    (_0, _1) ≤ Cθ_n^2 N_0^1/2≤ Cθ_n^2( N_0)^1/2.

Some simple algebra shows that N_0 ≤ Cn^2h_n, thus by choosing a sufficiently small c in the definition of h_n in (<ref>), we have 

    (_0, _1) ≤ Cθ_n^2nh_n^1/2≤ c.


To complete the proof, we now show that (_1,_1) ≤ c. Consider an arbitrary realization of {X_i}_i=1^n, and assume that based on their location indices {b_i}_i=1^n, {X_i}_i=1^n is partitioned into L clusters with corresponding cardinality s_ℓ so that the X_i's in the same cluster have the same value b_i. Apparently, we have the relations 1≤ L≤ n and ∑_ℓ=1^L s_ℓ = n. Let m_max be the maximum cluster size, and define the “good event" Ω_n{m_max≤ K}, where K2/(1-4α)) + 2. Then, it holds that

    (_1,_1)    = *1_Ω_n(_1(y|x), _1(y|x) + *1_Ω_n^c(_1(y|x),_1(y|x))
       ≤*1_Ω_n(_1(y|x), _1(y|x) + (Ω_n^c).

Under the choice of h_n in (<ref>), N is of the order n^2/(4α + 1), and

    λ_K lim_n→∞n^K/K!N^K-1 = 0.

Thus by Lemma <ref> (and continuity), it holds that Ω_n has asymptotic probability 1 under both _1 and _1. As a result, it suffices to upper bound (_1(y|x), _1(y|x)) for each realization x in Ω_n, where the maximum cluster size m_max is bounded by a fixed constant.

Denoting p_1,π_ℓ and p_1,π_ℓ for each ℓ∈[L] as the joint density of those y_i's in the ℓth cluster π_ℓ conditioning on the given realization of {X_i}_i=1^n under _1 and _1, we obtain that 

    p_1(y|x) - p_1(y|x) = ∏_ℓ=1^L p_1,π_ℓ - ∏_ℓ=1^L p_1,π_ℓ.

The above inequality further implies by telescoping that

    *p_1(y|x) - p_1(y|x)≤∑_ℓ=1^L*p_1,π_ℓ - p_1,π_ℓ.

For each ℓ∈[L], *p_1,π_ℓ - p_1,π_ℓ only depends on the ℓth cluster through its cardinality, which we now control for a general cluster size d≥ 1. Without loss of generality, we assume that ℓ=1 and the y_i's in this cluster are {y_1,…, y_d} with common location index b_i = 1 for i∈[d]. Then, under the choice of θ_n^2 in (<ref>), we clearly have Y_i  = θ_nr_1 + ε_i under _1 and Y_i = θ_nr_1 + ε_i under _1 for i∈[d], where the sequence {ε_i}_i=1^d follows the standard normal distribution under both _1 and _1. Therefore it holds that

    p_1,π_1(y_1,…,y_d)    = ∫_-∞^∞φ(y_1-θ_nv)…φ(y_d-θ_nv)(dv), 
    p_1,π_1(y_1,…,y_d)     =∫_-∞^∞φ(y_1-θ_nv)…φ(y_d-θ_nv)φ(v)dv,
 
where  is the distribution of {r_i}_i=1^N specified in (<ref>). Using the well-known equality φ(t-θ_nv) = φ(t)(∑_k=0^∞ v^kθ_n^kH_k(t)/k!) for any t,v, where H_k is the kth order Hermite polynomial, it holds that

    φ(y_1-θ_nv)…φ(y_d-θ_nv)
       = φ(y_1)…φ(y_d)∑_k_1,…,k_d=0^∞ v^∑_i=1^d k_iθ_n^∑_i=1^d k_iH_k_1(y_1)/k_1!…H_k_d(y_d)/k_d!
       = φ(y_1)…φ(y_d)∑_k=0^∞ v^kθ_n^k∑_k_1+…+k_d=kH_k_1(y_1)/k_1!…H_k_d(y_d)/k_d!

and therefore

    p_1,π_1(y_1,…,y_d) - p_1,π_1(y_1,…,y_d)
       = φ(y_1)…φ(y_d)∑_k=0^∞θ_n^k∑_k_1+…+k_d =kH_k_1(y_1)/k_1!…H_k_d(y_d)/k_d!∫ v^k ( - Φ)(dv)
       =φ(y_1)…φ(y_d)∑_k=p^∞θ_n^2k∑_k_1+…+k_d =2kH_k_1(y_1)/k_1!…H_k_d(y_d)/k_d!∫ v^2k ( - Φ)(dv),

where the second equality follows by the symmetry and moment matching property of  in (<ref>) and p  (q+1)/2 is a positive integer. This further yields

    *p_1,π_1(y_1,…,y_d) - p_1,π_1(y_1,…,y_d)
       ≤φ(y_1)…φ(y_d)∑_k=p^∞θ_n^2k∑_k_1+…+k_d=2k*H_k_1(y_1)/k_1!…*H_k_d(y_d)/k_d!∫ v^2k(dv) +
       φ(y_1)…φ(y_d)∑_k=p^∞θ_n^2k∑_k_1+…+k_d=2k*H_k_1(y_1)/k_1!…*H_k_d(y_d)/k_d!∫ v^2kφ(v) dv
        I + II.

For term I, since  is compactly supported on [-B,B], one clearly has

    I ≤φ(y_1)…φ(y_d)∑_k=p^∞θ_n^2kB^2k∑_k_1+…+k_d=2k*H_k_1(y_1)/k_1!…*H_k_d(y_d)/k_d!.

For term II, using the equality ∫φ(v)v^2kdv = (2k-1)!!, with (2k-1)!! (2k-1)(2k-3)… 1, we obtain

    II = φ(y_1)…φ(y_d)∑_k=p^∞θ_n^2k(2k-1)!!∑_k_1+…+k_d=2k*H_k_1(y_1)/k_1!…*H_k_d(y_d)/k_d!.

We now upper bound ∫_-∞^∞*H_k(t)φ(t)dt for an arbitrary positive integer k. When k is even, as has been calculated in <cit.> (cf. chain of inequality after Equation (19) on Page 662), ∫_-∞^∞*H_k(t)φ(t)dt ≤ 2^k/2(k-1)!!. When k is odd, set k = 2k + 1, then we have

    ∫_-∞^∞*H_k(t)φ(t)dt    = ∫_-∞^∞φ(t)*(2k+1)!∑_m=0^k(-1)^mt^2k+1-2m/m!(2k+1-2m)!2^mdt
       ≤∑_m=0^k(2k+1)!/m!(2k+1-2m)!2^m∫_-∞^∞ |t|^2k+1-2mφ(t)dt
       = √(2/π)∑_m=0^k(2k+1)!(2k-2m)!!/m!(2k+1-2m)!2^m
       = √(2/π)∑_m=0^k(2k+1)!(2m)!!/(k-m)!(2m+1)!2^k-m
       = √(2/π)(2k+1)!!∑_m=0^kk!/m!(k-m)!(m!)^22^2m/(2m+1)!
       ≤ (2k+1)!!∑_m=0^kk!/m!(k-m)!
       = (2k+1)!!2^k,

where in the third line we use the fact that ∫_-∞^∞ |t|^2ℓ+1φ(t)dt = √(2/π)(2ℓ)!! for any positive integer ℓ. Define for any positive integer k: [k]_1  k-1 if k is even and k if k is odd, and [k]_2  k/2 if k is even and (k-1)/2 if k is odd. Then, the above calculation implies that ∫_-∞^∞*H_k(t)φ(t)dt ≤ ([k]_1)!!2^[k]_2 for any k, and moreover, it can be readily checked that ([k]_1)!!/(k!) = 1/(2^[k]_2([k]_2)!). Therefore, for term I = I(y_1,…,y_d), we have

    ∫_^dI(y_1,…,y_d)dy_1… dy_d
       ≤∑_k=p^∞θ_n^2k(B^2)^k∑_k_1+…+k_d=2k1/(k_1)!…(k_d)!([k_1]_1)!!2^[k_1]_2…([k_d]_1)!!2^[k_d]_2
       = ∑_k=p^∞θ_n^2k(B^2)^k∑_k_1+…+k_d=2k1/([k_1]_2)!…([k_d]_2)!.

Now note that the number of d-tuple (k_1,…,k_d) such that k_1+…+k_d = 2k is upper bounded by (Ck)^d, which is further bounded by C^k for every k≥ 0 with some sufficiently large C that only depends on d, and for each such tuple, it holds that 

    k-d/2= ∑_i=1^d k_i-1/2≤∑_i=1^d [k_i]_2≤∑_i=1^d k_i/2 = k,

thus we have

    ∑_k_1+…+k_d=2k{([k_1]_2)!…([k_d]_2)!}^-1≤ C^k∑_k-d/2≤k̅_1+…+k̅_d≤ k{(k̅_1)!…(k̅_d)!}^-1.

For the latter quantity, we have by the multinomial identity

    ∑_x_1+…+x_d+1=kk!/(x_1!… x_d+1!)(d+1)^-k = 1

that

    (d+1)^k/k!   = ∑_k̅_1+…+k̅_d+1 = k1/(k̅_1)!…(k̅_d+1)!
        = ∑_k̅_1+…+k̅_d ≤ k1/(k̅_1)!…(k̅_d)!(k-(k̅_1+…+k̅_d))!
       ≥∑_k-d/2≤k̅_1+…+k̅_d ≤ k1/(k̅_1)!…(k̅_d)!(k-(k̅_1+…+k̅_d))!
       ≥**d/2!^-1∑_k-d/2≤k̅_1+…+k̅_d ≤ k1/(k̅_1)!…(k̅_d)!.

This concludes that

    ∫_^dI(y_1,…,y_d)dy_1… dy_d ≤θ_n^2p∑_k=p^∞(CB^2)^k/k!≤θ_n^2pe^CB^2.

Using a similar argument for II = II(y_1,…,y_d), we obtain

    ∫_^d II(y_1,…,y_d)dy_1… dy_d ≤∑_k=p^∞(2k-1)!!/k!θ_n^2kC^k     = ∑_k=p^∞(2k-1)!!/(2k)!!θ_n^2k(2C)^k  ≤θ_n^2pC^p

since θ_n^2 < 1/C for sufficiently large n. 

Putting together the pieces, we have for every realization x in Ω_n

    ∫_^n*p_1(y|x) - p_1(y|x) dy   ≤∑_ℓ=1^L∫_^|π_ℓ|*p_1,π_ℓ - p_1,π_ℓ≤ Lmax_1≤ d≤ Kθ_n^2p(e^CB^2 + C^p) 
       ≤ nθ_n^2p(e^CB^2 + C^p) ≤ c.

Here, the second inequality follows since every *p_1,π_ℓ - p_1,π_ℓ depends on the ℓth cluster only through its cardinality, the third inequality follows since L≤ n and K is a fixed absolute constant that only depends on α, and the last inequality follows due to the choice θ_n^2= h_n^2α = cn^-4α/(4α+1) and the value of p. This completes the proof.




For any fixed positive integer q, there exist a B < ∞ and a symmetric distribution  on [-B, B] such that  and the standard normal distribution have the same first q moments, that is,

    ∫_-B^B x^j(dx) = ∫_-∞^∞ x^jφ(x)dx,    j=1,…, q.





If μ∈^d and Σ_1 and Σ_2 are positive definite d× d matrices, then

    1/100≤TV*N_d(μ, Σ_1), N_d(μ, Σ_2)/min{1, Σ_1^-1Σ_2 - 𝐈_d_F}≤3/2.
 





For the following lemma, we first introduce some terminology regarding the multinomial distribution. Let m, M be two positive integers, and the random vector (f_1,…, f_M) be the multinomial count with total count m and equal probability (1/M,1/M,…, 1/M). Define ρ m/M. For any positive integer r≥ 2, define 
λλ_rlim_m→∞ m^r/(r!M^r-1). Following <cit.> (Chapter 2, Equation (11)), we will call the domain of variation m,M→∞, in which 

    ρ→ 0,    0 < λ_r < ∞

the left-hand r-domain. The following lemma characterizes the asymptotic behavior of the maximum frequency f_max defined as max_1≤ j≤ Mf_j.




Suppose the multinomial distribution with total count m and equal probability (1/M,…, 1/M) is in the left-hand r-domain for some positive integer r≥ 2 with limit λ_r, then it holds that

    (f_max = r - 1)→ e^-λ_r and   (f_max = r)→ 1 - e^-λ_r,

i.e., the maximum frequency converges asymptotically to a two-point distribution.


















apalike











§ APPENDIX


Throughout the supplement, we continue to use the notation introduced in the main paper. We also use the following new notation. For any positive integer d≥ 2, 𝕊^d-1 stands for the unit Euclidean sphere in ^d. For a real vector x, define x_0 as the number of nonzero coordinates of x. When writing the Hölder class Λ_α,I(C), we will omit the domain I in the subscript for simplicity. For two distributions  and  on , we will write * as their convolution.






§ PROOFS OF RESULTS IN SECTION <REF>




 §.§ Proof of Theorem <ref>




























Throughout the proof, we will use C, c to denote two generic fixed positive constants that only depend on M_K, M_K,α, C_F, C_σ, C_ε, C_0, c_0. C and c might have different values at each occurrence. We also use the notation W_ij W_i - W_j for a generic random variable W. 

Denote the two U-statistics on the numerator and denominator of σ̂^2 respectively as U_1, U_2, with corresponding mean values θ_1,θ_2. That is, with i≠ j,

    θ_1 *K_h(X_i-X_j)(Y_i-Y_j)^2/2   and   θ_2 K_h(X_i-X_j).

Define the “good" event E{U_2 ≥θ_2/2} and E^c as its complement, then it holds that

    *σ̂^2 - σ^2^2 = **U_1 - U_2σ^2/U_2^21{E} + **U_1 - U_2σ^2/U_2^21{E^c}.

By definition of E, the first term satisfies that

    **U_1 - U_2σ^2/U_2^21{E}≤4/θ_2^2*U_1 - U_2σ^2^2.

For θ_2, we have

    θ_2    =  K_h(X_i-X_j) = ∫1/hK*v/hp_(v)dv = ∫_-1^1 K(u)p_(uh)du
       ≥∫_U_h K(u)p_(uh)du≥inf_u∈U_hp_(uh)inf_u∈[-1,1]K(u)λ(U_h) ≥M_Kc_0^2.

Here, the third equality follows from the fact that K(·) is supported in [-1,1], and U_h starting from the first inequality is defined in Condition (c) in P_vf,(X,ε) (note that for any fixed δ_0 > 0 given therein, h ≤δ_0 for sufficiently large n). Moreover, it holds that

    *U_1 - U_2σ^2^2 ≤ 3**U_1 - θ_1^2 + σ^4*U_2 - θ_2^2 + *θ_1 - θ_2σ^2^2.

By Lemmas <ref> and <ref> and the fact that σ^4 ≤ C_σ^2, we have

    *U_1 - θ_1^2 + σ^4*U_2 - θ_2^2 ≤ C(n^-1 + n^-2h^-1).

For the third term *θ_1 - θ_2σ^2^2, we have 

    θ_1 = *K_h(X_i-X_j)(Y_i-Y_j)^2/2 = *K_h(X_i-X_j)(f(X_i)-f(X_j))^2/2 + θ_2σ^2
 
and

    *K_h(X_i-X_j)(f(X_i)-f(X_j))^2/2≤ C*1/hK*/h*^2(α∧ 1)
       = C∫1/hK*u/h*u^2(α∧ 1)p_(u)du = C∫ K(v)h^2(α∧ 1)|v|^2(α∧ 1)p_(vh)dv
       ≤ Ch^2(α∧ 1)sup_u∈p_(u)∫ K(v)|v|^2(α∧ 1)dv ≤ Ch^2(α∧ 1).

Here, the first inequality follows since f∈Λ_α(C_F), and the last inequality follows from Condition (b) in P_vf,(X,ε) and the convolution formula. 



Putting together the pieces, the choice of h in (<ref>) in the main paper yields

    **U_1 - U_2σ^2/U_2^21{E}≤ C(h^4(α∧ 1) + n^-1 + n^-2h^-1) ≤ C(n^-8α/(4α +1)∨ n^-1).

For the second term in (<ref>), we have 

    **U_1 - U_2σ^2/U_2^21{E^c}≤ 2σ^4*E^c + 2**U_1/U_2^21{E^c}.

Direct calculation shows that

    *U_1/U_2^2    = ∑_i<j,i^' < j^' K_h*X_ijK_h*X_i^' j^'(Y_i - Y_j)^2(Y_i^'-Y_j^')^2/4/∑_i<j, i^' < j^' K_h*X_ijK_h*X_i^' j^'
       ≤∑_i<j,i^' < j^' K_h*X_ijK_h*X_i^' j^'*(f(X_i) - f(X_j))^2 + σ^2ε_ij^2*(f(X_i^') - f(X_j^'))^2 + σ^2ε_i^' j^'^2/∑_i<j,i^' < j^' K_h*X_ijK_h*X_i^' j^'
       ≤ C∑_i<j,i^' < j^' K_h*X_ijK_h*X_i^' j^'**X_ij^2(α∧ 1) + σ^2ε_ij^2**X_i^' j^'^2(α∧ 1) + σ^2ε_i^' j^'^2/∑_i<j, i^' < j^' K_h*X_ijK_h*X_i^' j^'
       ≤ C*h^4(α∧ 1) + σ^2h^2(α∧ 1)max_i<jε_ij^2 + σ^4max_i<j,i^' < j^'ε_ij^2ε_i^' j^'^2,

where the last inequality follows by the support of K(·). By the condition σ^2≤ C_σ and the independence of {ε_i}_i=1^n and 1{E}, this implies that 

    **U_1 - U_2σ^2/U_2^21{E^c}≤ C*E^c*1 + max_i<jε_ij^2 + max_i<j,i^' <j^'ε_ij^2ε_i^' j^'^2.

Applying the first part of Lemma <ref> with v = nθ_2^2/16, u = n^2hθ_2^2/16 with the condition h = Ω(n^-(2-δ)) being satisfied with δ = 8α/(4α+1) and δ = 1 for α≤ 1/4 and α > 1/4 respectively, it holds that

    *E^c = (|U_2 - θ_2|≥θ_2/2) ≤ C*exp(-θ_2^2n/16) + exp(-θ_2^2n^2h/16).

Moreover, by Condition (d) in P_vf,(X,ε), there exists some fixed positive constant η such that

    1 + max_i<jε_ij^2 + max_i<j,i^' <j^'ε_ij^2ε_i^' j^'^2 ≤ Cn^η.

Putting together the pieces and using the fact that n^2h→∞ as n→∞, it yields

    **U_1 - U_2σ^2/U_2^21{E^c} = o(n^-8α/(4α +1)∨ n^-1).


This completes the proof.
 



 §.§ Supporting lemmas



Suppose f∈Λ_α(C_F) and σ^2≤ C_σ for some fixed constants C_F,C_σ and the joint distribution of (X,ε) satisfies Conditions (a), (b) and (d) in P_cv,(X,ε) with constants C_0,C_ε. Then, the U-statistic U_1 defined in the proof of Theorem <ref> satisfies 

    *U_1 - θ_1^2 ≤ C*n^-1∨ n^-2h^-1,

where C is some fixed positive constant that only depends on M_K, M_K, α, C_F, C_σ, C_ε, C_0.


Denote g as the kernel of U_1, that is,

    g(D_i, D_j)  K_h(X_i - X_j)(Y_i - Y_j)^2/2,   D_i  (X_i, ε_i)^⊤.

Recall that θ_1 =  g(D_i, D_j) for i≠ j. Then, it holds that 

    *U_1 - θ_1^2 = n 2^-2∑_i<j,i^'<j^'**g(D_i, D_j) - θ_1*g(D_i^', D_j^') - θ_1.

When i,j,i^',j^' take four different values, the expectation is zero. When they take three values, say, i = i^' < j < j^', by writing _ε as the conditional expectation given {X_i}_i=1^n, we have

    *g(D_i,D_j)g(D_i, D_j^')
       = 1/4*1/h^2K*X_i - X_j/hK*X_i - X_j^'/h(Y_i - Y_j)^2(Y_i - Y_j^')^2
       ≲*1/h^2K*X_i - X_j/hK*X_i - X_j^'/h((f(X_i) - f(X_j))^2 + σ^2ε_ij^2)((f(X_i) - f(X_j^'))^2 + σ^2ε_ij^'^2)
       ≲*1/h^2K*X_i - X_j/hK*X_i - X_j^'/h_ε***^2(α∧ 1) + σ^2ε_ij^2**X_ij^'^2(α∧ 1) + σ^2ε_ij^'^2
       ≲*1/h^2K*X_i - X_j/hK*X_i - X_j^'/h**X_ij^'^2(α∧ 1) + 2σ^2**^2(α∧1)+*X_ij^'^2(α∧ 1) + ε^4·σ^4
       ≲*1/h^2K*X_i - X_j/hK*X_i - X_j^'/h
       = ∫ K(v)K(w)p_X(u)p_X(u +hv)p_X(u + hw)dudvdw ≲ 1.

In the last line, we invoke Conditions (b) and (d) in P_cv,(X,ε). Moreover, it can be readily calculated that θ_1 = O(1). This concludes that the summand in (<ref>) is bounded by a fixed constant when i,j,i^',j^' take three different values. Lastly, performing a similar analysis, we obtain that **g(D_i, D_j) - θ_1*g(D_i^', D_j^') - θ_1 = O(1/h) when i = i^' and j = j^'. We therefore conclude that

    (U_1) ≲n^3 + n^2h^-1/n^4≍ n^-1 + n^-2h^-1.

This completes the proof.




Suppose h_n ≳ n^-(2-δ) for some 0 < δ < 2. Then, assuming Condition (b) in P_cv,(X,ε) with constant C_0, the U-statistic U_2 defined in the proof of Theorem <ref> satisfies 

    **U_2 - θ_2≥ C(v^1/2n^-1/2 + u^1/2n^-1h^-1/2)≤ C(exp(-u) + exp(-v))

for any u, v> 0, and

    *U_2 - θ_2^2 ≤ C(n^-1∨ n^-2h^-1),

where C is some fixed positive constant that only depends on M_K, M_K,α, C_0.


We first prove the concentration inequality by upper bounding the 5 quantities in Lemma <ref>. Denote g as the kernel of U_2 and g_1 as its linear part, that is, for some  i≠ j,

    g_1(X_i) *g(X_i, X_j)| X_i*K_h(X_i - X_j)| X_i.
  
For B_1, we have

    g_1(X_i) = ∫1/hK*u-X_i/hp_X(u)du = ∫ K(u)p_X(uh + X_i)du ≲ 1

due to Condition (b) in P_cv,(X,ε). Thus it also holds ν^2_1 ≲ 1. For B_2, we have

    B_2^2    = nsup_X_i*g^2(X_i, X_j)| X_i = nsup_X_i∫1/h^2K^2*u-X_i/hp_X(u)du
       ≲n/hsup_X_i∫ K(u)p_X(uh+X_i)du ≲ nh^-1,

where in the first inequality we use the condition that K(·) is bounded by M_K. Moreover, we clearly have B_3≲ h^-1. Lastly, for ν_2^2, it holds that 

    ν_2^2 = ∫1/h^2K^2*u/hp_(u)du ≲1/h∫ K(u)p_(uh)du ≲1/h,

where the last inequality follows by Condition (b) in P_cv,(X,ε) and the convolution formula

    p_(u)= ∫ p_X(t)p_X(t-u)dt ≤sup_u∈p_X(u)∫ p_X(t)dt = sup_u∈p_X(u).

Therefore, Lemma <ref> yields that

    (|U_2 - θ_2| ≥ a_1v^1/2 + a_2v + b_1u^1/2 + b_2u + b_3u^3/2 + b_4u^2) ≤ C(exp(-v) + exp(-u)),

where a_1 ≲ n^-1/2, a_2 ≲ n^-1, b_1 ≲ n^-1h^-1/2, b_2≲ n^-1, b_3 ≲ n^-3/2h^-1/2, b_4 ≲ n^-2h^-1. Under the condition that h ≳ n^-(2-δ) for some δ > 0 and n is sufficiently large, the dominant terms in the above inequality are a_1 and b_1, that is,

    n^-1/2∨ n^-1h^-1/2.

This proves the first part of the theorem. The expectation version follows by Lemma <ref>.




Suppose a random variable X satisfies the tail condition (|X|≥ a_1t^1/2 + a_2t + a_3t^3/2 + a_4t^2) ≤ C_1exp(-C_2t) for any t > 0 and some positive constants a_1,a_2,a_3,a_4, C_1,C_2. Then, for any positive integer p, it holds that

    (|X|^p)^1/p≤ C_3p^1/p(a_1 + a_2 + a_3 + a_4)

for some positive constant C_3 that only depends on C_1,C_2.


We use C_3 to denote a positive constant that only depends on C_1 and C_2, which might have different values at each occurrence. The tail condition in the assumption is equivalent to

    (|X| ≥ t) ≤ C_3exp*-C_3*t^2/a_1^2∧t/a_2∧t^2/3/a_3^2/3∧t^1/2/a_4^1/2.
 
Let I_1-I_4 be a partition of (0,+∞) such that for t∈ I_1, t^2/a_1^2 = min*t^2/a_1^2∧t/a_2∧t^2/3/a_3^2/3∧t^1/2/a_4^1/2, and similarly for I_2,I_3, I_4. Then, we have

    (|X|^p)
       = ∫_0^∞ pt^p-1(|X|≥ t)dt
       ≤ C_3{∫_I_1 pt^p-1exp*-C_3t^2/a_1^2dt + ∫_I_2 pt^p-1exp*-C_3t/a_2dt + ∫_I_3 pt^p-1exp*-C_3t^2/3/a_3^2/3dt +
       ∫_I_4 pt^p-1exp*-C_3t^1/2/a_4^1/2dt}
       ≤ C_3{∫_0^∞ pt^p-1exp*-C_3t^2/a_1^2dt + ∫_0^∞ pt^p-1exp*-C_3t/a_2dt + ∫_0^∞ pt^p-1exp*-C_3t^2/3/a_3^2/3dt + 
       ∫_0^∞ pt^p-1exp*-C_3t^1/2/a_4^1/2dt}
       = C_3p{a_1^p∫_0^∞ t^p-1exp*-C_3t^2dt + a_2^p∫_0^∞ t^p-1exp*-C_3tdt + a_3^p∫_0^∞ t^p-1exp*-C_3t^2/3dt + 
        a_4^p∫_0^∞ t^p-1exp*-C_3t^1/2dt}
       ≤ C_3p(a_1^p + a_2^p + a_3^p + a_4^p).

This completes the proof.




Recall the Condition (c) in the definition of P_cv,(X,ε) in the main paper, and Condition (c^') in the subsequent paragraph. We have (c^')⇒ (c).


Choose some δ_0 < 1/8 and fix any δ≤δ_0 and u∈[-1,1]. By the convolution formula, we have

    p_(uδ) = ∫_S p_X(s)p_X(s - uδ)ds ≥ c_0∫_S p_X(s - uδ)ds.

Therefore, it suffices to show that λ**S - uδ⋂[0,1]⋂ S is lower bounded by some fixed constant, say, 1/8. Assume this does not hold, then we have

    1 = λ([0,1])≥λ**S - uδ⋂[0,1] + λ(S) - 1/8 ≥λ(S) - δ + λ(S) - 1/8 ≥ 5/4,

which is a contradiction. 




The covariance Σ_1 defined in the proof of Theorem <ref> in the main paper is positive definite.

 
We will prove that for any a∈^n with a = 1, it holds that a^⊤Σ_1 a > 0. For each realization of {X_i}_i=1^n, partition the set [n] into L clusters for some positive integer 1≤ L ≤ n such that the Y_i's in each cluster fall into the same trapezoid in the lower bound construction. Denote these L clusters as M_1,…,M_L. Then, it follows that

    *∑_i=1^n a_iY_i|{X_i}_i=1^n = ∑_ℓ=1^L*∑_i∈M_ℓ a_iY_i|{X_i}_i=1^n.

Since a = 1, there exists some ℓ_0∈[L] such that A∑_i∈M_ℓ_0a_i^2 > 0. Partition {i:i∈M_ℓ_0} according to the sign:

    A_+ {i∈M_ℓ_0: a_i≥ 0},   A_-{i∈M_ℓ_0: a_i<0}

and define S_+ ∑_i∈A_+ a_i and A_+ ∑_i∈A_+ a_i^2, and S_- and A_- similarly. Then, A = A_+ + A_-. Moreover, it holds that

    *∑_i∈M_ℓ_0 a_iY_i|{X_i}_i=1^n
       = (1+h_n^2α)∑_i∈M_ℓ_0a_i^2 + h_n^2α*∑_i,j∈A_+;i≠ ja_ia_j + ∑_i,j∈A_-;i≠ ja_ia_j + 2∑_i∈A_+,j∈A_-a_ia_j
       = A(1+h_n^2α) + h_n^2α*S_+^2 - A_+ + S_-^2  - A_- + 2S_+S_-
       =A + h_n^2α*S_+ + S_-^2 ≥ A > 0.

This completes the proof.




Let Z_1,…,Z_n,Z∈Z be i.i.d., and g:Z^2→ be a symmetric measurable function with *g(Z_1,Z_2) < ∞. Write U_n(g) ∑_i<jg(Z_i,Z_j) and g_1(z) *g(Z,z). Define

    B_1 sup_Z_2**g(Z_1,Z_2)| Z_2,    B_2 *nsup_Z_2*g^2(Z_1,Z_2)| Z_2^1/2,
       B_3 g_∞

and

    ν_1^2 *g_1^2(Z_2),   ν_2^2 *g^2(Z_1,Z_2).
 
Then, it holds that

    **U_n(g) - *U_n(g)≥ t + C_1nν_2 u^1/2 + C_2nB_1u + C_3B_2u^3/2 + C_4B_3u^2
       ≤ 2exp*-t^2/n^2/8nν_1^2 + 4B_1· t/n + C_5e^-u,

where C_1-C_5 are absolute constants.















§ PROOFS OF RESULTS IN SECTION <REF>




 §.§ Proof of Theorem <ref>


Throughout the proof, C and c will denote two generic positive constants that do not depend on n and might have different values at each occurrence. We only prove the case for the pointwise error and the result for the integrated error will follow. Consider a fixed x^*∈(X). We will continue to use the notation ℓ,q(·),𝐁_n, X_ij,K_ij introduced in Section <ref> in the main paper. We will drop the subscript in V_(x^*) for notational simplicity. Recall the choice of (h_1,h_2) in (<ref>) in the main paper.



















Define 𝐁𝐁_n and the good event Ω_n {*𝐁_n - 𝐁≤ 1/(2𝐁^-1)}. Note that Ω_n is well-defined as we now prove 𝐁 is indeed invertible. For any a∈𝕊^ℓ, it holds that

    a^⊤𝐁a
       = ∫∫*a^⊤q*(u+v)/2 - x^*/h_2^2 1/h_1K*u-v/h_11/h_2K*(u+v)/2 - x^*/h_2p_X(u)p_X(v)dudv
       = ∫∫*a^⊤q(v)^2 K(u)K(v)p_X(x^* + h_2v + h_1u/2)p_X(x^* + h_2v - h_1u/2)dudv
       = ∫∫*a^⊤q(v - h_1u/(2h_2))^2K(u)K(v - h_1u/(2h_2))p_X(x^* + h_2v)p_X(x^* + h_2v - h_1u)dudv
       = ∫_-1^1∫_-1+h_1u/(2h_2)^1+h_1u/(2h_2)*a^⊤q*v - h_1u/2h_2^2K(u)K*v - h_1u/2h_2p_X(x^* + h_2v)p_X(x^* + h_2v - h_1u)dvdu
       ≥M_K^2∫_-1^1∫_-1+h_1u/(2h_2)^1+h_1u/(2h_2)*a^⊤q*v - h_1u/2h_2^2p_X(x^* + h_2v)p_X(x^* + h_2v - h_1u)dvdu,

where in the last inequality we use the lower bound M_K on K(·). Note that the first term of the integrand *a^⊤q*v - h_1u/2h_2^2 is a polynomial of variables u,v and thus only takes zero value with Lebesgue measure at most 0. By the second part of Condition (c) in P_vf,(X,ε) with δ = h_2 (note that h_2≤δ_0 for any fixed δ_0 >0 and sufficiently large n), for the given x^*, there exists a set A_x^*⊂ [-1,1] with Lebesgue measure at least c_0 such that for all v∈A_x^*, x^* + h_2v∈(X), and moreover, for each v∈A_x^*, the second part of Condition (c) with δ = h_1 again implies the existence of a set A_x^*,v⊂[-1,1] with Lebesgue measure at least c_0 such that for all u∈A_x^*,v, it holds that x^* + h_2v - h_1u∈(X). Therefore, by the first part of Condition (c)  and the fact that h_1/h_2→ 0, there exists a fixed positive constant c such that

    λ_mininf_a∈𝕊^ℓa^⊤𝐁a≥ c > 0.

This concludes that Ω_n is well-defined. By triangle inequality, we have

    *V(x^*) - V(x^*)^2
       ≲*V(x^*) - V(x^*)^2 +*V(x^*)^2
       ≲**(V(x^*) - V(x^*))^2 + *V(x^*) - *V(x^*)^2 + 
       **V^2(x^*) + V^2(x^*).

By Lemma <ref>, we have for the first term

    **(V(x^*) - V(x^*))^2 ≤ C*h_1^4(α∧ 1) + h_2^2β + h_1^2(β∧ 1) + τ_n^2.
 
By Lemma <ref> with conditions nh_2→∞ and n^2h_1h_2→∞ satisfied with the choices of (h_1,h_2) in (<ref>), we have for the second term

    *V(x^*) - *V(x^*)^2 ≤ C*n^-1h_2^-1 + n^-2(h_1h_2)^-1 + τ_n^2.

Plugging in the values of (h_1,h_2) as in (<ref>) and choosing τ_n ≍ n^-κ for some fixed κ≥ 1, we obtain that

    *V(x^*) - V(x^*)^2 ≤ C(n^-8αβ/4αβ + 2α + β + n^-2β/2β + 1) + **V^2(x^*) + V^2(x^*).

Lastly, note that V(x^*) = ∑_i<jD_ijw_ij / (∑_i<jw_ij + τ_n) is a linear estimator with weight w_ij = n 2^-1q^⊤(0)𝐁_n^*q((X_ij-x^*)/h_2)K_ij. By definition of 𝐁_n^*, w_ij is thus a weighted polynomial of K_ij(X_ij-x^*)/h_2 up to some order that only depends on ℓ. Therefore, in view of the choice of τ_n (decaying to 0 polynomially with n), h_1,h_2, there exists some sufficiently large constant η (only depending on α,β,κ) such that

    **V^2(x^*)|{X_i}_i=1^n^2 ≲  n^η,

and thus by Cauchy-Schwarz and the exponential inequality in Lemma <ref>, it holds that

    **V^2(x^*) + V^2(x^*)   = *(V^2(x^*)|{X_i}_i=1^n + V^2(x^*)) 
       ≲ ((V^2(x^*)|{X_i}_i=1^n)^2 + V^4(x^*))^1/2^1/2(Ω_n^c) 
       = o*n^-8αβ/4αβ + 2α + β + n^-2β/2β + 1.

This completes the proof.
























































































 §.§ Proof of Theorem <ref>


Note that the boundary of n^-8αβ/(4αβ + β+ 2α) and n^-2β/(2β + 1) lies at α = β/(4β + 2). When α≥β/(4β + 2), the statement can be proved using a slight variation of the proof of Theorem 4.2 in <cit.>, and we omit the details here. Next, we will focus on the case where α < β/(4β + 2). Consider a fixed point x^*∈(X). Throughout the proof, C and c represent two generic positive constants which only depend on α,β,C_F,C_V, C_σ,C_0,c_0,C_ε 

and might have different values at each occurrence, but like in the proof of Theorem <ref>, let c be always smaller than 1/4. Also, without loss of generality, assume that the sample size n and C_F,C_V, C_σ,C_ε,C_0 are sufficiently large, c_0 is sufficiently small, and [0,1]⊂ I.

We will make use of Le Cam's two point method. Introduce the constants

    θ_n^2  h_1^2α h_2^β cn^-4αβ/4αβ + β + 2α,    M h_2/(4h_1) - 1/2,    N  2M + 1 = h_2/(2h_1),

where we tune the constant c in θ_n^2 so that M is a positive integer. Note that under the above choice, h_2/h_1→∞ as n→∞. We now specify f(·), V(·), distribution of X and distribution of ε in the null and alternative hypotheses, H_0 and H_1, respectively.



   Choice of ε: Under both H_0 and H_1, let ε∼N(0,1).

   Choice of V(·): Under H_0, let V≡ 1. Under H_1, let V = 1 - θ_n^2H((x - x^*)/h_2), where H(·) is β-Hölder smooth, infinitely differentiable, compactly supported on [-2,2], and takes value 1 on [-1,1].

   Choice of f(·): Under H_0, let f≡ 0. Under H_1, let f be zero outside [x^* - h_2, x^* + h_2], and inside this interval, the linear interpolant of the function that takes value r_i on [x^* - h_2 + (4i-3)h_1, x^* - h_2 +  (4i-1)h_1] and zero at x^* - h_2 + 4(i-1)h_1 for all i∈[N], where {r_i}_i=1^N is an i.i.d. sequence of symmetric and compactly supported random variables with distribution  satisfying

    ∫_-∞^∞ x^j(dx) = ∫_-∞^∞ x^jφ(x)dx,    j=1,…,q,
 
where q is some fixed odd integer strictly larger than 1 + (β+2α)/(2αβ).

   Choice of X: Under both H_0 and H_1, let X be uniformly distributed on the union of the intervals 

    [0, 1] ⋂*[0,x^* - 2h_2] ⋃ [x^* + 2h_2, 1] ⋃_i=1^N [x^* - h_2 + (4i-3)h_1, x^* - h_2 +  (4i-1)h_1].


See Figure <ref> for an illustration. We now make a few remarks about the above construction. For the design of V(·) under H_1, one example of the smooth bump function H(·) is (1_[-3/2,3/2]∗φ_1/2)(·), where φ_ε(x)φ(x/ε)/ε with φ(x)exp(-1/(1-x^2))1{|x|≤ 1} being a smooth and compactly supported mollifier. The design of f(·) under H_1 is a “localized" version of f(·) in the proof of Theorem <ref>. The existence of {r_i}_i=1^N is again guaranteed by Lemma <ref>, and their range, which we denote as B, only depends on α and β and is thus fixed. Lastly, we indeed have x^*∈(X) since it is in the (M+1)th interval in the N intervals specified in the support of X. Moreover, under H_1, conditioning on the event that X_i ∈ [x^*-h_2, x^* + h_2] and any realization of {r_i}_i=1^N, f(X_i) is uniformly distributed over {h_1^αr_1,…,h_1^αr_N}. 












Clearly, under both the null and the alternative hypotheses, V(·) is β-Hölder smooth, and under H_1, f(·) is α-Hölder smooth for each realization of {r_i}_i=1^N due to their compact support. Next, we show that the joint distribution of (X,ε) satisfies the three conditions in P_vf,(X,ε). Condition (d) clearly holds and Condition (a) holds with I=[0,1]. Condition (b) holds as well since for any u in the support of X, p_X(u) = 1/(1-3h_2) for x^*∈(0,1) and p_X(u) = 2/(2-3h_2) for x^*∈{0,1}, both of which are smaller than 2 for sufficiently large n. Lastly, for Condition (c), the first part clearly holds since inf_u∈(X) p_X(u) ≥ 1. For the second part, define A_x,δ{u∈[-1,1]: x+δ u∈(X)}. Then, for any x^*∈(0,1) and any 0 < δ≤ 1/2, we have λ(A_x,δ) ≥ 1/2 if x∈(0,x^* - 2h_2]⋃[x^* + 2h_2, 1) and λ(A_x,δ) ≥ 1/4 if x∈⋃_i=1^N [x^* - h_2 + (4i-3)h_1, x^* - h_2 +  (4i-1)h_1]. A similar statement holds for x^*∈{0,1}. We therefore conclude that Condition (c) also holds.


Denote _0 and _1 to be the joint distributions of {X_i,Y_i}_i=1^n under H_0 and H_1, then the pointwise squared distance between _0 and _1 *V_0(x^*) - V_1(x^*)^2 ≍θ_n^4 is the desired minimax rate. Further define _1 as the corresponding joint distributions of {X_i,Y_i}_i=1^n under H_1 with {r_i}_i=1^N replaced by an i.i.d. standard normal sequence {r_i}_i=1^N. Then, following the same line of proof of Theorem <ref>, it suffices to show that (_0, _1) ≤ c and (_1, _1) ≤ c. 

For the first inequality, in view of (<ref>) in the proof of Theorem <ref>, it suffices to upper bound (_0(y|x), _1(y|x)) for each realization {x_i}_i=1^n.










Note that under _0, y|x∼N_n(0, Σ_0), with Σ_0 = 𝐈_n. Denote {b_i}_i=1^n as the location index sequence of {X_i}_i=1^n taking values in {0,1,…, N}, that is, b_i = 0 if X_i∉ [x^* - h_2, x^* + h_2] and b_i = j if X_i∈[x^* - h_2 + (4j-3)h_1, x^* - h_2 +  (4j-1)h_1] for j∈[N]. Then, due to the symmetry of {r_i}_i=1^N, the design of the nonparametric component f, and the fact that H(·) takes value 1 on [-1, 1], it holds that under _1, y|x∼N_n(0, Σ_1), with 

    (Σ_1)_ii = 1 - θ_n^21**X_i-x^*≤ h_2 +  h_1^2α1**X_i-x^*≤ h_2 = 1

and (Σ_1)_ij = h_1^2α1{b_i = b_j, b_i ≥ 1, b_j≥ 1} for i≠ j. Define N_0∑_i≠ j1{b_i = b_j, b_i ≥ 1, b_j≥ 1}. Then, we have by Lemma <ref> that

    (_0(y|x), _1(y|x)) ≤ C*h_1^4αN_0^1/2 = Cθ_n^2N_0^1/2.

Note that N_0 is a random variable that depends on {X_i}_i=1^n, and by (<ref>) in the proof of Theorem <ref>,

    (_0, _1) ≤ Cθ_n^2(N_0^1/2) ≤ Cθ_n^2( N_0)^1/2.

Since direct calculation implies that (N_0) ≤ Cn^2Nh_1^2 = Cn^2h_1h_2, we have (_0, _1) ≲ 1 under the given choice of h_1,h_2 and θ_n. 

Using a conditioning argument, the second part of proving (_1, _1) ≲ 1 follows similarly from that of Theorem <ref> by noting that n^2h_1 →∞ and nh_1→ 0 as n→∞ under the constraint α < β/(4β+2) so that Lemma <ref> can be similarly applied. The proof is complete. 





 §.§ Proof of Theorem <ref>


As in the proof of Theorem <ref>, we focus on the regime α < β/(4β + 2). We will couple the proof of Theorem <ref> with a standard technique via multiple hypotheses in the classic setting of mean function estimation. 

Introduce the following notation:

    θ_n^2  h_1^2α h_2^β  cn^-4αβ/4αβ + β + 2α,    N_2  1/(4h_2),    N_1  h_2/(2h_1),
       and    x^*_i  2h_2 + (i-1)4h_2,    i∈[N_2],

where we tune the constant c in θ_n^2 so that N_1 and N_2 are both positive integers. Note that under the above choice, h_2/h_1→∞ as n→∞. By the renowned Varshamov-Gilbert bound (cf. Lemma 2.8 in <cit.>), there exists a set of length-N_2 binary sequences {Δ_j}_j=0^M with M ≥ 2^N_2/8 such that Δ_0 = 0_N_2 and for any 0 ≤ k < ℓ≤ M, it holds that ρ(Δ_k, Δ_ℓ) ≥ N_2/8, where ρ is the Hamming distance. We now choose a number of M + 1 hypotheses with {Δ_j}_j=0^M satisfying the above property, which we denote as _0,_1,…,_M. We now specify f(·),V(·), distribution of X and distribution of ε under each hypothesis. 



   Choice of ε: Under _0 and _j for all j∈[M], let ε∼N(0,1).

   Choice of V(·): Under _0, let V_0≡ 1. Under _j, let V_j(x)  1 - ∑_i=1^N_2Δ_j,iθ_n^2H*(x-x^*_i)/h_2, where H(·) is infinitely differentiable, compactly supported on [-2,2] and takes value 1 on [-1,1].

   Choice of f(·): Under _0, let f_0≡ 0. Under _j, for all i∈[N_2] such that Δ_j,i = 1, let f be the linear interpolation of the function that takes value r^(j)_i,k on the interval [x^*_i - h_2 + (4k-3)h_1, x^*_i - h_2 + (4k-1)h_1] and value zero at x^*_i - h_2 + 4(k-1)h_1 for all k∈[N_1], where by denoting m_jΔ_j_0, {r^(j)_i,k}_j∈[M], i∈[m_j], k∈[N_1] is an i.i.d. sequence of symmetric and compactly supported random variables with distribution  satisfying

    ∫_-∞^∞ x^j(dx) = ∫_-∞^∞ x^jφ(x)dx,    j=1,…,q,

where q is some fixed odd integer that only depends on α and β.

   Choice of X: Under _0 and _j for all j∈[M], let X be uniformly distributed on the union of the disjoint intervals 

    ⋃_i=1^N_2⋃_k=1^N_1 [x^*_i - h_2 + (4k-3)h_1, x^*_i - h_2 + (4k-1)h_1].


The existence of H(·) in the design of V(·) and variables {r^(j)_i,k}_j∈[M], i∈[m_j], k∈[N_1] is as argued in the proof of Theorem <ref>. Moreover, one can readily check that for each 0≤ k<ℓ≤ M, the integrated squared distance between each _k and _ℓ satisfies

    d(_k, _ℓ)∫*V_k(x) - V_ℓ(x)^2p_X(x)dx ≳ h_2^2β≍ n^-8αβ/4αβ + β + 2α,

which is the desired lower bound.
 


















Clearly, under each _j, 0≤ j≤ M and for each realization of {r^(j)_i,k}_j∈[M], i∈[m_j], k∈[N_1], f_j(·) and V_j(·) are α- and β-Hölder smooth, respectively, due to the compact support of {r^(j)_i,k}_j∈[M], i∈[m_j], k∈[N_1]. Moreover, the joint distribution of (X,ε) (same in all hypothese) satisfies the conditions in P_vf,(X,ε) with a similar argument as in Theorem <ref>.

We now proceed with the proof. Note that under the above design, the support of X is segmented into N_3 N_1× N_2 intervals, and we let {b_i}_i=1^n be the location index of {X_i}_i=1^n, taking values in [N_2]× [N_1], that is, b_i = (k,ℓ) if X_i∈[x^*_k - h_2 + (4ℓ-3)h_1, x^*_k - h_2 + (4ℓ-1) h_1]. As in the proof of Theorem <ref>, define the event Ω_n {max_(k,ℓ)∈[N_2]×[N_1]#*b_i = (k,ℓ)≤ K}, where K is the smallest integer strictly larger than 2β/(β - 4αβ - 2α). Then, by Lemma <ref>, it holds that Ω_n has asymptotic probability 1 under all of _j and _j for 0≤ j≤ M. Now, by a standard reduction scheme with multiple hypotheses (cf. Chapter 2.2 in <cit.>) and Lemma <ref>, it suffices to show that

    1/M∑_j=1^M K(_j, _0; Ω_n) ≤ clog(M)

for some 0 < c < 1/8, where K(, ; Ω_n) is the “conditional" Kullback divergence between probability measures  and  defined as K(, ; E) ∫_Elog(d/d)d for any measurable set E. In order to show (<ref>), it further suffices to show that K(_j, _0; Ω_n) ≤log(M) for all j∈[M]. We now focus on a particular j∈[M]. For notational brevity, we will drop the superscript (j) in the sequence of variables {r^(j)_i,k}_i∈[m_j], k∈[N_1] for this particular j. Note that N_2/8 ≤ m_j ≤ N_2 by the property that ρ(Δ_0, Δ_j)≥ N_2/8. Moreover, by the design of f, there are a total of m_jN_1 trapezoids in the union of the intervals [x^*_i - h_2, x^*_i + h_2] for those i such that Δ_j,i = 1. Define _j as the joint distribution of {(X_i,Y_i)}_i=1^n under _j but with {r_i,k}_i∈[m_j], k∈[N_1] replaced by a sequence of i.i.d. standard normal variables denoted as {r_i,k}_i∈[m_j], k∈[N_1]. By definition, it holds that 

    K(_j, _0; Ω_n)    = ∫_Ω_n p_jlogp_j/p_0
        = ∫_Ω_n p_jlogp_j/p_j + ∫_Ω_np_jlogp_j/p_0 + ∫_Ω_n (p_j - p_j) logp_j/p_0
       = K(_j, _j;Ω_n) + K(_j, _0; Ω_n) + ∫_Ω_n (p_j - p_j) logp_j/p_0

for density functions with respect to some common dominating measure. Next, we will show respectively that, by matching the moments of {r_i,k}_i∈[m_j],k∈[N_1] and the standard Gaussian random variable up to some sufficiently high order, it holds that

    K(_j, _j; Ω_n)≲ 1,   K(_j, _0;Ω_n)≲log(M),   and   ∫_Ω_n (p_j - p_j) log*p_j/p_0≲ 1.


First note that, by denoting x (x_1,…,x_n), dx dx_1… dx_n and similarly for y and dy, we have

    K(_j, _j; Ω_n)    = ∫1*Ω_np(x)dx∫log*d_j(y|x)/d_j(y|x)_j(dy|x)
       = *1*Ω_nK*_j(y|x), _j(y|x)
       ≤*1*Ω_nχ^2*_j(y|x), _j(y|x)
       χ^2*_j, _j; Ω_n,

where the inequality follows by Lemma 2.7 in <cit.>. Therefore the first inequality K(_j, _j; Ω_n)≲ 1 holds by Lemma <ref>. 

Next we prove K(_j, _0;Ω_n)≲log(M). Again, it suffices to prove that for any realization of {X_i}_i=1^n in Ω_n, it holds that K(_j(y|x), _0(y|x)) ≲log(M) ≍ N_2. Note that under _0, y|x∼N_n(0, Σ_0), with Σ_0 = 𝐈_n. Recall that the location index sequence {b_i}_i=1^n = {(k_i,ℓ_i)}_i=1^n takes the value (k_i,ℓ_i) = (k,ℓ) if X_i∈[x^*_k - h_2 + (4ℓ-3)h_1, x^*_k - h_2 + (4ℓ-1) h_1]. Then, due to the symmetry of {r_i,k}_i∈[m_j],k∈[N_1], the design of the nonparametric component f, and the fact that K(·) takes value 1 on [-1, 1], it holds that under _j, y|x∼N_n(0, Σ_1), with 

    (Σ_1)_ii = 1 - θ_n^21*Δ_j,k_i = 1 +  h_1^2α1*Δ_j,k_i = 1 = 1

and (Σ_1)_i_1i_2 = h_1^2α1*Δ_j,k_i_1 = 1, (k_i_1,ℓ_i_1) = (k_i_2,ℓ_i_2) for i_1≠ i_2. Define 

    N_0∑_i_1≠ i_21*Δ_j,k_i_1 = 1, (k_i_1,ℓ_i_1) = (k_i_2,ℓ_i_2).

Then, by the proof of Lemma 3.6 in <cit.>, it holds that

    K(_j(y|x), _0(y|x)) ≤ Ch_1^4αN_0 = Cθ_n^4N_0.

Note that N_0 is a random variable that depends on {X_i}_i=1^n, and by direct calculation we have 

    (N_0) ≤ n^2m_jh_2h_1 ≍ n^2N_2h_1h_2.

Putting together the pieces, we have K(_j, _0; Ω_n) ≤θ_n^4n^2h_1h_2N_2 ≲ N_2. This completes the proof of the second inequality.

Lastly, we show that ∫_Ω_n (p_j - p_j) log*p_j/p_0≲ 1. First note that

    ∫_Ω_n (p_j - p_j) log*p_j/p_0
       ≤∫_Ω_n*p_j - p_j*log(p_j/p_0)
       ≤*∫_Ω_n*p_j - p_j^1/2*∫*p_j - p_jlog^2(p_j/p_0)^1/2
       ≤*∫_Ω_n*p_j - p_j^1/2**∫ p_j log^2(p_j/p_0)^1/2 + *∫p_j log^2(p_j/p_0)^1/2.

By Lemmas <ref> and <ref>, by matching moments up to some sufficiently high order, the first term above can be upper bounded (up to some constant) by n^-η for any η > 0, therefore it suffices to show that both ∫ p_j log^2(p_j/p_0) and ∫p_j log^2(p_j/p_0) can be upper bounded by some polynomial of n of fixed order. Consider any realization of {X_i}_i=1^n in Ω_n, and assume that based on their location indices {b_i}_i=1^n, the n data points are partitioned into L_1 + L_2 clusters with cardinality s_ℓ such that the X_i's in the same cluster have the same value b_i. Moreover, for each data point in the first L_1 clusters, the location index b_i = (k_i,ℓ_i) satisfies that Δ_j,k_i = 1 while for the data points in the last L_2 clusters, it holds that Δ_j,k_i = 0. Apparently, we have the relations 1≤ L_1 + L_2 ≤ n, ∑_ℓ=1^L_1+L_2s_ℓ = n and 1≤ s_ℓ≤ K for ℓ∈[L_1 + L_2]. Moreover, denoting _j,π_ℓ and _0,π_ℓ (resp. p_j,π_ℓ and p_0,π_ℓ) for each ℓ∈[L_1 + L_2] as the joint distribution (resp. density) of those Y_i's in the ℓth cluster conditioning on the given realization {X_i}_i=1^n under _j and _0, we have 

    p_j = ∏_ℓ=1^L_1 + L_2p_j,π_ℓ   and     p_0 = ∏_ℓ=1^L_1 + L_2 p_0,π_ℓ.

Moreover, for any L_1 + 1≤ℓ≤ L_1 + L_2, it holds that p_j,π_ℓ = p_0,π_ℓ, therefore it holds that

    log^2*p_j/p_0 = *∑_ℓ=1^L_1log(p_j,π_ℓ) - log*p_0,π_ℓ^2 ≲ n∑_ℓ=1^L_1log^2(p_j,π_ℓ/p_0,π_ℓ).

Now consider any ℓ∈[L_1] and assume that s_ℓ = d for some positive integer d. Without loss of generality, assume the y_i's in this cluster are {y_1,…,y_d}, and they take the form Y_i = h_1^αr_1,1 + (1-h_2^β)^1/2ε_i = θ_nr_1,1 + (1-h_2^β)^1/2ε_i under _j and Y_i = ε_i under _0. Define σ^2  (1 - h_2^β) which is positive for large enough n. Then, the previous equalities imply that

    p_0,π_ℓ = φ(y_1)…φ(y_d) = (2π)^-d/2exp*-∑_i=1^d y_i^2/2

and

    p_j,π_ℓ   = ∫1/σφ*y_1-θ_nv/σ…1/σφ*y_d-θ_nv/σφ(v)dv
       = (2π)^-d/21/σ^d-1(dθ_n^2 + σ^2)^1/2exp*-∑_i=1^d y_i^2/2σ^2 + (∑_i=1^d y_iθ_n)^2/2σ^2(dθ_n^2 + σ^2)^2.
 
Putting together the pieces, we obtain that

    log^2(p_j,π_ℓ/p_0,π_ℓ) ≲ d^2log^2(1/σ) + *∑_i=1^d y_i^2^2 + *∑_i=1^d y_iθ_n^4 ≲ 1 + ∑_i=1^d y_i^4.

Therefore we have

    ∫ p_j log^2(p_j/p_0) ≲ n∑_ℓ=1^L_1∫ p_j(1 + ∑_i=1^d y_i^4) ≲ n∑_ℓ=1^L_1∑_i=1^d ∫ y_i^4_j(dy_i) ≲ n^2,

where we use the fact that L_1 ≤ n. Similarly, we have ∫p_j log^2(p_j/p_0)≲ n^2. The proof is thus complete.




 §.§ Supporting lemmas












































































































































Suppose f∈Λ_α(C_F), V∈Λ_β(C_V), σ^2≤ C_σ for some fixed constants C_F,C_V,C_σ, and the joint distribution of (X,ε) belongs to P_vf,(X,ε). Then, with Ω_n defined in the proof of Theorem <ref>, it holds that

    ***V(x^*) - V(x^*)≤ C*h_1^2(α∧ 1) + h_2^β + h_1^β∧ 1 + τ_n

for some fixed positive constant C that only depends on α,β,C_F,C_V,C_σ,C_0,C_ε.



We adopt the notation ℓ,q(·), 𝐁_n, X_ij, and K_ij from the proof of Theorem <ref>. Also recall the definition of 𝐁_n^*, D_ij, w_ij, and w_ij from the definition of V(x^*). Writing _ε as the conditional expectation given {X_i}_i=1^n, it holds that

    ***V(x^*) - V(x^*) = **_ε*V(x^*) - V(x^*).





Then V(x^*) = ∑_i<jw_ijD_ij and

    **_ε*V(x^*) - V(x^*)
       ≤**∑_i<jw_ij*_ε D_ij - V(x^*) + V(x^*)τ_n** (|𝐁_n| + τ_n)^-1.

By definition, it holds on Ω_n that 𝐁_n - 𝐁≤λ_min(𝐁)/2, where λ_min(𝐁) is the smallest eigenvalue of 𝐁. Thus by Weyl's inequality, it holds that λ_min(𝐁_n)≥λ_min(𝐁)/2 ≥ c for some fixed positive constant c due to the invertibility of 𝐁 as proved in Theorem <ref> under Condition (c) in P_vf,(X,ε). Then, using the fact that |𝐁_n| ≥λ^ℓ+1_min(𝐁_n) and the boundedness of V(·), it holds that

    **_ε*V(x^*) - V(x^*)≤**∑_i<jw_ij*_ε D_ij - V(x^*) + Cτ_n.

Direct calculation shows that _ε D_ij - V(x^*)= (f(X_i) - f(X_j))^2/2 + (V(X_i) - V(x^*))/2 + (V(X_j) - V(x^*))/2. Due to symmetry, we only need to control the first two terms. For the first term, using the fact 𝐁_n^* = |𝐁_n|𝐁_n^-1 on Ω_n (𝐁_n invertible on Ω_n), we have

    *∑_i<j (|𝐁_n| +τ_n)^-1w_ij(f(X_i) - f(X_j))^2
       = *n 2^-1∑_i<j (|𝐁_n|+τ_n)^-1q^⊤(0)𝐁_n^*q*X_ij-x^*/h_2K_ij(f(X_i)-f(X_j))^2
       = *(|𝐁_n|+τ_n)^-1q^⊤(0)𝐁_n^*q*X_ij - x^*/h_2K_ij(f(X_i)-f(X_j))^2
       ≲ h_1^2(α∧ 1)*(|𝐁_n|+τ_n)^-1𝐁_n^*K_ij = h_1^2(α∧1)*|𝐁_n|/|𝐁_n| + τ_n𝐁_n^-1K_ij
       ≲ Ch_1^2(α∧ 1) K_ij≲ Ch_1^2(α∧ 1).

Here, the second inequality follows from the fact that 𝐁_n^-1 = λ_min(𝐁_n)^-1 is bounded by some fixed constant on Ω_n, and in the last inequality we use the fact that K_ij is bounded Condition (b) in P_vf,(X,ε). For the second term, it holds that

    *∑_i<jw_ij(V(X_i) - V(x^*))
       ≤*∑_i<jw_ij(V(X_ij) - V(x^*)) + *∑_i<jw_ij(V(X_i) - V(X_ij)) I + II.

For I, using the Hölder property of V(·) and the reproducing property of local polynomial estimators (see (<ref>) in the main paper), that is,

    ∑_i<j w_ij(X_ij - x^*)^k = ∑_i<jw_ij(X_ij - x^*)^k = 0,    k=1,2,…, ℓ,

it holds that

    **∑_i<jw_ij(V(X_ij) - V(x^*))
       = **∑_i<jw_ij*∑_k=1^ℓ-1V^(k)(x^*)/k!(X_ij - x^*)^k + V^(ℓ)(x^* + τ(X_ij - x^*))/ℓ!(X_ij - x^*)^ℓ
       = **∑_i<jw_ij(X_ij-x^*)^ℓ/ℓ!*V^(ℓ)(x^* + τ(X_ij - x^*)) - V^(ℓ)(x^*)
       = ** (|𝐁_n|+τ_n)^-1q^⊤(0)𝐁_n^*q*X_ij-x^*/h_2K_ij(X_ij-x^*)^ℓ/ℓ!*V^(ℓ)(x^* + τ(X_ij - x^*)) - V^(ℓ)(x^*)
       ≲* |X_ij -x^*|^β (|𝐁_n|+τ_n)^-1𝐁_n^*K_ij
       ≲ h_2^β*𝐁_n^-1K_ij≲ h_2^β,

where in the second line we use the Taylor expansion of V(·) around x^* with some τ∈[0,1], in the fifth line we use the fact that q((X_ij-x^*)/h_2) has bounded ℓ_2 norm due to the compact support of K(·), and in the last inequality we use again the fact that 𝐁_n^-1 is bounded by some fixed constant on Ω_n. With a similar calculation and the fact that β-smooth functions are Lipschitz for β≥ 1, we have II ≲ h_1^β∧ 1. Therefore, putting together the pieces, we obtain

    ***V(x^*) - V(x^*)≤ C*h_1^2(α∧ 1) + h_2^β + h_1^β∧ 1 + τ_n.

This completes the proof.




Suppose f∈Λ_α(C_F), V∈Λ_β(C_V), σ^2≤ C_σ for some fixed constants C_F,C_V,C_σ, and the joint distribution of (X,ε) belongs to P_vf,(X,ε). Assume that nh_2→∞ and n^2h_1h_2→∞ as n→∞. Then, with Ω_n defined in the proof of Theorem <ref>, it holds that

    *V(x^*)≤ C*n^-1h_2^-1 + n^-2(h_1h_2)^-1 + τ_n^2

for some fixed positive constant C that only depends on α,β,C_F,C_V,C_σ,C_0,C_ε.


We adopt the notation ℓ,q(·), 𝐁_n, X_ij, and K_ij from the proof of Theorem <ref>. Also recall the definition of 𝐁_n^* and D_ij from the definition of V(x^*). Define the vector-valued U-statistic 

    U_nn 2^-1g(X_i,X_j) n 2^-1∑_i<j D_ijq*X_ij-x^*/h_2K_ij

and let θg(X_i,X_j)∈^ℓ+1. Then, for each j∈[ℓ+1], writing _ε as the conditional expectation given {X_i}_i=1^n, we have

    θ_j    = *D_ij((X_ij-x^*)/h_2)^j/j!K_ij = *((X_ij-x^*)/h_2)^j/j!K_ij_ε D_ij
       = *((X_ij-x^*)/h_2)^j/2j!K_ij((f(X_i) -f(X_j))^2 + V(X_i) + V(X_j))≲ K_ij≲ 1,

where we have used Condition (b) in P_vf,(X,ε), boundedness of V(·) and compact support of K(·). Therefore we have θ = O(1). With the above notation, and using the fact that 𝐁_n is invertible and satisfies 𝐁_n^* = |𝐁_n|𝐁_n^-1 on Ω_n, we have

    V(x^*) = |𝐁_n|/|𝐁_n|+τ_nq^⊤(0)𝐁_n^-1U_n = q^⊤(0)𝐁_n^-1U_n - τ_n/|𝐁_n|+τ_nq^⊤(0)𝐁_n^-1U_n.

Thus in order to upper bound *V(x^*), it suffices to upper bound the variances of the two terms in the above display. For the second term, using the fact that |𝐁_n| is bounded away from zero on Ω_n under Condition (c) in P_vf,(X,ε), we have

    *τ_n/|𝐁_n|+τ_nq^⊤(0)𝐁_n^-1U_n   ≲τ_n^2*q^⊤(0)𝐁_n^-1U_n^2 ≲τ_n^2*𝐁_n^-1^2U_n^2≲τ_n^2,

where the last inequality follows since θ = O(1) and U_n concentrates to θ by Lemma <ref> since nh_2→∞ and n^2h_1h_2→∞. The first term can be decomposed as

    q^⊤(0)𝐁_n^-1U_n = q^⊤(0)(𝐁_n^-1 - 𝐁^-1)U_n + q^⊤(0)𝐁^-1U_n I + II.

For the first term, it holds on the event Ω_n that

    𝐁_n^-1 - 𝐁^-1 = 𝐁^-1𝐁_n^-1𝐁_n - 𝐁≤𝐁^-1(𝐁_n^-1 - 𝐁^-1 + 𝐁^-1)𝐁_n - 𝐁.

Thus on the event Ω_n, it holds that 𝐁_n^-1 - 𝐁^-1≤𝐁^-1^2𝐁_n - 𝐁/(1 - 𝐁^-1𝐁_n - 𝐁)≤ 2𝐁^-1^2𝐁_n-𝐁. This implies that

    *I   ≤*q^⊤(0)(𝐁_n^-1 - 𝐁^-1)U_n^2 ≤*(𝐁_n^-1 - 𝐁^-1)^2U_n^2
       ≲*𝐁_n -𝐁^2U_n^2≤ (𝐁_n - 𝐁^4)^1/2(U_n^4)^1/2.

Clearly, (U_n^4)^1/2 = O(1) since U_n concentrates to θ and θ = O(1), and by Lemmas <ref> and <ref>, it holds that (𝐁_n - 𝐁^4)^1/2≲ n^-1h_2^-1 + n^-2h_1^-1h_2^-1. This concludes that

    (I) ≲ n^-1h_2^-1 + n^-2h_1^-1h_2^-1.

Lastly, for II, writing Zq^⊤(0)𝐁^-1n 2^-1∑_i<jg(X_i,X_j), we have

    (II) = (Z - (Z))^2 = ((Z- Z) + (Z) - Z)^2 ≲(Z)  + (Z)^2.
 
By Lemma <ref>, it holds that

    (Z) ≲ n^-1h_2^-1 + n^-2h_1^-1h_2^-1.

Lastly, by Cauchy's inequality and Lemma <ref>, 

    (Z)^2 ≤ ( Z^4)^1/2^1/2(Ω_n^c) ≲*𝐁^-1^4U_n^4^1/2^1/2(Ω_n^c) = o(n^-1h_2^-1 + n^-2h_1^-1h_2^-1).

Thus we conclude that

    (II) ≲ n^-1h_2^-1 + n^-2h_1^-1h_2^-1.

Putting together the pieces, we have proved 

    *V(x^*) - *V(x^*)^2 ≲ n^-1h_2^-1 + n^-2(h_1h_2)^-1 + τ_n^2.

This completes the proof.




Consider the term Z defined in the proof of Lemma <ref>:

    Z = q^⊤(0)𝐁^-1n 2^-1∑_i<j D_ijq*X_ij-x^*/h_2K_ij.

Then, under the same conditions of Lemma <ref>, it holds that 

    *Z≤ C*n^-1h_2^-1 + n^-2h_1^-1h_2^-1

for some fixed positive constant C that only depends on α,β,C_0.


Denote g(X_i, X_j) D_ijq*X_ij-x^*/h_2K_ij and θg(X_i,X_j). Then, it holds that

    (Z)    = *q^⊤(0)𝐁^-1n 2^-1∑_i<j (g(X_i, X_j) - θ)^2
       ≤𝐁^-1^2*n 2^-1∑_i<j(g(X_i,X_j)-θ)^2
       ≲ n^-4∑_i<j,i^' < j^'*(g(X_i,X_j) - θ)^⊤(g(X_i^',X_j^') - θ),

where the last inequality follows since 𝐁^-1 is invertible under Condition (c) in P_vf,(X,ε) as proved in Theorem <ref>. Apparently, when i,j,i^', j^' are all different, the summand equals to zero. When i,j,i^', j^' take three different values, say, i = i^' < j < j^', we have

    *(g(X_i,X_j) - θ)^⊤(g(X_i,X_j^') - θ)
       = *D_ijD_ij^'q^⊤*X_ij-x^*/h_2q^⊤*X_ij^'-x^*/h_2K_ijK_ij^' - θ^2.

Let Z_1g(X_i, X_j) and Z_2 g(X_i, X_j^'). Then, for any k∈[ℓ+1], it holds that *Z_1,k≲ D_ijK_h_1(X_i-X_j)K_h_2(X_ij-x^*) and similarly for Z_2,k. Therefore, using the finite fourth moment of ε in Condition (d) of P_vf,(X,ε) in the calculation of _ε (D_ijD_ij^') and the fact that both f(·) and V(·) are bounded, we have

    (|Z_1,kZ_2,k|)
       ≲*D_ijK_h_1(X_i-X_j)K_h_2(X_ij-x^*)D_ij^'K_h_1(X_i - X_j^')K_h_2(X_ij^'-x^*)
       ≲*K_h_1(X_i-X_j)K_h_2(X_ij-x^*)K_h_1(X_i - X_j^')K_h_2(X_ij^'-x^*)
       = ∫_^31/h_1^2h_2^2K*v-u/h_1K*u+v/2 - x^*/h_2K*w-u/h_1K*u+w/2 - x^*/h_2p_X(u)p_X(v)p_X(w)dudvdw
       = 1/h_2∫_^3 K(u)K(v)K(w)K*uh_2 + h_1(w-v)/2/h_2p_X(s_1)p_X(s_2)p_X(s_3)dudvdw
       ≲ 1/h_2,

where we again invoke Condition (b) in P_vf,(X,ε) and the compact support of K(·), and

    s_1 = uh_2+x^*-h_1v/2,    s_2 = uh_2+x^*+h_1v/2,    s_3 = uh_2 + x^* + wh_1 - vh_1/2.

This, along with the fact that θ = O(1), concludes that 

    *(g(X_i,X_j) - θ)^⊤(g(X_i^',X_j^') - θ)≲ h_2^-1

when i,j,i^', j^' take three different values. Similarly, one can prove that when i,j,i^', j^' take two different values, that is i = i^' and j = j^',

    *(g(X_i,X_j) - θ)^⊤(g(X_i^',X_j^') - θ)≲ (h_1h_2)^-1.

Putting together the pieces, we obtain that 

    (Z) ≲n^3h_2^-1 + n^2(h_1h_2)^-1/n^4 = n^-1h_2^-1 + n^-2h_1^-1h_2^-1.

This completes the proof.




Suppose Condition (b) in P_vf,(X,ε) holds. Assume n^2h_1h_2→∞, nh_2→∞, and h_1/h_2→ 0. Then, for any u,v > 0, the matrices 𝐁_n and 𝐁 defined in the proof of Theorem <ref> satisfy

    *𝐁_n - 𝐁≥ C(v^1/2n^-1/2h_2^-1/2 + u^1/2n^-1h_1^-1/2h_2^-1/2)≤ C*exp(-u) + exp(-v)

for some fixed positive constant C.


Using a standard entropy argument (see, for example, Lemma 5.3 in <cit.>), it holds that for any t>0,

    (𝐁_n - 𝐁≥ t) ≤ Nmax_1≤ i≤ N**a_i^⊤(𝐁_n - 𝐁)a_i≥ t/2,
 
where N  5^ℓ+1 and {a_i}_i=1^N is a 1/2-net on the unit sphere 𝕊^ℓ. We now upper bound **a^⊤(𝐁_n - 𝐁)a≥ t/2 for an arbitrary a∈𝕊^ℓ with the help of Lemma <ref>. For this, we upper bound the five quantities B_1,B_2,B_3,ν_1^2,ν_2^2 therein. Denote the kernel of a^⊤𝐁_na as g and its linear part as g_1, that is,

    g(X_i, X_j) *a^⊤q*X_ij-x^*/h_2^2 K_ij and  g_1(X_i) = *g(X_i,X_j)| X_i.

Then, we have

    g_1(X_i)    = ∫*a^⊤q*(u+X_i)/2 - x^*/h_2^21/h_1K*u - X_i/h_11/h_2K*(u+X_i)/2 - x^*/h_2p_X(u)du
       ≲1/h_2∫ K(u)K*X_i + h_1u/2 - x^*/h_2p_X(X_i + uh_1)du≲ h_2^-1,

where in the first inequality we use the fact that q*(u+X_i)/2 - x^*/h_2 has bounded ℓ_2 norm due to the compact support of K(·), and in the last inequality we apply Condition (b) in P_vf,(X,ε). Therefore B_1≲ h_2^-1 and one can similarly show that ν_1^2≲ h_2^-1.

For B_2, we have

    B_2^2    = nsup_X_i**a^⊤q*X_ij-x^*/h_2^4 1/h_1^2K^2*X_i-X_j/h_11/h_2^2K^2*X_ij-x^*/h_2| X_i
       ≲n/h_1h_2^2sup_X_i*1/h_1K*X_i - X_j/h_1| X_i
       = n/h_1h_2^2sup_X_i∫ K(u)p_X(X_i + uh_1)du ≲ nh_1^-1h_2^-2.

This concludes that B_2≲ n^1/2h_1^-1/2h_2^-1. Moreover, one can easily show that ν_2^2 ≲ (h_1h_2)^-1 and B_3≲ (h_1h_2)^-1. Putting together the pieces and applying Lemma <ref>, we obtain that for any u,v>0,

    **a^⊤(𝐁_n - 𝐁)a≥ a_1v^1/2 + a_2v + b_1u^1/2 + b_2u + b_3u^3/2 + b_4u^2≤ C*exp(-v) + exp(-u),

where a_1≲ n^-1/2h_2^-1/2, a_2≲ n^-1h_2^-1 and b_1≲ n^-1h_1^-1/2h_2^-1/2, b_2≲ n^-1h_2^-1, b_3≲ n^-3/2h_1^-1/2h_2^-1, b_4≲ n^-2h_1^-1h_2^-1. Under the conditions n^2h_1h_2→∞, nh_2→∞ and h_1/h_2→ 0 as n→∞, the dominant terms are a_1 and b_1, that is,

    n^-1/2h_2^-1/2∨ n^-1h_1^-1/2h_2^-1/2.

This completes the proof.












































Under the setting and conditions of Theorem <ref>, for any positive η > 0, there exists an i.i.d. sequence {r_i,k}_i∈[m_j],k∈[N_1] (with m_j, N_1 defined in Theorem <ref>) with range contained in [-B, B] for some B only depending on α,β,η, such that the probability measures _j and _j defined therein satisfy that

    χ^2*_j, _j; Ω_n≲ n^-η,

where, for any measurable subset E and two probability measures  and , χ^2(,; E) is the conditional χ^2-distance defined as

    χ^2(,; E) ∫_E(p-q)^2/q

with p,q being the densities of  and  with respect to some common dominating measure.


Without loss of generality, let j = 1. First note that the conditional χ^2-distance can be written as

    χ^2*_1, _1; Ω_n   = ∫1{Ω_n}*p_1(x, y) - p_1(x, y)^2/p_1(x, y) dxdy
       = ∫1{Ω_n}p(x)dx∫*p_1(y|x) - p_1(y|x)^2/p_1(y|x)dy
       = *1{Ω_n}χ^2*_1(y|x), _1(y|x),

where we have used p(·) to represent the density of {X_i}_i=1^n under both _1 and _1.

Recall the definition of the location index sequence {b_i}_i=1^n in Theorem <ref>. Consider any realization of {X_i}_i=1^n in Ω_n, and assume that based on their location indices {b_i}_i=1^n, {X_i}_i=1^n is partitioned into L_1 + L_2 clusters with cardinality s_ℓ such that the X_i's in the same cluster have the same value b_i. Moreover, for those data points in the first L_1 clusters, the location index b_i = (k_i,ℓ_i) satisfies that Δ_1,k_i = 1 while for the data points in the last L_2 clusters, it holds that Δ_1,k_i = 0 (recall the definition of Δ_1 = (Δ_1,1,…, Δ_1,N_2) in the lower bound design of Theorem <ref>. Apparently, we have the relations 1≤ L_1 + L_2 ≤ n, ∑_ℓ=1^L_1+L_2s_ℓ = n and 1≤ s_ℓ≤ K for ℓ∈[L_1 + L_2] (recall the definition of K in Theorem <ref>). Moreover, denoting _1,π_ℓ and _1,π_ℓ (resp. p_1,π_ℓ and p_1,π_ℓ) for each ℓ∈[L_1 + L_2] as the joint distribution (resp. density) of those Y_i's in the ℓth cluster conditioning on the given realization {X_i}_i=1^n under _1 and _1, we have 

    χ^2*_1(y|x), _1(y|x)   = ∏_ℓ=1^L_1+L_2*1 + χ^2*_1,π_ℓ, _1,π_ℓ - 1
       =  ∏_ℓ=1^L_1*1 + χ^2*_1,π_ℓ, _1,π_ℓ - 1
       ≤exp*∑_ℓ=1^L_1χ^2*_1,π_ℓ, _1,π_ℓ - 1
  
where the first equality follows by mutual independence of data points in each cluster, the second inequality follows from the fact that for those data points Y_i's in the latter L_2 clusters, each Y_i| X_i∼N(0,1) under both _1 and _1. Since L_1≤ n, it suffices to show that for any realization of {X_i}_i=1^n in Ω_n, by matching enough moments, we have χ^2*_1,π_ℓ, _1,π_ℓ≤ n^-η for any η > 0.

For each ℓ∈[L_1], *p_1,π_ℓ - p_1,π_ℓ only depends on the ℓth cluster via its cardinality, which we now control for a general cluster size 1≤ d≤ K. Without loss of generality, we assume that ℓ=1 and the y_i's in this cluster are {y_1,…, y_d} with common location index b_i = (1,1). In this case, under the choice of θ_n^2 and h_1 given in Theorem <ref>, we clearly have Y_i = θ_nr_1,1 + (1-h_2^β)^1/2ε_i under _1 and Y_i = θ_nr_1,1 + (1-h_2^β)^1/2ε_i under _1 for i∈[d], where the sequence {ε_i}_i=1^d follows the standard normal distribution under both _1 and _1. Define σ^2  1-h_2^β = 1-θ_n^2. Then, it holds that

    p_1,π_ℓ(y_1,…,y_d)    = ∫_-∞^∞1/σφ*y_1-θ_nv/σ…1/σφ*y_d - θ_nv/σ(dv),
    p_1,π_ℓ(y_1,…,y_d)    = ∫_-∞^∞1/σφ*y_1-θ_nv/σ…1/σφ*y_d - θ_nv/σφ(v)dv,

where  is the distribution of {r_i,k}_i∈[m_j],k∈[N_1]. Using the well-known equality φ(t - θ_nv) = φ(t)*∑_k=0^∞ v^kθ_n^kH_k(t)/k! for any t,v, where H_k is the kth order Hermite polynomial, it holds that

    φ*y_1 -θ_nv/σ…φ*y_d -θ_nv/σ
       = φ*y_1/σ…φ*y_d/σ∑_k_1,…,k_d = 0^∞*θ_nv/σ^k_1+…+k_dH_k_1(y_1/σ)/k_1!…H_k_d(y_d/σ)/k_d!
       = φ*y_1/σ…φ*y_d/σ∑_k=0^∞*θ_nv/σ^k∑_k_1+…+k_d = k^∞*θ_nv/σ^kH_k_1(y_1/σ)/k_1!…H_k_d(y_d/σ)/k_d!,

and therefore by the symmetry of  and by matching the moments of  and the standard normal distribution up to order 2p for some positive integer p to be chosen later, we obtain

    p_1,π_ℓ(y_1,…,y_d) - p_1,π_ℓ(y_1,…,y_d)
       = φ*y_1/σ…φ*y_d/σ∑_k=0^∞*θ_n/σ^k ∑_k_1+…+k_d=kH_k_1(y_1/σ)/k_1!…H_k_d(y_d/σ)/k_d!∫ v^k( - Φ)(dv)
       = φ*y_1/σ…φ*y_d/σ∑_k=p^∞*θ_n/σ^2k∑_k_1+…+k_d=2kH_k_1(y_1/σ)/k_1!…H_k_d(y_d/σ)/k_d!∫ v^2k( - Φ)(dv).

Define δ_2k∫ v^2k( - Φ)(dv). Then, the above inequality further implies that

    *p_1,π_ℓ(y_1,…,y_d) - p_1,π_ℓ(y_1,…,y_d)^2
        = φ^2*y_1/σ…φ^2*y_d/σ∑_k,ℓ=p^∞*θ_n/σ^2k+2ℓ∑_k_1+…+k_d=2k
    ℓ_1+…+ℓ_d = 2ℓH_k_1(y_1/σ)/k_1!H_ℓ_1(y_1/σ)/ℓ_1!…H_k_d(y_d/σ)/k_d!H_ℓ_d(y_d/σ)/ℓ_d!

On the other hand, letting Z∼, we have

    p_1,π_ℓ(y_1,…,y_d)    = ∫1/σφ*y_1-θ_nv/σ…1/σφ*y_d-θ_nv/σ(dv)
       = 1/σ^dφ*y_1/σ…φ*y_d/σ∫exp*-d/2σ^2(θ_nv)^2 +∑_i=1^d y_iθ_nv/σ^2(dv)
       = 1/σ^dφ*y_1/σ…φ*y_d/σ*exp*-d/2σ^2θ_n^2Z^2+ ∑_i=1^d y_iθ_n/σ^2Z
       ≥1/σ^dφ*y_1/σ…φ*y_d/σexp*-dθ_n^2/2σ^2,

where the last inequality follows from Jensen's inequality and the fact Z = 0, Z^2 = 1 from moment matching. Combining (<ref>) and (<ref>), we obtain that

    χ^2(_1,π_ℓ, _1,π_ℓ)    = ∫*p_1,π_ℓ(y_1,…,y_d) - p_1,π_ℓ(y_1,…,y_d)^2/p_1,π_ℓ(y_1,…,y_d)dy_1… dy_d
       ≤σ^d∑_k,ℓ=p^∞δ_2kδ_2ℓ∑_k_1+…+k_d=2k
    ℓ_1+…+ℓ_d = 2ℓ*θ_n/σ^2k+2ℓ∏_j=1^d∫φ*y_j/σH_k_j(y_j/σ)/k_j!H_ℓ_j(y_j/σ)/ℓ_j!dy_j
       = σ^2d∑_k,ℓ=p^∞δ_2kδ_2ℓ∑_k_1+…+k_d=2k
    ℓ_1+…+ℓ_d = 2ℓ*θ_n/σ^2k+2ℓ∏_j=1^d∫φ*y_jH_k_j(y_j)/k_j!H_ℓ_j(y_j)/ℓ_j!dy_j
       ≤∑_k=p^∞*θ_n/σ^4kδ_2k^2∑_k_1+…+k_d=2k1/k_1!…1/k_d!,
 
where the last inequality follows from the fact that σ≤ 1 and ∫φ(t)H_k(t)H_ℓ(t)dt = k!1*k = ℓ. Now, using the multinomial identity 

    ∑_k_1+…+k_d = 2k(2k)!/k_1!… k_d!*1/d^2k = 1,

we obtain that ∑_k_1+…+k_d=2k1/(k_1!… k_d!) = d^2k/(2k)!, therefore it holds that

    χ^2(_1,π_ℓ, _1,π_ℓ) ≤∑_k=p^∞*θ_n√(d)/σ^4k1/(2k)!δ_2k^2.

Now, by Lemma <ref>, for any positive integer p,  we can find a symmetric distribution  that has the same first p moments as the standard normal distribution and is compactly supported on [-B,B] for some B that only depends on p. This combined with the fact that ∫ t^2kφ(t)dt = (2k-1)!! implies that δ_2k^2≲ B^4k + (2k)!. We therefore obtain

    χ^2(_1,π_ℓ, _1,π_ℓ) ≲∑_k=p^∞*θ_n√(d)/σ^4k≲ n^-η

by choosing a sufficiently large p that only depends on α,β and η, where we also use the fact that d is bounded by an absolute constant and for sufficiently large n, it holds that σ > 1/2. This completes the proof.




For some M≥ 2, let _0,_1,…, _M be M+1 hypotheses on some measurable space (X, A) such that for each 0 ≤ i≠ j≤ M, _i and _j are mutually absolutely continuous. Let Ω be a measurable subset of X such that _j(Ω) is identical for all 0≤ j≤ M. Define the “conditional" version of Kullback divergence as

    K(,; Ω) ∫_Ωlog*d/dd.

Then, if

    1/M∑_j=1^MK(_j, _0; Ω) ≤ c^*log(M)

for some 0 < c^* < 1/8, the following statement holds

    inf_ψp_e,M(ψ) ≥√(M)/1+√(M)*_0(Ω) - 2c^* - √(2c^*/log(M)),
  
where the infimum ranges over all tests taking values in {0,1,…, M} and p_e,M(ψ)max_0≤ j≤ M_j(ψ≠ j).


This is exactly the conditional version of Theorem 2.5 in <cit.>. We first show that subject to the condition

    1/M∑_j=1^MK(_j, _0; Ω) ≤ c

for some c > 0, for all 0 < τ < 1, it holds that

    1/M∑_j=1^M_j*d_0/d_j≥τ≥_0(Ω)- c^',

where c^' -(c+√(c/2))/log(τ). For this, we have for each j∈[M]

    _j*d_0/d_j≥τ   = _j*d_j/d_0≤1/τ
       =1 - *_j**d_j/d_0≥1/τ⋂Ω + _j**d_j/d_0≥1/τ⋂Ω^c
       ≥_0(Ω) - _j**d_j/d_0≥1/τ⋂Ω
       = _0(Ω) - _j**log*d_j/d_0≥log*1/τ⋂Ω
       ≥_0(Ω) - (log(1/τ))^-1__j*log*d_j/d_0_+1*Ω,

where in the third line we use the fact that _j(Ω) = _0(Ω), and for a real number a, a_+ max{0,a}. Let p_0 and p_j be the densities of _0 and _j with respect to some common dominating measure. Then, by definition of the conditional Kullback divergence, Lemma <ref>, and Lemma <ref>, it holds that

    __j*log*d_j/d_0_+1*Ω   = ∫_Ω p_j*logp_j/p_0_+
       = K(_j, _0; Ω) + ∫_Ω p_j*logp_j/p_0_-
       ≤ K(_j, _0; Ω) + (_j, _0; Ω)
       ≤ K(_j, _0; Ω) + √(K(_j, _0; Ω)/2).

Now, by the condition ∑_j=1^MK(_j, _0; Ω)/M ≤ c and Cauchy's inequality, it holds that

    1/M∑_j=1^M √(K(_j, _0; Ω))≤*1/M∑_j=1^M K(_j, _0; Ω)^1/2≤√(c).

We therefore conclude that (<ref>) is true.



Next, by Proposition 2.2 in <cit.>, we obtain that

    inf_ψp_e,M(ψ) ≥sup_0<τ<1τ M/τ M+ 1*1/M∑_j=1^M _j*d_0/d_j≥τ≥sup_0<τ<1τ M/τ M+ 1*_0(Ω) + c + √(c/2)/logτ.

Lastly, by choosing c = c^*log M and τ = 1/√(M), we obtain

    inf_ψp_e,M(ψ)    ≥sup_0<τ<1τ M/τ M+ 1*1/M∑_j=1^M _j*d_0/d_j≥τ
       ≥√(M)/1+√(M)*_0(Ω)-2c^* - √(2c^*/log M).


This completes the proof.




Let  and  be two probability measures on a measurable space (X,A) such that ≪, and Ω be a measurable subset of X. Define the conditional version of the total variation distance as follows

    (, ; Ω) sup_A∈A*(A⋂Ω) - (A⋂Ω).

Then, it holds that

    ∫_Ω*log*d/d_- ≤(,; Ω),

where a_- max{0, -a}.


Let p and q be the densities of  and  with respect to some common dominating measure, and define A{q≥ p > 0}. Then, it holds that

    ∫_Ω*log*d/d_-d   = ∫_Ω⋂{p,q>0} p*logp/q_- = ∫_A⋂Ω plogq/p≤∫_A⋂Ω q - p ≤(,; Ω).

This completes the proof.




Let  and  be two probability measures on a measurable space (X,A) such that ≪, and Ω be a measurable subset of X such that (Ω) = (Ω). For the conditional version of the Kullback divergence (defined in (<ref>)) and total variation distance (defined in (<ref>)),



it holds that

    (,; Ω) ≤√(K(,;Ω)/2).



Firstly, using the condition (Ω) = (Ω), it can be readily verified that the conditional total variation distance can be equivalently written as

    (,; Ω) = 1/2∫_Ω|p - q|,

where p and q are the densities of  and  with respect to some common dominating measure (cf. Lemma 2.1 in <cit.>). Then, following the proof of the first Pinsker's inequality in Lemma 2.5 in <cit.>, it holds that

    (,;Ω)    = 1/2∫_Ω*p-q
       = 1/2∫_Ω⋂{q>0}*p/q-1q
       ≤1/2∫_Ω⋂{q>0} q√(*4/3 + 2p/3qψ(p/q))
       ≤1/2*∫_Ω*4q/3+2p/3^1/2*∫_Ω q*p/qlogp/q - p/q + 1^1/2
       = √((Ω)/2)K^1/2(,; Ω)
       ≤√(K(,; Ω)/2),

where ψ(x)  xlog x - x + 1. This completes the proof.





§ PROOFS OF RESULTS IN SECTION <REF>

We only provide the proofs for Propositions <ref>, <ref>, <ref>, <ref>, <ref>-<ref>. The proofs of Propositions <ref>, <ref>, <ref>, and <ref> are straightforward, and the proof of Proposition <ref> is similar to that of Proposition <ref>.



 §.§ Proof of Proposition <ref>























Given the proof of Theorem <ref> and its supporting lemmas, the proof here is relatively straightforward. We only provide here a sketched version for completeness. For simplicity, we only prove the case with d = 2, and we will show that the desired upper bound can be achieved with the bandwidths choices 

    h_1≍ n^-2α_2/(4α_1α_2 + α_1 + α_2) and  h_2≍ n^-2α_1/(4α_1α_2 + α_1 + α_2).

C and c still represent two generic fixed positive constants whose values may change at each occurrence. 

Define U_1,U_2, θ_1, θ_2 and the “good" event E the same way as in Theorem <ref>. Following its proof, we now lower bound θ_2 and upper bound the term *θ_1 - θ_2σ^2. For θ_2, we have
















    θ_2    = *K_h_1(X_i,1-X_j,1)K_h_2(X_i,2-X_j,2)
       = ∫_^21/h_1h_2K*u/h_1K*v/h_2p_(u,v)dudv
       = ∫_^2 K(u)K(v)p_(uh_1, vh_2)dudv
       = ∫_-1^1 ∫_-1^1 K(u)K(v)p_(uh_1, vh_2)dudv
       ≥∫_U_(h_1,h_2) K(u)K(v)p_(uh_1, vh_2)dudv
       ≥λ(U_(h_1,h_2))inf_u∈U_(h_1,h_2)p_(u_1h_1,u_2h_2)M_K
       ≥ c_0^2M_K.

Here, the fourth equality follows from the kernel condition that K(·) is supported in [-1,1], and the set U_(h_1,h_2) starting from the first inequality follows from Condition (b) in P_mcv,(X,ε) since for any fixed δ_0 > 0 chosen therein, δ_∞(h_1, h_2)_∞≤δ_0 for sufficiently large n. For *θ_1 - θ_2σ^2, using the condition α_i∈(0,1], i=1,2, it holds that













    *θ_1 - θ_2σ^2
       = *K_h_1(X_ij,1)K_h_2(X_ij,2)(f(X_i)-f(X_j))^2/2
       ≲*1/h_1h_2K*X_ij,1/h_1K*X_ij,2/h_2**X_ij,1^2α_1 + *X_ij,2^2α_2
       = ∫1/h_1h_2K*u/h_1K*v/h_1**u^2α_1 + *v^2α_2p_(u,v)dudv
       = h_1^2α_1∫ K(u)K(v)|u|^2α_1p_(uh_1,vh_2)dudv + 
        h_2^2α_2∫ K(u)K(v)|v|^2α_2p_(uh_1,vh_2)dudv
       ≤ C_0*h_1^2α_1∫ K(u)K(v)|u|^2α_1dudv + h_2^2α_2∫ K(u)K(v)|v|^2α_2dudv
       ≤ C(h_1^2α_1 + h_2^2α_2),





where we have applied Condition (a) in P_mcv,(X,ε) and the compact support of K(·). Therefore, Lemmas <ref> and <ref> and the above estimates imply that 

    **U_1 - U_2σ^2/U_2^21{E}≲ (h_1^4α_1 + h_2^4α_2 + n^-1 + n^-2(h_1h_2)^-1) ≍ (n^-8α_1α_2/4α_1α_2 + α_1 + α_2 + n^-1).


Moreover, using the same argument as in the proof of Theorem <ref>, it holds that

    **U_1 - U_2σ^2/U_2^21{E^c} = o(n^-8α_1α_2/4α_1α_2 + α_1 + α_2 + n^-1).

This completes the proof for d = 2.



























In the case of general dimension d with heterogeneous smoothness index α = (α_1, …, α_d)^⊤, the upper bound takes the form

    *σ̂_d^2 - σ^2^2≲ n^-1 + n^-2*∏_k=1^dh_k^-1 + ∑_k=1^d h_k^4α_k

and we choose h_k ≍ n^-2α/(α_k(4α+d)). This completes the proof.




 §.§ Proof of Proposition <ref>


Given the proof of Theorem <ref>, the proof here is relatively straightforward. We will thus only present the construction of the hardest sub-problem. For simplicity, we will only prove the case for d = 2. We also only consider the regime of (α_1,α_2) in which the lower bound is sub-parametric: 4α_1α_2  < α_1+ α_2. Throughout the proof, C represents some generic positive constant and does not depend on n, and c represents a generic sufficiently small positive constant which also does not depend on n. In particular, c is always taken to be smaller than 1. Both C and c might have different values for each occurrence.

Introduce the following constants:

    θ_n^2  h_1^2α_1 h_2^2α_2 cn^-4α_1α_2/4α_1α_2 + α_1 + α_2,    N_1  1/(6h_1),    N_2  1/(6h_2),

where we tune the constant c in θ_n^2 so that N_1 and N_2 are both positive integers. We now specify f(·), distribution of X, σ^2 and distribution of ε in the null and alternative hypotheses, H_0 and H_1, respectively.



  Choice of σ^2: Under H_0, let σ^2 = 1 + θ_n^2. Under H_1, let σ^2 = 1.

  Choice of ε: Under both H_0 and H_1, let ε∼N(0,1).

  Choice of X: Under both H_0 and H_1, let X be uniformly distributed on the union of the rectangles [(6i_1-5)h_1, (6i_1-1)h_1]× [(6i_2-5)h_2, (6i_2-1)h_2] for i_1∈[N_1] and i_2∈[N_2].

  Choice of f(·): Under H_0, let f≡ 0. Under H_1, let f be a smooth bump function that takes value θ_nr_i_1,i_2 on the rectangle [(6i_1-5)h_1, (6i_1-1)h_1]× [(6i_2-5)h_2, (6i_2-1)h_2], and then smoothly decays to 0 on the union of the segments {x_1 = 6(i_1-1)h_1, 0 ≤ x_2≤ 1} for i_1∈[N_1] and {0≤ x_1≤ 1, x_2 = 6(i_2-1)h_2} for i_2∈[N_2]. Here, the double indexed sequence {r_i_1,i_2}_i_1∈[N_1], i_2∈[N_2] are N_1× N_2 i.i.d. symmetric and compactly supported random variables with distribution  satisfying

    ∫_-∞^∞ x^j(dx) = ∫_-∞^∞ x^jφ(x)dx,    j=1,…, q,

where q is some odd integer strictly larger than 1+(α_1+α_2)/(2α_1α_2).

In the definition of f(·) under H_1, the existence of the distribution G is guaranteed by Lemma <ref> and its range only depends on α_1 and α_2. The smoothness property of f(·) can be achieved by mollifying an indicator function.














We only verify Condition (c) in P_mcv,(X,ε), which holds by the convolution formula that for any 0≤ u≤ 1/2 and 0 ≤ v≤ 1/2

    p_(u,v)    = ∫_u^1∫_v^1 p_X(t_1,t_2)p_X(t_1-u, t_2-v)dt_1dt_2
       ≥∑_i_1=u/(6h_1)+1^N_1∑_i_2=v/(6h_2)+1^N_2∫_(6i_1-5)h_1^(6i_1-1)h_1∫_(6i_2-5)h_2^(6i_2-1)h_2 p_X(t_1,t_2)p_X(t_1-u, t_2-v)dt_1dt_2
       ≥∑_i_1=u/(6h_1)+1^N_1∑_i_2=v/(6h_2)+1^N_2(2h_1)(2h_2)·9/4·9/4
       ≥9/256

for sufficiently large n. A similar calculation holds for all |u|≤ 1/2 and |v|≤ 1/2. Therefore, Condition (c) holds with δ_0 = 1/2 and U_δ≡[-1,1]^2.












































 §.§ Proof of Proposition <ref>


We employ an iterative usage of pairwise difference. Under the regular design, we have Y_i_1,…,i_d = ∑_k=1^d f_k(i_k / n^1/d) + σε_i_1,…,i_d for (i_1,…,i_d)∈ [n^1/d]×…×[n^1/d], where we assume without loss of generality that n^1/d is an even integer. Let m  n^1/d/2 and define I{(1,2),…, (2m-1, 2m)} with cardinality m. For all index pairs (i^(1)_k, i^(2)_k)∈I, k∈[d], we have

    Y_(i^(1)_1,i^(2)_1),…,(i^(1)_d,i^(2)_d)   ∑_j_k∈{i_k^(1),i_k^(2)},k∈[d] Y_j_1,…, j_d(-1)^∑_k=1^d 1{j_k = i_k^(1)}
        = ∑_j_k∈{i_k^(1),i_k^(2)},k∈[d]σε_j_1,…, j_d(-1)^∑_k=1^d 1{j_k = i_k^(1)}.

Clearly, we have *Y_(i^(1)_1,i^(2)_1),…,(i^(1)_d,i^(2)_d) = 0 and *Y_(i^(1)_1,i^(2)_1),…,(i^(1)_d,i^(2)_d) = 2^dσ^2. More importantly, the newly formed data sequence {Y_(i^(1)_1,i^(2)_1),…,(i^(1)_d,i^(2)_d)}_(i_k^(1),i^(2)_k)∈I, k∈[d] with cardinality m^d = n/2^d is i.i.d. with mean 0 and variance 2^dσ^2. Therefore, by defining Y̅ to be the average of this newly formed data sequence and σ^2 to be 

    σ^2 1/n∑_(i^(1)_k, i^(2)_k)∈I, k∈[d]*Y_(i^(1)_1,i^(2)_1),…,(i^(1)_d,i^(2)_d) - Y̅^2,

we clearly have *σ^2 - σ^2^2 ≲ n^-1 for some absolute positive constant C under a finite fourth moment assumption, which is clearly not improvable. Thus the proof is complete.




 §.§ Proof of Proposition <ref>


Throughout the proof, C represents a positive constant that only depends α and C_0. Following the argument before the statement of Proposition <ref>, define 

    ε_i^(ℓ)∑_k∈[d], k≠ℓf_k(X_i,k) + ε_i

and its variance

    σ^2_(ℓ)∑_k∈[d],k≠ℓ f^2_k(X_i,k) + σ^2

for all ℓ∈[d]. Clearly, under the mutual independence of the components of (X_i,1,…, X_i,d), it holds that ε_i^(ℓ) = 0 and ε_i^(ℓ) is independent of f_ℓ(X_i,ℓ). For each ℓ∈[d], by viewing the model equivalently as Y_i = f_ℓ(X_i,ℓ) + ε^(ℓ)_i for i∈[n] and then applying the univariate kernel smoother defined in (<ref>), which renders an estimator which we denote as σ^2_(ℓ), we obtain by Theorem <ref>

    *σ^2_(ℓ) - σ^2_(ℓ)^2 = *σ^2_(ℓ) - ∑_k∈[d],k≠ℓ f_k^2(W_i,k) - σ^2^2 ≤ Cn^-8α_ℓ/4α_ℓ+1.

Moreover, letting Y̅ be the average of {Y_i}_i=1^n and σ^2_Y be the sample variance estimator, that is, σ^2_Y ∑_i=1^n (Y_i - Y̅)^2/n, it holds that

    *σ^2_Y - (Y)^2 = *σ^2_Y - ∑_k∈[d] f_k^2(X_i,k) - σ^2^2 ≤ Cn^-1.

Since σ^2 = ∑_ℓ=1^d σ^2_(ℓ) - (d-1)(Y), thus by defining σ̂^2 ∑_ℓ=1^d σ^2_(ℓ) - (d-1)σ^2_Y, we obtain that

    *σ^2 - σ^2^2    = **σ^2_(1) - σ^2_(1) + … + *σ^2_(d) - σ^2_(d) + (d-1)*σ^2_Y - (Y)^2
       ≤ C*n^-8α_min/4α_min+1 + n^-1.

This completes the proof.














 §.§ Proof of Proposition <ref>


For simplicity, we only prove the case with two additive components f(X) and g(W) which are α- and β-Hölder smooth, respectively. Throughout the proof, C represents a generic fixed positive constant that only depends on α,β, C_0 and the joint distribution of (X,W). Denote the marginal distribution of X and W as F_X and F_W. Since the transition boundary for both α and β is 1/4, we may assume without of loss of generality that 0 < α, β < 1. As a result, since F^-1_X and F^-1_W are both Lipschitz with fixed positive constants, f̅ f∘ F^-1_X and g̅ g∘ F^-1_W are still α- and β-Hölder smooth. With a standard wavelet expansion (cf. Proposition 2.5 in <cit.>), we can write the model equivalently as

    Y_i = f̅_1(U_1,i) + ∑_j=1^2^J_1ψ_j(U_1,i)γ_1,j + g̅_1(U_2,i) + ∑_j=1^2^J_2φ_j(U_2,i)γ_2,j + σε_i,

where {ψ_j}_j=1^∞ and {φ_j}_j=1^∞ are two sets of orthonormal wavelet basis (with respect to the Lebesgue measure on [0,1]), {U_1,i}_i=1^n = {F_X(X_i)}_i=1^n and {U_2,i}_i-=1^n = {F_W(W_i)}_i=1^n are two uniform [0,1] sequences, and f̅_1_∞≤ C(2^-α J_1) and g̅_1_∞≤ C(2^-β J_2). Define U_i  (ψ_1(U_1,i), …, ψ_2^J_1(U_1,i), φ_1(U_2,i),…, φ_2^J_2(U_2,i)) as the new feature vector of length 2^J_1 + 2^J_2. Without loss of generality, we assume U_i = 0 (a mean shift does not affect the estimation of variance) and let Σ(U_i). Without loss of generality, we can assume Σ is strictly positive definite (otherwise we can orthogonalize with respect to the linear span of (ψ_1(U_1,i), …, ψ_2^J_1(U_2,i), φ_1(U_1,i),…, φ_2^J_2(U_2,i)) in (<ref>) below), and thus it holds that 

    Y_i = V_i^⊤γ + f̅_1(U_1,i) + g̅_1(U_2,i) + σε_i,

where γΣ^1/2(γ_1,1,…,γ_1,2^J_1,γ_2,1,…,γ_2,2^J_2), and V_iΣ^-1/2U_i.

We now calculate the bias and variance of the estimator σ^2_proj,add defined as

    σ^2_proj,add1/n-1∑_i=1^n(Y_i - Y̅)^2 - n 2^-1∑_i<jY_iY_jV_i^⊤V_i.

Direct calculation shows that

    *σ^2_proj,add - σ^2   = *f̅_1^2(U_1) + g̅_1^2(U_2) + 2*f̅_1(U_1)g̅_1(U_2) - **(f̅_1(U_1) - g̅_1(U_2))V_2^2
       ≲f̅_1^2(U_1) + g̅_1^2(U_2) + **f̅_1(U_1)V_2^2 + **g̅_1(U_2)V_2^2
       ≲ 2^-2α J_1 + 2^-2β J_2 +  **f̅_1(U_1)V_2^2 + **g̅_1(U_2)V_2^2.
 













Moreover, we have

    **f̅_1(U_1)V_2^2 = sup_a≤ 1**f̅_1(U_1)V^⊤a^2 ≤f̅_1^2(U_1)·sup_a≤ 1*a^⊤(VV^⊤)a≲ 2^-2α J_1,

where the last inequality again follows by the identity covariance of V. Similarly, it holds that **g̅_1(U_2)V_2^2 ≲ 2^-2β J_2. We therefore conclude that the bias of σ^2_proj,add is smaller than the order 2^-2α J_1 + 2^-2β J_2. 

Next, we calculate the variance of σ^2_proj,add. For this, it suffices to upper bound the variance of ∑_i=1^n (Y_i - Y̅)^2/(n-1) and n 2^-1∑_i<jY_iY_jV_i^⊤V_j. The first variance is clearly of the order n^-1 under the boundedness of f(·) and g(·) and the fact that ε_i^4≤ C_ε, so we focus on the second variance. Direct calculation shows that

    *∑_i<jY_iY_jV_i^⊤V_j
       = ∑_i<j,i^'<j^'**Y_iY_jY_i^'Y_j^'(V_i^⊤V_j)(V_i^'^⊤V_j^') - *Y_iY_jV_i^⊤V_j*Y_i^'Y_j^'V_i^'^⊤V_j^'.
 
When i,j,i^',j^' take four different values, the above summand is clearly 0. When they take 3 values (i = i^' < j < j^'), denoting z_i f̅_1(U_1,i) + g̅_1(U_2,i), i∈[n], we have

    *Y_i^2Y_jY_j^'(V_i^⊤V_j)(V_i^⊤V_j^')
        = *(z_i + V_i^⊤γ + σε_i)^2(z_j + V_j^⊤γ + σε_j)(z_j^' + V_j^'^⊤γ + σε_j^')(V_i^⊤V_j)(V_i^⊤V_j^').

Next, we expand the above display and upper bound each term individually. For simplicity, we only show the calculation for the following dominating term and the other terms follow similarly. 

    *(V_i^⊤γ)^2(V_j^⊤γ)(V_j^'^⊤γ)(V_i^⊤V_j)(V_i^⊤V_j^')
       = *(V_j^⊤γ)V_j^⊤*(V_i^⊤γ)^2V_iV_i^⊤V_j^'(V_j^'^⊤γ)
       = *γ^⊤V_jV_j^⊤*(V_i^⊤γ)^2V_iV_i^⊤*V_j^'V_j^'^⊤γ
       = γ^⊤*(V_i^⊤γ)^2V_iV_i^⊤γ
       =  (V_i^⊤γ)^4,

where in the second equality we use the independence of V_i,V_j,V_j^' and in the third equality we use the fact that (V_i^⊤V_i) = (V_i) = 𝐈_L since V_i = Σ^-1/2U_i = 0, where L 2^J_1 + 2^J_2. Moreover, by definition, it holds that *V_i^⊤γ = *f(X_i) + g(W_i) - z_i≤f_∞ + g_∞ + |z_i| ≤ C due to the boundedness of f(·),g(·) and |z_i|. This concludes that when i,j,i^',j^' take three different values, the summand in (<ref>) can be upper bounded by a fixed constant.
When they take two different values (i = i^', j = j^'), we have

    *Y_i^2Y_j^2(V_i^⊤V_j)^2 = *(z_i + V_i^⊤γ + σε_i)^2(z_j + V_j^⊤γ + σε_j)^2(V_i^⊤V_j)^2.

We again upper bound the dominating term in the expansion of the above display.

    *(V_i^⊤γ)^2(V_j^⊤γ)^2(V_i^⊤V_j)^2≤ C(V_i^⊤V_j)^2 = C(𝐈_L) = CL,

where we again use the fact that (V_i^⊤γ) and (V_j^⊤γ) are bounded by a fixed constant. Putting together the pieces, we obtain that

    *σ^2_proj,add≤ C(n^-1 + n^-4(n^3 + n^2L)) = Cn + 2^J_1 + 2^J_2/n^2.

Optimal choice of 2^J_1≍ n^2/(4α+1) and 2^J_2≍ n^2/(4β +1) then gives the desired error bound.




 §.§ Proof of Proposition <ref>


We use H_* as a shorthand for H_δ^*. Fix any α≥α_*. Define the oracle bandwidth 

    h^*
    n^-1,    α > 1/4,
    max*h∈H_*: h^2α≤ c√(log n/(n^2h)),     0<α≤ 1/4,

 for some positive constant c to be specified later. When 0<α≤ 1/4, h^* is taken to be n^-2/(4α+1) if the set being maximized is empty. If not, then it holds that (2h^*)^2α> c√(log n/(2n^2h^*)), and thus (h^*)^2α≍√(log n/(n^2h^*)), or h^*≍ (log n/n^2)^1/(4α+1).

We first prove that with high probability, we have h_δ_*≥ h^*. For this, we have

    *h_δ_*<h^*   ≤*∃ h∈H_*, h≤ h^*, *σ^2(h) - σ^2(h^*)≥τ√(log n/(n^2h))
       ≤∑_h∈H_*,h≤ h^***σ^2(h) - σ^2(h^*)≥τ√(log n/(n^2h))
       ≤∑_h∈H_*,h≤ h^***σ^2(h) - σ^2≥τ/2√(log n/(n^2h)) + *H_*·**σ^2(h^*) - σ^2≥τ/2√(log n/(n^2h^*)).

We now upper bound each probability in the above summation for any h≤ h^*. As in the proof of Theorem <ref>, denote the two U-statistics on the numerator and denominator of σ^2(h) as U_1,U_2, with corresponding mean values θ_1,θ_2. That is, 

    θ_1*K_h(X_i - X_j)(Y_i-Y_j)^2/2 and θ_2  K_h(X_i-X_j).

Define the “good" event E*U_2≥θ_2/2 and E^c as its complement, then it holds that

    **σ^2(h) - σ^2≥τ/2√(log n/(n^2h))⋂E
       ≤**U_1 - U_2σ^2≥τ/θ_2√(log n/(n^2h))
       = **(U_1-θ_1) + (θ_1 - θ_2σ^2) + (U_2-θ_2)≥τ/θ_2√(log n/(n^2h))
       ≤**(U_1 - θ_1)≥τ/4θ_2√(log n/(n^2h)) + **(U_2 - θ_2)≥τ/4θ_2√(log n/(n^2h)),

where the last inequality follows from the fact h≤ h^* and the bound *θ_1 - θ_2σ^2≲ h^2(α∧ 1) calculated in the proof Theorem <ref>. By choose u≍log n and v≍log n/(nh) in Lemma <ref> and τ to be sufficiently large, it holds that

    **U_2 - θ_2≥τ/4θ_2√(log n/(n^2h))≲ n^-C

for arbitrarily large C. Furthermore, for sufficiently large τ and η in Lemma <ref> below, choosing the same u and v yields that

    **U_1 - θ_1≥τ/4θ_2√(log n/(n^2h))≲ n^-C

for arbitrarily large C. This, combined with the calculation 

    *E^c≲exp(-θ_2^2n/16)+exp(-θ_2^2n^2h/16)

in the proof of Theorem <ref>, concludes that (E^c)≲ n^-C, where E*h_δ_*<h^*. Therefore, we have

    *σ_adapt^2 - σ^2^2    ≲**σ_adapt^2 - σ^2^21*E^c + **σ_adapt^2 - σ^2(h^*)^21*E + **σ^2(h^*) - σ^2^2
       ≲ n^-C +  log n/n^2h^* + *n^-1 + (h^*)^4(α∧ 1)+ (n^2h^*)^-1
       ≲*log n/n^2^4α/(4α+1) + n^-1.

This completes the proof.




 §.§ Proof of Proposition <ref>




Note that the desired result is equivalent to the following statement:

    inf_σ^2max*sup_f∈Λ_α_1(C_F)
    σ^2≤ C_σ,_(X,ε)∈P^adapt_cv, (X,ε)*(σ^2 - σ^2)/ϕ_n,α_1^2,sup_f∈Λ_α_2(C_F)
    σ^2≤ C_σ,_(X,ε)∈P^adapt_cv, (X,ε)*(σ^2 - σ^2)/ϕ_n,α_2^2≥ c

for sufficiently large n and sufficiently small c. By applying Lemma <ref> with A={α_1,α_2} with α_*≤α_1<α_2, it suffices to lower bound the adaptive minimax rate under measure  defined therein. More precisely, we will prove that for n≥ n_0 with some sufficiently large n_0,

    inf_σ^2max*sup_f∈Λ_α_1(C_F)
    σ^2≤ C_σ,_(X,ε)∈P^adapt_cv, (X,ε)_*(σ^2 - σ^2)/ϕ_n,α_1^2,sup_f∈Λ_α_2(C_F)
    σ^2≤ C_σ,_(X,ε)∈P^adapt_cv, (X,ε)_*(σ^2 - σ^2)/ϕ_n,α_2^2 > c
 
for some sufficiently small positive c. In order to show this, we will prove that, for any n≥ n_0 and any estimator σ^2, if 

    sup_σ^2≤ C_σsup_f∈Λ_α_2(C_F)sup__(X,ε)∈P^adapt_cv, (X,ε)_*(σ^2 - σ^2)/ϕ_n,α_2^2 ≤ c,

then it holds that

    sup_σ^2≤ C_σsup_f∈Λ_α_1(C_F)sup__(X,ε)∈P^adapt_cv, (X,ε)_*(σ^2 - σ^2)/ϕ_n,α_1^2 > c.

If α_2>1/4, then ϕ_n,α_2≍ n^-1/2, and we can choose a sufficiently small c such that (<ref>) never holds. Therefore, in what follows, we will assume α_*≤α_1<α_2≤ 1/4, in which case ϕ_n,α_i≍ (log n/n^2)^2α_i/(4α_i+1) for i=1,2. 

We will now apply Lemma <ref>. To this end, we adopt a two-point method and introduce the two probability measures _0 and _1 (conditioning on a fixed realization of m∼Poi(2n)) as follows. Introduce

    h_n = c*log n/n^2^1/(4α_1+1) and  N N_n h_n^-1,

where c is some sufficiently small constant tuned such that N is a positive integer. Let Q be a discrete distribution that takes value 0 with probability 1/2, -1 with probability 1/4 and 1 with probability 1/4. It then can be readily checked that m_1(Q) = 0 and m_2(Q) = 1/2.



   Choice of ε: Under H_0, let ε∼ (1+h_n^2α_1/2)^-1/2((h_n^α_1 Q)*N(0,1)). Under H_1, let ε∼N(0,1). 

   Choice of σ^2: Under H_0, let σ^2 = 1 + h_n^2α_1/2. Under H_1, let σ^2 = 1.

  Choice of X: Under both H_0 and H_1, let X be uniformly distributed on the union of the intervals [(6i-5)h_n, (6i-1)h_n] for i∈[N].

   Choice of f(·): Under H_0, let f≡ 0. Under H_1, let f take the value h_n^α_1 r_i on [(6i-5)h_n, (6i-1)h_n], where {r_i}_i=1^N are N i.i.d. variables with law Q. 

Clearly, by the boundedness of Q, H_0 belongs to the model class indexed by the smoothness index α_2, and H_1 belongs to the model class indexed by α_1. Moreover, the absolute difference in σ^2 under H_0 and H_1 is lower bounded by the order h_n^2α_1≍ (log n/n^2)^2α_1/(4α_1+1).

Denote p_0 and p_1 as the densities of _0 and _1. Define f_max2/(1-4α_1)+1 and let d_i be the number of X's that fall into [(6i-5)h_n, (6i-1)h_n] for each i∈ [N]. Consider the following event 

    E*{m≤ 3n}⋂{max_1≤ i≤ Nd_i ≤ f_max}.

Note that under both _0 and _1, the sequence {d_i}_i=1^N are i.i.d. Poisson variables with mean 2n/N. Thus, a standard Poisson tail estimate and Lemma <ref> imply that the event E has asymptotic probability 1 under both _0 and _1. Next, we calculate the χ^2 distance between _0 and _1 conditioning on the event E. First, we have

    ∫p_1^2/p_01*E   = ∫ p(m)p(x_1,…,x_m)1*E∫∏_j=1^Np_1,j^2/p_0,j
       = ∫ p(m)p(x_1,…,x_m)1*E∫∏_j=1^N *1 + χ^2*p_1,j, p_0,j,

where p(m) and p(x_1,…,x_m) are the pmf and pdf of m and {X_i}_i=1^m under both _0 and _1, and p_0,j is the conditional density of all those Y_i's with corresponding X_i∈[(6j-5)h_n, (6j-1)h_n] and similarly for p_1,j. 

We now upper bound each χ^2(p_0,j, p_1,j), where we assume that there are d_j X_i's that belong to [(6j-5)h_n, (6j-1)h_n]. Write d instead of d_j for short. Here, d is a random variable that only depends on m and {X_i}_i=1^m, and on the event E, we have d≤ f_max. Clearly, if d = 0 or 1, then χ^2(p_0,j, p_1,j) = 0. Assume d≥ 2. Assume for simplicity that the d data points are y_1,…,y_d. Then, by definition of _0, we have

    p_0,j≥ (1/2)^f_maxφ(y_1)…φ(y_d)≥ cφ(y_1)…φ(y_d).

On the other hand, we have by direct calculation 

    ∫p_0,j^2/(φ(y_1)…φ(y_d))    = _s_1,…,s_d∼ h^α_1 Q
    s_1,…,s_d∼ h^α_1 Qexp*∑_i=1^d s_is_i,
    ∫p_1,j^2/(φ(y_1)…φ(y_d))    = _t,t∼ h^α_1 Qexp(dtt),
    ∫p_0,jp_1,j/(φ(y_1)…φ(y_d))    = _t,s_1,…,s_d∼ h^α_1 Qexp*∑_i=1^d ts_i.

We therefore conclude that

    χ^2(p_0,j, p_1,j)    ≲∑_k=1^∞1/k!∑_1≤ i_1,…,i_k≤ d*_s_1,…,s_d∼ h^α_1 Q*s_i_1… s_i_k - _t∼ h^α_1 Qt^k^2
       ∑_k=1^∞1/k!∑_1≤ i_1,…,i_k≤ dΔ_i_1,…,i_k^2.

For any k≥ 1, if (i_1,…,i_k) are all identical, then Δ_i_1,…,i_k = 0. More generally, if (i_1,…,i_k) take ℓ different values, where ℓ≤ d≤ f_max on the event E, there are dℓ ways of choosing ℓ different values among [d], and there are a total of k-1ℓ-1 ways to distribute ℓ values in (i_1,…,i_k), thus we obtain the estimate

    ∑_k=1^∞1/k!∑_1≤ i_1,…,i_k≤ dΔ_i_1,…,i_k^2    ≤∑_k=2^∞1/k!∑_ℓ=2^k h^2α_1 kdℓk-1ℓ-1
       = ∑_ℓ=2^f_maxdℓ∑_k=ℓ^∞1/k!h^2α_1 kk-1ℓ-1.

For each 2≤ℓ≤ f_max, by the Stirling's formula, we have

    ∑_k=ℓ^∞1/k!h^2α_1 kk-1ℓ-1≲∑_k=ℓ^∞h^2α_1 k/k^k+1/2e^-k(ek/ℓ)^k ≲∑_k=ℓ^∞(e^2h^2α_1)^k/k^1/2ℓ^k≲ h^2α_1ℓ.

Next, using the trivial bound dℓ≤ d^2f_max^ℓ-2, we obtain that χ^2*p_0,j, p_1,j≲ d^2h^4α_1 if d≥ 2. This implies that

    ∫p_1^2/p_01*E≲∫ p(m)p(x_1,…,x_m)∏_j=1^N *1+χ^2_j,

where χ^2_j = 0 for d_j = 0,1 and χ^2_j≤ d_j^2h^4α_1 for d_j≥ 2. Next, 

    ∏_j=1^N *1+χ^2_j = 1 + ∑_j=1^N χ^2_j + ∑_1≤ i< j≤ Nχ^2_iχ^2_j + … + ∑_1≤ i_1<…<i_N≤ Nχ^2_i_1…χ^2_i_N.

Consider the kth term in the above display. Note that on the event E, we have ∑_j=1^N d_j≤ 3n, thus there are at most (3n/2) j's with d_j≥ 2. Then, using the estimate n k≤ n^k/k!, we have

    ∑_1≤ i_1<…<i_k≤ Nχ^2_i_1…χ^2_i_k≤(3n/2) kh^4kα_1d^2_i_1… d^2_i_k≤C^kn^k/k!h^4kα_1d^2_i_1… d^2_i_k.

We therefore conclude that

    ∫p_1^2/p_01*E≤_d_1,…,d_N∼Poi(2n/N)*1+∑_k=1^NC^kn^kh^4α_1 k/k!d_1^2… d_k^2≤ 1 + ∑_k=1^∞(clog n)^k/k! = n^c

for some sufficiently small constant c. Then, by Lemma <ref> and the calculation

    ε√(I_E)/Δ≲n^c/2(log n/n^2)^4α_2/(4α_2+1)/(log n/n^2)^4α_1/(4α_1+1)→ 0,

where Δ,I,ε are defined as in Lemma <ref>, we conclude that for sufficiently large n and any considered estimator σ^2, if (<ref>) holds, then (<ref>) will follow. This shows that, even over two smooth classes α∈{α_1,α_2}, the adaptive minimax rate can be no faster than ϕ_n,α.









 §.§ Supporting Lemmas



Suppose f∈Λ_α(C_F) and σ^2≤ C_σfor some fixed constants C_F,C_σ, and the joint distribution of (X,ε) satisfies the conditions in P_mcv,(X,ε). Then, the U-statistic U_1 defined in the proof of Proposition <ref> satisfies 

    *U_1 - θ_1^2 ≤ C(n^-1∨ n^-2(h_1h_2)^-1)

for some positive constant C that only depends on M_K,M_K,α,C_F,C_σ,C_0,C_ε.


Denote g as the kernel of U_1, that is,

    g(D_i, D_j)  K_h_1(X_i,1- X_j,1)K_h_2(X_i,2-X_j,2)(Y_i - Y_j)^2/2,   D_i  (X_i, ε_i)^⊤.

Then, it holds that

    (U_1) = n 2^-1∑_i<j,i^' < j^'**g(D_i, D_j) - θ_1*g(D_i^', D_j^') - θ_1.

When i,j,i^', j^' take four different values, the expectation is zero. Using a similar argument as in the proof of Lemma <ref>, when i,j,i^', j^' take three different values, it holds that

    **g(D_i, D_j) - θ_1*g(D_i^', D_j^') - θ_1 = O(1).

When they take two different values,

    **g(D_i, D_j) - θ_1*g(D_i^', D_j^') - θ_1 = O((h_1h_2)^-1).

We therefore conclude that

    (U_1) ≲n^3 + n^2(h_1h_2)^-1/n^4≍ n^-1 + n^-2(h_1h_2)^-1.

This completes the proof.









































































































Suppose h_1h_2 ≳ n^-(2-δ) for some 0 < δ < 2, and the joint distribution of (X,ε) satisfies the conditions in P_mcv,(X,ε) Then, for any u, v> 0, the U-statistic U_2 defined in the proof of Theorem <ref> satisfies 

    **U_2 - θ_2≥ C(v^1/2n^-1/2 + u^1/2n^-1(h_1h_2)^-1/2)≤ C(exp(-u) + exp(-v))

for sufficiently large n and

    *U_2 - θ_2^2 ≤ C(n^-1∨ n^-2(h_1h_2)^-1),

where C is some positive constant that only depends on M_K,M_K,α,C_0. 



The proof is similar to that of Lemma <ref>. In the application of Lemma <ref>, the five quantities are of the order B_1 ≲ 1, B_2≲ (h_1h_2)^-1, B_3≤ n^1/2(h_1h_2)^-1/2, ν_1^2≲ 1, ν_2^2 ≲ (h_1h_2)^-1. Therefore, for any u,v> 0, it holds that

    (|U_2 - θ_2| ≥ a_1v^1/2 + a_2v + b_1u^1/2 + b_2u + b_3u^3/2 + b_4u^2) ≤ C(exp(-v) + exp(-u)),

where a_1 ≲ n^-1/2, a_2 ≲ n^-1, b_1 ≲ n^-1(h_1h_2)^-1/2, b_2≲ n^-1, b_3 ≲ n^-3/2(h_1h_2)^-1/2, b_4 ≲ n^-2(h_1h_2)^-1. Under the condition that h_1h_2 = Ω(n^-(2-δ)) for some δ > 0 and n is sufficiently large, the dominant terms in the above inequality are a_1 and b_1, that is,

    n^-1/2∨ n^-1(h_1h_2)^-1/2.

This proves the first part of the theorem. The expectation version follows by Lemma <ref>. 




Suppose the conditions of Proposition <ref> hold. For the U-statistic U_1 defined therein, suppose h≳ n^-(2-δ) for some 0≤δ≤ 2. Then, for any η>0, there exists some positive constant C = C(M_K,M_K, α, η) such that

    **U_1 - θ_1≥ C(v^1/2n^-1/2 + u^1/2n^-1h^-1/2)≤ C*exp(-u) + exp(-v) + n^-η.



Denote g_1 as the kernel of U_1, and consider its truncated version of defined as

    g̅_1(D_i,D_j) 1/2hK*/h*(f(X_i) - f(X_j)) + ε_ij^21*|ε_i|≤κ_n1*|ε_j|≤κ_n,

where D_i = (X_i,ε_i) and κ_n is some truncation parameter satisfying κ_n↑∞ as n→∞ to be specified later. We first consider the concentration of g̅_1 around its mean value θ̅_1*g̅(D_i,D_j). For this, we will make use of Lemma <ref> by upper bounding the 5 quantities B_1,B_2,B_3, ν_1^2, ν_2^2 therein. 

For B_1, denoting g_1(D) = *g̅_1(D,D_j)| D, it holds that

    g_1(D_i)    = *1/2hK*/h((f(X_i) - f(X_j)) + ε_ij)^21*|ε_i|≤κ_n1*|ε_j|≤κ_n| D_i
       ≲*1/hK*/h(|X_i-X_j|^2α + ε^2_ij)1*|ε_i|≤κ_n1*|ε_j|≤κ_n| D_i
       ≲*1/hK*/h||^2α| X_i + *1/hK*/hε_ij^21*|ε_i|≤κ_n| X_i,ε_i.

For the first term, we have

    *1/hK*/h||^2α| X_i = ∫1/hK*u-X_i/h|u-X_i|^2αp_X(u)du
       = ∫ K(v)|vh|^2αp_X(X_i + vh)dv≤sup_u∈p_X(u)h^2α∫ K(v)|v|^2αdv ≲ h^2α.

For the second term, we obtain similarly that 

    *1/hK*/hε_ij^21*|ε_i|≤κ_n| X_i,ε_i   ≲ (ε_i^2+σ^2)1*|ε_i|≤κ_n≲κ_n^2 + σ^2.

Putting together the pieces and using the fact that κ_n↑∞ as n→∞, we obtain that 

    B_1 = g_1_∞≲κ_n^2.

Moreover, with similar analysis, it can be readily checked that ν_1^2≲ 1.

For B_2, it holds that 

    *g̅_1(D_i,D_j)   ≲1/hK*/h*(f(X_i)-f(X_j))^2 + ε_ij^21*|ε_i|≤κ_n1*|ε_j|≤κ_n
       ≲1/hK*/h*||^2α + ε_ij^21*|ε_i|≤κ_n1*|ε_j|≤κ_n
       ≲1/hK*/h*/h^2αh^2α + 1/hK*/hκ_n^2 ≲1/hκ_n^2.

We therefore conclude that B_2 = g̅_1_∞≲ h^-1κ_n^2.

For B_3, we have

    B_3^2    = nsup_D_i*g̅_1^2(D_i,D_j)| D_i
       = nsup_D_i*1/h^2K^2*/h(f(X_i) - f(X_j) + ε_ij)^41**ε_i≤κ_n1**ε_j≤κ_n| X_i,ε_i
       ≲nM_K/h*1/hK*/h*(f(X_i) - f(X_j))^4 + ε_ij^41**ε_i≤κ_n1**ε_j≤κ_n| X_i,ε_i
       ≲nM_K/h*1/hK*/h*^4α| X_i,ε_i + nM_K/h*1/hK*/hε_ij^41**ε_i≤κ_n1**ε_j≤κ_n| X_i,ε_i.

Now, using similar calculation in the analysis of B_1, it holds that

    *1/hK*/h*^4α| X_i,ε_i≲ h^4α,
       *1/hK*/hε_ij^41**ε_i≤κ_n1**ε_j≤κ_n| X_i,ε_i≲κ_n^4.

Putting together the pieces, we conclude that B_3 ≲κ_n^4nh^-1. 

Lastly, for ν_2^2, we have

    ν_2^2    = *g̅_1^2(D_i,D_j)≲*1/h^2K^2*/h*f(X_i) - f(X_j) + ε_ij^4
       ≲ h^-1*1/hK*/h*f(X_i) - f(X_j)^4 + M_K/h*1/hK*/hε_ij^4
       ≲ h^-1(h^4α + (ε_i^4))≲ h^-1.


Define U̅_1 to be the U-statistic generated by the kernel g̅_1, that is, U̅_1 n 2^-1∑_i<jg̅_1(D_i,D_j), and define θ̅_1 to be the mean value *g̅_1(D_i,D_j). Define the event E**ε_i≤κ_n  for all i∈[n]. Then, we have for any t≥ 0,

    **U_1 - θ_1≥ t   = **U̅_1 - θ_1≥ t ⋂E + *E^c
       ≤**U̅_1 - θ̅_1≥ t - *θ_1 - θ̅_1 + *E^c.

For the first term, we have by Lemma <ref> that, for any given u,v> 0, it holds that

    **U̅_1 - θ̅_1≥ a_1v^1/2 + a_2v + b_1u^1/2 + b_2u + b_3u^3/2 + b_4u^2≤ C(exp(-u) + exp(-v)),

where a_1 ≲ n^-1/2, a_2 ≲κ_n^2/n, b_1 ≲ n^-1h^-1/2, b_2≲ n^-1κ_n^2, b_3≲κ_n^2n^-3/2h^-1/2, b_4≲ n^-2h^-1κ_n^2. Choosing κ_n = κ√(log n) for some sufficiently large constant κ, then as long as h = Ω(n^-(2-δ)) for some δ > 0, then the dominant terms in the above inequality are a_1 and b_1, that is, 

    n^-1/2∨ n^-1h^-1/2.

Therefore, Lemma <ref> implies that

    **U̅_1 - θ̅_1≤ C(n^-1/2∨ n^-1h^-1/2).

Now we calculate the difference between θ_1 and θ̅_1. By definition, we have

    *θ_1 - θ̅_1   = **1/hK*/h(f(X_i) - f(X_j) + ε_ij)^21*|ε_i|≥κ_n⋃ |ε_j|≥κ_n
       ≲**1/hK*/h(f(X_i) - f(X_j) + ε_ij)^21*|ε_i|≥κ_n
       ≲**1/hK*/h(f(X_i) - f(X_j))^21*|ε_i|≥κ_n+**1/hK*/hε_ij^21*|ε_i|≥κ_n
       ≲ h^2α(|ε_i|≥κ_n) + *ε_i^21*|ε_i|≥κ_n
       ≲ h^2αn^-κ^2/(2κ_ε^2) + κ_ε^2n^-κ^2/(4κ_ε^2),

where the second line is by symmetry, and in the last line we use the sub-Gaussianity of ε_i. 	Therefore, as long as h = Ω(n^-(2-δ)) for some δ > 0, by choosing κ large enough (depending only on δ,κ_ε,u,v), it holds that 

    *θ_1 - θ̅_1 = o(a_1v^1/2 + a_2v + b_1u^1/2 + b_2u + b_3u^3/2 + b_4u^2).

Lastly, by the sub-Gaussanity of ε_i, it holds that

    (E^c) ≤ n(ε_i≥κ_n) ≲ n^-η

for sufficiently large η by choosing κ correspondingly large enough. This completes the proof.


The following Poissonization lemma reduces the original problem of Proposition <ref> into the case with a random sample size, which facilitates the calculation of χ^2 distance. We introduce some notation. Consider the following experiment: for any given positive integer n, f(·), σ, distribution p_X(·) of X, and distribution p_ε(·) of ε,  


  * generate m∼Poi(2n);

  * generate X_1,…, X_m∼ p_X and ε_1,…,ε_m∼ p_ε;

  * generate Y_i = f(X_i) + σε_i for each i∈[m].

Denote the original experiment and the above experiment as  and , respectively, where we omit the dependence on n, f(·), σ, distribution p_X(·) of X, and distribution p_ε(·) of ε.



Let ϕ_n,α be defined as in Proposition <ref>.
For any fixed α_*>0 and set A⊂[α_*, ∞), the following inequality holds:

    inf_σ^2sup_α∈Asup_σ^2≤ C_σsup_f∈Λ_α(C_F)sup__(X,ε)∈P_cv,(X,ε)_*(σ^2 - σ^2)/ϕ_n,α^2 
       ≥inf_σ^2sup_α∈Asup_σ^2≤ C_σsup_f∈Λ_α(C_F)sup__(X,ε)∈P_cv,(X,ε)_*(σ^2 - σ^2)/ϕ_n,α^2  - 4nC_σ^2exp(-n/6).



Define the event E{m≥ n}, where m∼Poi(2n). Then, a standard tail estimate has (E)≥ 1 - e^-n/6. For the adaptive minimax rate under , we have

    inf_σ^2sup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_*(σ^2 - σ^2)/ϕ_n,α^2
       = inf_σ^2sup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_**(σ^2 - σ^2)/ϕ_n,α^2*1*E + 1*E^c
       = inf_σ^2≤ C_σsup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_**(σ^2 - σ^2)/ϕ_n,α^2*1*E + 1*E^c
       ≤inf_σ^2≤ C_σsup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_**(σ^2 - σ^2)/ϕ_n,α^21*E + 4nC_σ^2·sup_σ^2≤ C_σsup_f(·), _(X,ε)(E^c)
       ≤inf_σ^2sup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_**(σ^2 - σ^2)/ϕ_n,α^21*E + 4nC_σ^2exp(-n/6)
       ≤inf_σ^2=σ^2*{(X_i,Y_i)}_i=1^nsup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_**(σ^2 - σ^2)/ϕ_n,α^21*E + 4nC_σ^2exp(-n/6)
       ≤inf_σ^2sup_α∈Asup_σ^2≤ C_σsup_f(·),_(X,ε)_**(σ^2 - σ^2)/ϕ_n,α^2 + 4nC_σ^2exp(-n/6).

This completes the proof.




The following lemma will also be used in the proof of Proposition <ref>, and is a slight variation of the constrained risk inequality derived in <cit.> (see Theorem 1 therein). We first introduce some notation. Consider some measurable space equipped with a class of probability measures {_θ}_θ∈Θ, where (Θ, d) is a metric space. For each θ∈Θ, let f_θ be the density of _θ with respect to some common dominating measure ν, and denote by _θ the expectation under the measure _θ. For any estimator T of θ, define its risk as

    R_θ R(θ, T) _θ(T - θ)^2 = ∫ (T(x) - θ)^2f_θ(x)ν(dx).

Now, fix two measures _θ_1 and _θ_2, and let E be a measurable set.  Define

    I_E I(θ_1,θ_2, E) = _θ_1(q^2(X)1{X∈E}),

where q(x) f_θ_2(x)/f_θ_1(x).



Let Δ d(θ_1,θ_2) and assume that, for certain estimator T, R(θ_1,T)≤ε^2 and 0 < ε < Δ((_θ_2(E)/√(I_E))∧ 1). Then, 

    R(θ_2,T) ≥Δ^2_θ_2(E)^2*1 - 2ε√(I_E)/_θ_2(E)Δ.



We follow the proof of Theorem 1 in <cit.> by considering the same estimator T therein which minimizes R(θ_2,T) subject to the condition R(θ_1,T)≤ε^2. Then, with ρ defined therein, we have

    T(x) = ρΔ q(x)/1+ρ q(x)

and R(θ_1,T) = ε^2, so that (see Equation (2.6) in the proof of Theorem 1 in <cit.>)

    ε^2 = Δ^2∫*ρ q(x)/1+ρ q(x)^2 f_θ_1(x)ν(dx)

under the condition ε<Δ. Then, by Cauchy-Schwarz, 

    ε√(I_E)   = Δ*∫*ρ q(x)/1+ρ q(x)^2 f_θ_1(x)ν(dx)^1/2*∫ q^2(x)f_θ_1(x)1_E(x)ν(dx)^1/2
       ≥Δ∫ρ q(x)/1+ρ q(x)f_θ_2(x)1_E(x)ν(dx).

Then, under the condition ε√(I_E)≤Δ_θ_2(E), we have

    *Δ_θ_2(E) - ε√(I_E)^2    ≤Δ^2*_θ_2(E) - ∫ρ q(x)/1+ρ q(x)f_θ_2(x)1_E(x)ν(dx)^2
       = Δ^2*∫1/1+ρ q(x)f_θ_2(x)1_E(x)ν(dx)^2
       ≤Δ^2*∫1/1+ρ q(x)f_θ_2(x)ν(dx)^2
       ≤Δ^2∫*1/1+ρ q(x)^2f_θ_2(x)ν(dx)
       = R_θ_2,

where the last equality is true due to Equation (2.5) in the proof of Theorem 1 in <cit.>. The statement then follows from (a-b)^2≥ a^2(1-2b/a) for a,b>0.




Suppose X_1,…,X_N are i.i.d. Poisson variables with mean value n/N, where N = N_n ≥ Cn^1+δ for some positive constant δ and absolute constant C. Then, there exists some positive integer f_max that only depends on δ such that

    *max_1≤ i≤ NX_i≤ f_max→ 1
  
as n→∞.


Let λ n/N. We will show that the above statement holds for f_max = (1+δ)/δ. For each variable X_i and positive integer k, we have

    *X_i ≤ k = ∑_ℓ=0^k λ^ℓ/ℓ!e^-λ = e^-λ(1+∑_ℓ=1^k (n/N)^ℓ/ℓ!).

Therefore, it holds that

    *max_1≤ i≤ NX_i≤ k   = exp(-n)exp*Nlog*1+∑_ℓ=1^k (n/N)^ℓ/ℓ!
       = exp(-n)exp*N∑_m=1^∞1/m(-1)^m-1∑_ℓ_1,…,ℓ_m=1^k 1/ℓ_1!…ℓ_m!(n/N)^ℓ_1+…+ℓ_m.

The exponent in the above display is a polynomial function of (n/N), and clearly the coefficient for (n/N)^1 is 1. Next, we will show next that the coefficients corresponding to (n/N)^ℓ for ℓ=2,…, k are all zero. For simplicity, we will show this for ℓ = k. By Lemma <ref>, the coefficient for (n/N)^k is

    1/k!∑_m=1^k (-1)^m-11/m∑_ℓ=1^m ℓ^k(-1)^m-ℓmℓ = 1/k!∑_ℓ=1^k (-1)^ℓℓ^k-1∑_m=ℓ^k m-1ℓ-1 = 1/k!∑_ℓ=1^k(-1)^ℓℓ^k-1mℓ = 0,

where the first identity is by direct calculation, the second is the Hockey-Stick identity, and the third is proved in <cit.>.

Next, we consider the coefficient of (n/N)^p for some general p≥ k+1, which takes the form

    ∑_m=1^∞1/m(-1)^m-1∑_1≤ℓ_1,…,ℓ_m≤ k
    ℓ_1+…+ℓ_m=p1/ℓ_1!…ℓ_m!= ∑_m=1^p 1/m(-1)^m-1∑_1≤ℓ_1,…,ℓ_m≤ k
    ℓ_1+…+ℓ_m=p1/ℓ_1!…ℓ_m!≤∑_m=1^p 1/mm^p/p!
       ≲1/pp^p/p!≲ e^pp^-3/2.

Thus, in view of the fact that N≥ Cn^1+δ, we have

    ∑_p=k+1^∞ (n/N)^p∑_m=1^∞1/m(-1)^m-1∑_1≤ℓ_1,…,ℓ_m≤ k
    ℓ_1+…+ℓ_m=p1/ℓ_1!…ℓ_m!≲∑_p=k+1^∞*en/N^pp^-3/2≲ (en/N)^k+1.

By definition of f_max, we have (n/N)^f_max+1→ 0 as n→∞, thus *max_1≤ i≤ NX_i≤ f_max→ 1 as n→∞.




Fix any positive integer k. Then, for any positive integer 1≤ m≤ k, the following identity holds:

    ∑_1≤ℓ_1,…,ℓ_m≤ k
    ℓ_1+…+ℓ_m=k1/ℓ_1!…ℓ_m! = 1/k!∑_ℓ=1^m (-1)^m-ℓℓ^kmℓ.



We will prove by induction. Denote the LHS by C(m). Suppose the statement holds up to m. Using the identity

    ∑_0≤ℓ_1,…,ℓ_m+1≤ k
    ℓ_1+…+ℓ_m+1 = kk!/ℓ_1!…ℓ_m+1!(m+1)^-k = 1,

we obtain the recursive equation

    C(m+1) + m+1 1C(m) + … + m+1 mC(1) = (m+1)^k/k!.

Therefore, plugging in the equation for C(ℓ), ℓ=1,…,m, we obtain that

    C(m+1)    = (m+1)^k/k! - ∑_ℓ=1^m m+1ℓ/k!∑_j=0^ℓ-1(ℓ-j)^k(-1)^jℓ j
       = 1/k!*(m+1)^k - ∑_ℓ=1^m m+1ℓ∑_j=1^ℓ j^k(-1)^ℓ-jℓ j
       = 1/k!*(m+1)^k - ∑_j=1^m j^k∑_ℓ=j^m (-1)^ℓ-jℓ jm+1ℓ.

Thus it suffices to show that 

    ∑_ℓ=j^m (-1)^ℓℓ jm+1ℓ = (-1)^mm+1 j.

This is indeed true since the LHS equals

    m+1 j∑_ℓ=j^m m+1-jℓ-j(-1)^ℓ = ∑_ℓ=0^m-jm+1 jm-j+1ℓ(-1)^ℓ+j
       = m+1 j(-1)^j*∑_ℓ=0^m-j+1m-j+1ℓ(-1)^ℓ + (-1)^m-j = m+1 j(-1)^m,

which completes the proof.









