Conceptio › Archive › arXiv CS

arXiv CSopen access

Prediction-powered Inference by Mixture of Experts

2026 · arxiv_cs

arXiv CS · Papers · License: Open Access · 2026

Open Source ↗Direct PDF ↓

neural-networks

machine learning, deep learning, neural networks

Prediction-powered Inference by Mixture of Experts Yanwu Gu1 , Linglong Kong2 , Dong Xia1∗ 1

Department of Mathematics, Hong Kong University of Science and Technology

arXiv:2604.27892v1 [stat.ML] 30 Apr 2026

Department of Mathematical and Statistical Sciences, University of Alberta

(May 1, 2026)

Abstract The rapidly expanding artificial intelligence (AI) industry has produced diverse yet powerful prediction tools, each with its own network architecture, training strategy, data-processing pipeline, and domain-specific strengths. These tools create new opportunities for semi-supervised inference, in which labeled data are limited and expensive to obtain, whereas unlabeled data are abundant and widely available. Given a collection of predictors, we treat them as a mixture of experts (MOE) and introduce an MOE-powered semi-supervised inference framework built upon prediction-powered inference (PPI). Motivated by the variance reduction principle underlying PPI, the proposed framework seeks the mixture of experts that achieves the smallest possible variance. Compared with standard PPI, the MOE-powered inference framework adapts to the unknown performance of individual predictors, benefits from their collective predictive power, and enjoys a best-expert guarantee. The framework is flexible and applies to mean estimation, linear regression, quantile estimation, and general M-estimation. We develop non-asymptotic theory for the MOE-powered inference framework and establish upper bounds on the coverage error of the resulting confidence intervals. Numerical experiments demonstrate the practical effectiveness of MOE-powered inference and corroborate our theoretical findings.

Keywords: Semi-supervised Learning, Model Averaging, Ensemble Learning, Variance Reduction.

∗

Xia’s research is partially supported by Hong Kong RGC grant GRF 16303224.

Introduction

Semi-supervised learning (Zhu, 2005; Chapelle et al., 2009) is a machine learning paradigm that aims to integrate both labelled and unlabelled data to improve overall prediction accuracy and inferential efficiency. Suppose we observe i.i.d. labelled data L = {(Xi , Yi )}ni=1 ⊂ X × Y, drawn from a joint ei }N , where Xi and X ei share the same marginal distribution PX,Y , and i.i.d. unlabelled data U = {X i=1 distribution PX . Since unlabelled data are often much cheaper to obtain than labelled data in many applications, the unlabelled sample size typically far exceeds the labelled sample size; that is, N ≫ n. Because the unlabelled data provide abundant information about PX , inferential efficiency can potentially be improved whenever the target parameter depends on the covariate distribution (Chakrabortty and Cai, 2018). More specifically, let θ∗ = θ(PX,Y ) = θ PX , PY |X denote the target parameter. Then, whenever θ depends on PX , one may potentially improve inferential efficiency by leveraging the unlabelled data. For example, Zhang et al. (2019) studied semi-supervised inference for the population mean EY and showed that even a simple linear-model-based correction can outperform the labelled-only sample mean, especially when Y and X are strongly linearly related. When the covariate X is high-dimensional, inferential efficiency for EY can still be improved by incorporating unlabelled data through sparse linear regression methods (Zhang and Bradic, 2022). See also Tony Cai and Guo (2020); Liu et al. (2025); Xu et al. (2025). The effectiveness of semi-supervised inference also depends on the postulated model for PY |X . For instance, Zhang et al. (2019) and Zhang and Bradic (2022) studied linear models and showed that their proposed semi-supervised inference methods achieve lower variance when the linear relationship between Y and X is strong. This observation suggests that the efficiency of semisupervised inference may be further improved when more powerful predictive models for Y |X are available. Deep learning methods, or artificial intelligence more broadly, provide a rich class of such predictive tools. Prediction-powered inference (PPI; Angelopoulos et al., 2023a) is a general semi-supervised inference framework that leverages both unlabelled data and a powerful predictor f : X → Y. The predictor is regarded as powerful in the sense that Y ≈ f (X). For example, this can be formalized by requiring that Var Y − f (X) be much smaller than Var(Y ) in mean value inference; that is, f (X) explains a substantial fraction of the variation in Y . There are two key components in PPI: the imputation of outcomes for the unlabelled data, yielding pseudo-pairs ei , Yei )}N with Yei = f (X ei ), and the rectifier, constructed from the labelled data through the {(X i=1 2

ei ) may be biased, the rectifier serves residuals {Yi − f (Xi )}ni=1 . Because the imputed values f (X to debias the resulting imputation-assisted estimators. PPI is a flexible framework that applies to a variety of classical statistical problems, including mean estimation, linear regression, logistic regression, and general M -estimation. See also Angelopoulos et al. (2023b); Zrnic and Candès (2024). It is often the case that the predictor f (·) is obtained from an off-the-shelf machine learning method trained on a massive dataset whose distribution does not necessarily match PX,Y . Because of the potential distribution shift, there is no guarantee that a predictor f (·) that performs well on the training data will still be powerful on the target labelled data L. This issue is especially pronounced in deep learning methods, where networks are often highly overparameterized (Allen Zhu et al., 2019). Although a preliminary check of Var Y − f (X) may be convenient before applying PPI in some settings (e.g., mean estimation), it can be prohibitively time-consuming in others (e.g., M-estimation). Moreover, the rapidly evolving and expanding AI industry has produced a diverse array of powerful prediction tools, many of which are publicly or commercially accessible. These tools may be built on distinct network architectures tailored to specific data types, practical applications, or domain knowledge, and they often produce different predictions for the same covariate input. Large language models (LLMs), for example, illustrate this diversity well: GPT (Achiam et al., 2023), Claude, Gemini (Team et al., 2023), Qwen (Yang et al., 2025), DeepSeek, and Kimi each have their own strengths, i.e., being an expert in different domains. This raises natural and crucial questions: is it possible to enhance the inference efficiency when multiple predictors {fk (·)}K k=1 are available whose performance on the target labelled data L is unknown or expensive to evaluate? Is it possible to generalize the PPI framework so that it can robustly incorporate additional predictors whenever they are available? Treating each predictor fk (·) as a potentially domain-specific expert, we propose a mixture-ofexperts (MOE) PPI framework. MOE (Jacobs et al., 1991; Dai et al., 2024) is a general framework for predictive modeling in which multiple specialized predictors {fk (·)}K k=1 are combined adaptively to produce stronger predictions. For example, a linear mixture of experts takes the form P ⊤ Fβ (x) = K k=1 βk (x)fk (x), where the weight vector β(x) = (β1 (x), · · · , βK (x)) varies with the input covariate. The central idea of MOE is to decompose a complex prediction problem into simpler subproblems, each handled by an expert with localized competence.

Motivated by the variance reduction principle underlying PPI, our MOE-powered inference framework seeks a mixture of experts Fβ such that PPI achieves the smallest possible variance, Var Y − Fβ (X) . This MOE-PPI framework enjoys several advantages including collective prediction power, best-expert guarantee, and safe expert expansion, which we shall detail in Section 2. Learning the optimal weight function β(x) is challenging and may itself introduce bias and substantial variability, especially when the dimension of X is high. For simplicity, we therefore focus primarily on global weighting, that is, β(x) ≡ β. We discuss localized weighting in Section 5. Note that global weighting reduces MOE to classical model averaging in the statistical literature (Claeskens and Hjort, 2008; Raftery et al., 1997). Our contributions are as follows. We propose a flexible MOE-powered inference framework for prediction-powered inference that leverages the collective predictive power of multiple predictors. We define the oracle mixture of experts as the combination that minimizes the variance of PPI-based estimators. Notably, unlike classical model averaging, which typically focuses on mean squared error (Claeskens and Hjort, 2008), the oracle MOE is designed specifically to minimize variance. This distinction is natural in the PPI framework, since the rectification step already removes bias. The oracle MOE generally achieves smaller variance than PPI based on any single predictor. We estimate the oracle mixture weights by minimizing the sample (co)variance on the labeled data and then construct confidence sets using the resulting sample MOE. Unlike standard PPI-based estimators, which are unbiased, the MOE-powered estimator is generally biased because the oracle MOE must be estimated from the labelled data. Nevertheless, we show that this bias is negligible and, in many applications, is dominated by the standard deviation. Our MOE-powered inference framework is highly flexible and applies to a broad range of problems, including M-estimation, mean estimation, linear regression, and quantile estimation. We develop a non-asymptotic theory for the proposed framework. Under mild conditions, we show that the bias of the MOE-powered estimator converges at the rate O(n−1 ). We also establish a Berry–Esseen bound for the normal approximation of the MOE-powered estimator, which yields non-asymptotic guarantees for the coverage probabilities of the resulting confidence sets. Moreover, the constructed confidence intervals have widths comparable to those based on the oracle MOE. Finally, we present comprehensive numerical experiments demonstrating that the proposed MOE framework generally outperforms standard PPI in terms of inferential efficiency.

The rest of the paper is organized as follows. In Section 2, we introduce the MOE-powered inference framework for general M-estimation, describe the construction of the confidence set, and provide theoretical guarantees for its coverage probability. Section 3 presents specific applications to mean estimation, quantile estimation, linear regression, and logistic regression, accompanied with corresponding theoretical guarantees. Comprehensive simulation studies and real-world data experiments are presented in Section 4. We discuss possible extensions in Section 5. All proofs and technical lemmas are deferred to the Appendix.

MOE-powered Inference for M-Estimation

2.1

Variance reduction by PPI

M-estimation (Geer, 2000) is a unified statistical framework for parametric estimation. Let (X, Y ) ∈ X × Y ⊂ Rd × R be a pair of random variables and ℓθ : X × Y 7→ R be a loss function in θ ∈ Θ ⊂ Rp . Denote the population minimizer θ∗ = arg min Eℓθ (X, Y ),

(1)

θ∈Θ

where we assume Θ is a convex set. Note that θ∗ may represent a set if the minimizer is not unique. For simplicity, we assume that ℓθ (·, ) is continuously differentiable with respect to θ and denote its gradient by gθ (·) := ∇θ ℓθ (·). Under the optimality, we can equivalently view θ∗ as the solution to E gθ (X, Y ) = 0.

(2)

Let L := {(X, Y ), (X1 , Y1 ), · · · , (Xn , Yn )} be a collection of labelled i.i.d. observations and U := e X e1 , · · · , X eN } be a collection of unlabelled i.i.d. observations. Our goal is to construct valid {X, 100(1 − α)% confidence sets for θ∗ . Following the practice in Angelopoulos et al. (2023b), we search for the candidate solutions from a fine-grained grid Θgrid ⊂ Θ. Let f : X 7→ Y be a given “powerful” expert or predictor. Motivated by (2) and inspired by Angelopoulos et al. (2023a), we define the sample gradient and imputed gradient by n

bθ := g

1X gθ (Xi , Yi ) n

eθ,f = and g

i=1

N 1 X e i , f (X ei ) , gθ X N i=1

bθ is unbiased, the imputed one may have a non-zero respectively. While the sample gradient g (negative) bias ∆θ,f := E gθ (X, Y ) − E gθ (X, f (X)), which is referred to as the rectifier in 5

Angelopoulos et al. (2023b). The empirical rectifier is thus defined by b θ,f := 1 ∆ n

n X

gθ (Xi , Yi ) − gθ (Xi , f (Xi )) .

i=1

b θ,f are unbiased, yet their bθ and the rectified imputed gradient g eθ,f + ∆ Both the sample gradient g covariances can be strikingly different. Denote Wθ,Y = Cov gθ (X, Y ) , Wθ,f = Cov gθ (X, f (X)) , and Wθ,Y −f = Cov gθ (X, Y ) − gθ (X, f (X)) , respectively. Then, Wθ,Y bθ = Cov g n

and

b θ,f = Wθ,f + Wθ,Y −f , eθ,f + ∆ Cov g N n

implying that the rectified imputed gradient achieves a smaller covariance if N ≫ n and Wθ,Y −f ≤ Wθ,Y . For brevity, we abuse the notation “≤” for matrices which, depending on the context, may represent comparison w.r.t. their traces or diagonal entries one-by-one. b θ,f , the conventional and PPI-based bθ and n1/2 ∆ Motivated by the asymptotic normality of n1/2 g approaches construct the confidence sets by q n o c θ,Y )jj /n, ∀j ∈ [p] ; Cαcon := θ ∈ Θgrid : (b gθ )j ≤ zα/(2p) (W q n o PPI b c θ,Y −f )jj /n, ∀j ∈ [p] , eθ,f + ∆θ,f j ≤ zα/(2p) (W Cα,f := θ ∈ Θgrid : g

(3)

c θ,Y and W c θ,Y −f represent the sample estimates where the z-critical value zα := Φ−1 (1 − α) and W of Wθ,Y and Wθ,Y −f , respectively. See Section 2.2 for more details. It was shown by AngelopouPPI is at least 1 − α los et al. (2023b) that, under mild conditions, the coverage probability of Cα,f

asymptotically. The advantage of PPI-based approach is reliant on the assumption that the expert is powerful in predicting Y so that Wθ,Y −f ≤ Wθ,Y . This assumption may not hold in practice, or it may be costly to verify. Fortunately, it is common that multiple experts are available, e.g., AI models released by different companies. Our MOE-powered inference framework can leverage the prediction power of multiple experts and outperform the conventional one as long as one expert provides accurate prediction.

2.2

Oracle MOE and its estimate

Let f1 (·), · · · , fK (·) : X 7→ Y be a collection of available experts. We consider the mixture of experts (MOE) through an ensemble operator G (·; β), parameterized by β ∈ B ⊆ Rq , FB := Fβ = G (f1 , · · · , fK ; β) : Fβ (x) ∈ C 2 (B ◦ ), ∀x ∈ X . 6

(4)

Here we assume the set B is convex and G (f1 , · · · , fK , β) is twice continuously differentiable w.r.t. the weight vector β ∈ B. Note that B ◦ denotes the interior if B is closed, and q is a constant much smaller than n, indicating the complexity of the candidate space is limited. The function family FB covers many important examples when choosing q = K: 1). Linear mixture: Fβ (·) =

K k=1 βk fk (·) for β ∈ B = R ;

2). ℓ1 -compressed sensing: Fβ (·) = 3). Soft max: Fβ (·) = log

K k=1 βk fk (·) for B := {β ∈ R : ∥β∥ℓ1 ≤ L};

k=1 exp(βk fk (·))

4). Logistic linear mixture: Fβ (·) = ψ

for β ∈ B = RK ;

k=1 βk fk (·)

for β ∈ B = RK with a sigmoid function

ψ(u) = 1/(1 + e−u ). Inspired by the variance-reduction principle, we call the weight vector that minimizes the variance of PPI-based estimators the oracle MOE weight vector. More precisely, we define β∗ := arg min tr Wθ,Y −Fβ .

(5)

β∈B

For notational brevity, we suppress the dependence of β∗ on θ. Correspondingly, we refer to Fβ∗ as the oracle mixture of experts, which enjoys several benign properties: 1). Collective prediction power: By definition, the oracle MOE attains the smallest variance among all possible mixtures in the class. In particular, its variance is never larger than that of any PPI-based estimator that uses a single expert; that is Wθ,Y −Fβ∗ ≤ Wθ,Y −fk for ∀k ∈ [K]. This variance-reduction guarantee highlights the key advantage of MOE-powered methods: they can exploit the collective predictive power of multiple experts. This advantage is especially appealing when different experts capture domain knowledge that is not widely shared. See the numerical experiments in Section 4. 2). Best-expert guarantee: The oracle MOE automatically achieves variance reduction at least as good as that of the PPI-based estimator built on the best single expert. Unlike standard PPI, the MOE approach does not require prior knowledge of the experts’ relative predictive performance. Put differently, the oracle MOE is robust to expert misspecification in the sense that it can still deliver effective variance reduction as long as at least one expert provides strong predictions. 7

3). Safe expert expansion: The variance of the oracle MOE does not increase as the expert pool expands. As the rapidly evolving AI industry continues to produce new expert models, these models can be safely incorporated into the MOE framework without requiring confidence check in their trustworthiness or predictive quality. We define the unbiased sample covariance matrix by c θ,Y −f := W

n ⊤ 1 X n b θ,f ∆ b⊤ , ·∆ gθ (Xi , Yi ) − gθ (Xi , f (Xi )) gθ (Xi , Yi ) − gθ (Xi , f (Xi )) − θ,f n−1 n−1 i=1

c θ,Y −F ). It remains to establish and estimate the oracle MOE weight by βbn := arg minβ∈B tr(W β non-asymptotic error rate for βbn − β∗ . b θ,f = For any fixed θ, let mθ (β, X, Y ) := gθ (X, Y ) − gθ (X, Fβ (X)). As a result, we have ∆ P c θ,Y −f = Wθ,Y −f . Then we can write Wθ,Y −F = Cov mθ (β, X, Y ) n−1 ni=1 mθ (β, Xi , Yi ) and EW β and define Q(β) = tr Wθ,Y −Fβ , where we suppress the dependence of Q(β) on θ for notational c θ,Y −F . brevity. Similarly, we write Qn (β) := tr SCov mθ (β, X, Y ) := tr W β Assumption 1. Suppose that B ⊂ Rq is compact with a fixed parameter dimension q and a diameter D > 0, and, for any θ ∈ Θ, (X, Y ) ∈ X × Y, the function mθ (β, X, Y ) satisfies: (a) mθ (β, X, Y ) ∈ C 2 (B ◦ ) and ∥mθ (β, X, Y )∥ ≤ U ; (b)

mθ (β1 , X, Y ) − mθ (β2 , X, Y ) ≤ τ1 ∥β1 − β2 ∥, for ∀β1 , β2 ∈ B;

(c) Q(β) − Q(β∗ ) ≥ τ0 ∥β − β∗ ∥2 for ∀β ∈ B; where τ0 , τ1 , U > 0 are absolute constants. Assumption 1 is standard for obtaining high-probability upper bound in empirical risk minimization. The uniform boundedness condition on mθ (β, X, Y ) can be relaxed, e.g., one may instead assume that mθ (β, X, Y ) is sub-Gaussian for any fixed θ and β. The assumption of global strong convexity can also be weakened to a local condition, provided additional requirements hold, for instance, that β∗ is the unique minimizer of Q(β). See Assumption 3 in Section 3.2. We shall explore these refinements in future work. Lemma 1. Suppose Assumption 1 holds and n ≥ C0 (τ12 /τ0 )p q + log(Dn) for a large C0 > 1. There exists a constant C1 > 0 such that −1 τ12 p(q + log(Dn)) 2 2 τ0 log n b ∥βn − β∗ ∥ ≤ C1 U + 2· , n n τ0 8

which holds with probability at least 1 − n−9 .

2.3

MOE-powered confidence sets

Based on the estimated oracle MOE Fβbn (·), we construct the confidence set following the PPI b b denote g b θ,F , respectively. eθ,βbn and ∆ eθ,Fβb and ∆ approach (3). For notational simplicity, let g b θ,βn β n

We then define the MOE-powered confidence set as q n o b b c θ,Y −F eθ,βbn + ∆ ≤ z CαMOE := θ ∈ Θgrid : g W /n, ∀j ∈ [p] , α/(2p) b θ,βn j jj β n

where Θgrid denotes a fine-grained grid of Θ. We note that the estimated MOE weight vector βbn , as well as1 the oracle β∗ , depends on θ. Therefore, β∗ (θ) must be re-estimated as θ ranges over Θgrid . The detailed steps are summarized as in Algorithm 1. If a point estimate of θ∗ is desired, we propose θbMOE := arg min

θ∈Θgrid

b b , eθ,βbn + ∆ g θ,βn

although we do not pursue the theoretical analysis of θbMOE here. Algorithm 1 MOE-powered Inference for M-Estimation ei }[N ] , MOE family FB := {Fβ , β ∈ B}, 1: Input: labeled and unlabelled dataset {(Xi , Yi )}[n] , {X the gradient of loss function gθ (·, ·), initial set CαMOE = ∅; 2: for θ ∈ Θgrid do 3: 4: 5: 6: 7:

solve βbn = arg minβ∈B tr SCov gθ (X, Y ) − gθ (X, Fβ (X)) ; P e e eθ,βbn ← N −1 N imputed gradient: g i=1 gθ (Xi , Fβbn (Xi )); b b ← n−1 Pn g (X , Y ) − g (X , F (X )) ; rectifier: ∆ i i i i θ θ i=1 βbn θ,βn c θ,Y −F ← SCov gθ (X, Y ) − gθ (X, F b (X)) ; sample covariance: W bn β β n q ⊤W b c e if e⊤ g + ∆ ≤ z e e /n for ∀s ∈ [p] then θ,Y −Fβb s α/(2p) s s θ,βbn θ,βbn n

8: 9:

update: CαMOE ← CαMOE ∪ {θ}; end if

10: end for MOE 11: Output: confidence set Cα

For a sequence of random variables {Xn } and a deterministic sequence {αn }, we write Xn = ep (αn ) if P |Xn | ≥ Cαn = O(n−2 ) for some constant C > 0. Recall that the classical notation O 1

As noted above, we suppress the dependence of βbn on θ for notational simplicity.

Xn = Op (αn ) means that, ∀ε > 0, there exists an nε > 0 such that P |Xn | ≥ Cαn ≤ ε for all n ≥ nε . Theorem 1 (M-estimation). Suppose Assumption 1 holds and n ≥ C1 p2 q 2 + log2 (Dn) and N ≥ C1 pn log n for a large constant C1 > 0. Then, for any θ ∈ Θ, r √ p(q + log(Dn)) pn log n b e √ eθ,βbn + ∆θ,βbn − Egθ (X, Y ) = Zn,θ,β∗ + Op , (6) n g + N n P where Zn,θ,β∗ = n−1/2 ni=1 gθ (Xi , Yi ) − gθ (Xi , Fβ∗ (Xi )) − Egθ (X, Y ) − Egθ (X, Fβ∗ (X)) . Given θ∗ ∈ Θgrid , if n ≥ C2 p4 q 2 + log2 (Dn) and N ≥ C2 p3 n log n for a large constant C2 > 0, then the coverage probability satisfies r 2 p (q + log(Dn)) p3 n log(n) MOE √ P(θ∗ ∈ Cα ) ≥ 1 − α + O . + N n

(7)

We remark that the global uniqueness of θ∗ is not required. The estimating equation E[gθ (X, Y )] = 0 may have multiple solutions, the asymptotic normality and coverage statement in Theorem 1 remain valid for all regular solutions. Moreover, the complexity of FB , as measured by its dimension q, should be kept under control. In particular, q should not much larger than log n; otherwise, the resulting procedure may suffer from increased estimation error and reduced statistical efficiency. By taking into consideration of the variance contributed by the unlabelled data, we can construct p a slightly larger confidence set, which improves the remainder term O( (pn/N ) log n) in (6) and allows a wider range of the size of unlabelled data. Theorem 2 (M-estimation+). Suppose Assumption 1 holds and n ≥ C1 p2 q 2 + log2 (Dn) and N ≥ C1 p log n for a large C1 > 0. Then, for any θ ∈ Θ, r √ p q + log(Dn) n b b − Egθ (X, Y ) = Zn,θ,β + e N,θ,β + O ep √ eθ,βbn + ∆ n g Z . (8) ∗ ∗ θ,βn N n h i e N,θ,β = N −1/2 PN gθ (X c θ,F = SCov gθ (X, F b (X)) ei , Fβ (X ei ))−Egθ (Xi , Fβ (Xi )) . Let W where Z ∗ ∗ ∗ i=1 b βn β n

and construct the confidence set as b b eθ,βbn + ∆ e⊤ s g θ,βn MOE grid e Cα = θ∗ ∈ Θ : q ≤ zα/(2p) , ∀s ∈ [p] . −1 c −1 W c θ,F es e⊤ s n Wθ,Y −Fβbn + N b βn Given θ∗ ∈ Θgrid , if n ≥ C2 p4 q 2 + log2 (Dn) and N ≥ C2 p2 log n for a large C2 > 0, then the coverage probability satisfies r 2 p q + log(Dn) p2 log(n) MOE e √ + P(θ∗ ∈ Cα ) ≥ 1 − α + O . N n 10

Remark 1. The oracle weight vector β∗ and its estimator βbn minimize the trace of WY −Fβ and c Y −F , respectively. From the perspective of Theorem 2, a more refined choice is to minimize the W β trace of n−1 WY −Fβ + N −1 WFβ

and

cF , c Y −F + N −1 W n−1 W β β

including the covariance contributed by the unlabelled data. The corresponding oracle and empirical minimizers are denoted by β∗+ and βbn+ , respectively. Theorem 2 continues to hold if these refined empirical estimates are employed.

Applications

For ease of exposition, we focus on the linear mixture of experts throughout this section. More specifically, given K experts f1 (·), · · · , fK (·), the candidate MOE resides in F :=

Fβ :=

k X

βk fk , β ∈ R

k=1

Let L := {(X, Y ), (X1 , Y1 ), · · · , (Xn , Yn )} ⊂ X × Y be a collection of labelled i.i.d. observations and e X e1 , · · · , X eN } be a collection of unlabelled i.i.d. observations, which is independent from U := {X, e and X are identically distributed and X ⊂ Rd . L. Moreover, all covariates X

3.1

Mean Value Inference

Mean value estimation is perhaps the simplest yet most classical application of semi-supervised learning (Angelopoulos et al., 2023a). The goal is to estimate the mean response θ∗ := EY , by leveraging the unlabelled observations and multiple experts f1 , · · · , fK : X 7→ Y which provide “reliable” prediction of the response. Without loss of generality, we focus on the case of univariate response where Y ⊂ R and is compact. For any X, denote f (X) = (f1 (X), · · · , fK (X))⊤ , and we write f for brevity. Assumption 2. There exist constants c0 , C0 > 0 such that |Y | ≤ C0 and |fk (X)| ≤ C0 almost surely for all k = 1, · · · , K, and Cov(f ) ≽ c0 · IK . The uniformly boundedness condition in Assumption 2 is imposed for technical simplicity, which is relaxable. We assume that the K experts are not overly correlated with each other, otherwise, the optimal mixture may not be unique and our MOE procedure becomes unstable. 11

The conventional estimator θbcon and the PPI-based (Angelopoulos et al., 2023a) estimator θbPPI , given an expert f (·), are defined by n X

1 θbcon : = n 1 θbfPPI : = N

Yi ,

i=1 N X i=1

(9)

ei ) − 1 f (X n

n X

(f (Xi ) − Yi ) ,

(10)

i=1

respectively. These two estimators are both unbiased, and as shown in Angelopoulos et al. (2023a), their covariances are Var(θbcon ) ≍

1 · Var(Y ) n

and

Var(θbPPI ) ≍

1 1 · Var Y − f (X) + · Var f (X) . n N

The PPI-based estimator achieves a smaller variance than the conventional one as long as the expert is accurate in that Var Y − f (X) ≪ Var(Y ) and the unlabelled data is abundant in that N ≫ n. Given that multiple experts f1 (·), · · · , fK (·) are available, our MOE-powered approach aims to exploit a mixture of these experts to achieve a robust PPI-based estimator. More specifically, given β ∈ RK , define the MOE-powered estimator by 1 θbβMOE = N

N X

1X (Fβ (Xi ) − Yi ) , n

ei ) − Fβ (X

i=1

(11)

i=1

where the linear mixture of experts Fβ (x) :=

k βk fk (x).

For any fixed β, the MOE-powered

= θbfPPI where ek denotes the k-th canonical basis vector. estimator θbβMOE is unbiased. Moreover, θbeMOE k k To this end, the MOE-powered approach targets the linear mixture which achieves the smallest variance and define the oracle MOE weight vector by β∗ := arg minβ Var Fβ (X) − Y . Its sample version based on labelled data is given by n 2 1X βbn := arg min SVar Fβ (X) − Y := Fβ (Xi ) − Yi − n β i=1

X n 2 1 Fβ (Xi ) − Yi , n

(12)

i=1

which admits a closed-form solution. ei , define the K-dimensional predictor vectors by fi = f1 (Xi ), · · · , fK (Xi ) ⊤ For each Xi and X ei ), · · · , fK (X ei ) ⊤ , respectively. We assume, without loss of generality, that E(f1 f ⊤ ) and e fi = f1 (X 1 is invertible. Collectively, we define the predictor matrices by F = (f1 , · · · , fn )⊤ ∈ Rn×K and e = (e F f1 , · · · , e fN )⊤ ∈ RN ×K . As a result, the sample MOE weight vector βbn from (12) actually 12

minimizes (Fβ − y)⊤ Pn (Fβ − y), where y = (Y1 , · · · , Yn )⊤ and Pn = In − n−1 1n 1⊤ n . Here, 1n denotes the n-dimensional all-one vector. The solution to (12) is βbn = (F⊤ Pn F)−1 F⊤ Pn y. −1 Similarly, the oracle MOE weight vector is β∗ = Cov(f ) Cov(f , Y ) and , = θbfPPI θbβMOE ∗ ∗

where f∗ (X) = ⟨β∗ , f ⟩ =

K X

β∗,k fk (X).

(13)

k=1

The function f∗ is the oracle mixture of experts. Finally, our MOE-powered mean value estimator can be written as −1 ⊤ e b b θbβMOE 1N Fβn − n−1 1⊤ n Fβn − y bn =N ⊤ −1 ⊤ −1 ⊤ −1 ⊤ e = N −1 1⊤ N F − n 1n F (F Pn F) F Pn y + n 1n y.

(14)

Hereafter, we shall write θbMOE for brevity whenever its dependence on βbn is clear from the context. Unlike the conventional and PPI-based estimator, the MOE-powered estimator can be biased. The following lemma provides an upper bound of the bias of θbMOE . Lemma 2. Suppose that Assumption 2 holds and n ≥ C1 K 2 log n for a large absolute constant C1 > 0. The MOE-powered estimator by (14) is biased and 2 K log n MOE b Eθ − θ∗ = O n Lemma 2 shows that the bias of θbMOE is negligible. The MOE-powered inference procedure b βn

is detailed in Algorithm 2. For any α ∈ (0, 1), we denote zα := Φ−1 (1 − α) where Φ(·) is the cumulative distribution function of a standard normal random variable. Algorithm 2 MOE-powered Inference on Mean Value 1: Input:

e = (e labelled and unlabelled predictor matrices F = (f1 , · · · , fn )⊤ and F f1 , · · · , e fN )⊤ ,

the observed response vector y = (Y1 , · · · , Yn )⊤ . ⊤ −1 ⊤ 2: Mean response and predictors: Ȳ ← n−1 1⊤ n y and f̄ ← n 1n F.

bn ← (n−1 F⊤ F − f̄ f̄ ⊤ )−1 (n−1 F⊤ y − f̄ Ȳ ) 3: Estimate MOE weight vector: β −1 ⊤ bMOE ← N −1 1⊤ F b eb 4: MOE-powered point estimate: θ N βn − n 1n Fβn − y

2 b n−1 1⊤ n (Fβn − y) bMOE and its 100(1 − α)% confidence interval Cα := θbMOE ± zα/2 n−1/2 σ bY −f∗ 6: Output: θ 5: Estimated the variance: σ bY2 −f∗ ← n−1 ∥Fβbn − y∥2 −

Theorem 3 (mean value inference). Suppose that Assumption 2 holds, n ≫ K 4 log2 n, and N ≫ n log n. The MOE-powered estimator θbMOE produced by Algorithm 2 satisfies r 2 √ K log n n log n MOE b e √ , (15) n θ − θ∗ = −Zn,f∗ + Op + N n P where n1/2 Zn,f∗ := ni=1 f∗ (Xi ) − Yi − E f (X) − Y . Moreover, the coverage probability of Cα output by Algorithm 2 satisfies P θ∗ ∈ Cα

K 2 log n √ =1−α+O + n

n log n . N

(16)

Remark 2. The confidence interval constructed by Algorithm 2 has a width characterized by σY2 −f∗ = Var Y − f∗ (X) . It is determined by the accuracy of the optimal mixture of experts. Since Var Y − f∗ (X) ≤ mink Var Y − fk (X) , our MOE-powered confidence interval is at least as short as the PPI-based confidence interval by any individual expert. This property holds for all the other algorithms in this section. The confidence interval constructed in Algorithm 2 ignores the variance contributed by the unlabelled data, which also introduces the remainder term (n/N ) log n in (15) and (16). While this remainder term is negligible for large N , we can still eliminate it by slightly modifying Algorithm 2 and using a slightly wider confidence interval. Theorem 4 (mean value inference+). Suppose that Assumption 2 holds and n ≫ K 4 log2 n. The MOE-powered estimator θbMOE produced by Algorithm 2 satisfies 2 √ K log n MOE b e e √ n θ − θ∗ = −Zn,f∗ + ZN + Op , n P eN := n1/2 N −1 PN f∗ (X ei )−Ef∗ (X) . where n1/2 Zn,f∗ := ni=1 f∗ (Xi )−Yi −E f (X)−Y and Z i=1 2 e βbn ∥2 − (N −1 1⊤ F eb 2 Let σ bf2∗ := N −1 ∥F N βn ) be an estimate of σf∗ = Var f∗ (X) and construct s s ! 2 2 2 2 σ b σ b σ b σ b Y −f f Y −f f ∗ ∗ Ceα := θbMOE − zα/2 + ∗ , θbMOE + zα/2 + ∗ . n N n N Then, P θ∗ ∈ Ceα = 1 − α + O K 2 n−1/2 log n + N −1/2 log n .

3.2

Quantile inference

Quantiles play a central role in statistical inference because they provide a complete description of the underlying distribution. For any fixed q ∈ (0, 1), we define the q-quantile of Y as θ∗ := 14

inf θ ∈ R : P(Y ≤ θ) ≥ q , where we suppress the dependence of θ∗ on q for brevity. Similarly as in Section 2, we search for the candidate solutions from a fine-grained grid Θgrid ⊂ Θ. Under mild continuity assumption on the distribution function of Y , we have E I{Y ≤ θ∗ } = q. Motivated by this equation (see Angelopoulos et al. (2023a) for more details), the conventional and PPI-based approaches construct the confidence sets by studying the following estimating functions: m b con (θ) := n−1

n X I Yi ≤ θ − q, i=1

−1 m b PPI f (θ) := n

n X

N X ei ) ≤ θ − q, I Yi ≤ θ − I f (Xi ) ≤ θ + N −1 I f (X

i=1

where f (·) is a given predictor. Both of them are unbiased in that Em b con (θ) = Em b PPI f (θ) = E I(Y ≤ θ) − q. Moreover, if θ = θ∗ , by central limit theorem, these two estimating functions can be normally approximated by Var I Y ≤ θ∗ and N 0, n

e ≤ θ∗ Var I Y ≤ θ∗ − I f (X) ≤ θ∗ Var I f (X) N 0, + , n N respectively. If N ≫ n and f (·) is a powerful predictor so that Var I{Y ≤ θ} − I{f (X) ≤ θ} ≤ b con (θ) and the PPI-based approach may yield a tighter Var I{Y ≤ θ} , then Var m b PPI f (θ) ≤ Var m confidence set. Note that Var I Y ≤ θ∗ = q − q 2 and 2 Var I Y ≤ θ∗ − I f (X) ≤ θ∗ = P f (X) − θ∗ Y − θ∗ < 0 − q − P(f (X) ≤ θ∗ ) . By the shape of the function t−t2 for t ∈ [0, 1], a necessary condition for Var I{Y ≤ θ∗ }−I{f (X) ≤ θ∗ } ≤ Var I{Y ≤ θ∗ } is that q − P(f (X) ≤ θ∗ ) ≤ min(q, 1 − q). However, due to the non-smoothness of indicator function, the general M-estimation framework in Section 2 is not immediately applicable. There are several common approaches to smooth the indicator function, such as linearization at the jump point (Koltchinskii and Xia, 2016). Here, we adopt the sigmoid function Sh (t) = (1 + exp{−t/h})−1 , ∀t ∈ R for its higher-order smoothness, where h > 0 is a bandwidth parameter. Note that Sh′ (t) ≤ h−1 . The bandwidth controls the tradeoff between bias and variance. Based on the smoothed indicator function, we define the conventional, PPI-based, and MOE-

powered estimating functions by n

1X m b (θ) := Sh θ − Yi − q, n con

m b PPI f (θ) :=

1 n

i=1 n X i=1

N 1 X ei ) − q, Sh θ − Yi − Sh θ − f (Xi ) + Sh θ − f (X N i=1

n N 1 X 1 X MOE ⊤ Sh θ − Y i − S h θ − f i β + Sh θ − e m b β (θ) := fi⊤ β − q, n N i=1

(17)

i=1

where we suppress the dependence of hn on n for notational brevity and fi , e fi are as defined in Section 3.1. Hereafter, we take h := n−γ ,

γ∈

1 1 , . 4 2

(18)

The oracle MOE and its sample counterpart are defined as the minimizers of the population and sample variances, respectively. Denote Q(θ, β) := Var Sh (θ − Y ) − Sh (θ − f ⊤ β) and Qn (θ, β) := SVar Sh (θ − Y ) − Sh (θ − f ⊤ β) . More specifically, we have n n 2 1 X 2 1 X ⊤ ⊤ Sh (θ − Yi ) − Sh (θ − fi β) − Sh (θ − Yi ) − Sh (θ − fi β) Qn (θ, β) : = . n n i=1

i=1

Finally, the oracle MOE weight vector and its sample counterpart are defined by β∗ (θ) := arg min Q(θ, β)

and

βbn (θ) := arg min Qn (θ, β),

β∈B

(19)

β∈B

where B ⊂ RK is a compact set. Hereafter, we will suppress the dependence of β∗ and βbn on θ for simplicity. Assumption 3. There exists a constant U > 0 such that ∥f (X)∥ ≤ U almost surely. Moreover, let β∗ = β∗ (θ∗ ) and (A1) there exist δ > 0 and C0 > 0 such that supβ∈B∗ (δ) pβ (·) ≤ C0 , where B∗ (δ) := {β ∈ B : ∥β − β∗ ∥ ≤ δ} and pβ is the probability density function of Fβ (X); (A2) there exist δ0 > 0, c0 > 0 and η0 > 0 such that inf β∈B∗ (δ0 )

λmin h∇2 Q(θ∗ , β) ≥ c0

and

inf

β ∈B / ∗ (δ0 )

Q(θ∗ , β) − Q(θ∗ , β∗ ) ≥ η0 .

(A3) there exists a neighborhood N0 of θ∗ such that the probability density function fY is differentiable on N0 , and fY′ is continuous at θ∗ . 16

Note that β∗ (θ) changes with respect to θ, but the conditions (A1) and (A2) in Assumption 3 are imposed only on β∗ (θ∗ ). The condition (A1) characterizes the smoothness of Q(θ∗ , β) with = O(∥β1 − β2 ∥) respect to β ∈ B∗ (δ). It implies that E Sh (θ∗ − f ⊤ β1 ) − E Sh (θ∗ − f ⊤ β2 ) for all β1 , β2 ∈ B∗ (δ). The condition (A2) imposes locally strong convexity and global separation condition around β∗ , which ensure the uniqueness of β∗ . This condition is weaker than the global strongly convexity in Assumption 1. And the condition (A3) ensures the feasibility of smoothing approach. Lemma 3. Suppose Assumption 3 holds and n ≥ C1 (K + log n) for a large constant C1 > 0. Then, there exists a constant C2 > 0 and an event E with P(E) ≥ 1 − n−10 on which r K + log n βbn − β∗ ≤ C2 , n

(20)

where βbn = βbn (θ∗ ) and β∗ = β∗ (θ∗ ) are as defined in (19). For any fixed β, the MOE-powered estimating function m b MOE β (θ∗ ) is biased due to the kernel 2 smoothing. Indeed, Em b MOE β (θ∗ ) = ESh (θ∗ − Y ) − q = O(h ), provided that Y has a smooth density

in a neighborhood of θ∗ . The following lemma shows that a similar bias remains even if the sample MOE weight βbn is employed. Lemma 4. Suppose that Assumption 3 holds and n ≥ C1 (K + log n) for a large constant C1 > 0. The MOE-powered estimating function by (17) satisfies K + log n MOE 2 . Em b βb (θ∗ ) = O h + n nh −1 Since the variance of m b MOE b (θ∗ ) is at the order of O(n ), Lemma 4 implies that a valid inference βn

of θ∗ requires 1/4 < γ < 1/2 if h ≍ n−γ . The MOE-powered inference procedure is detailed in Algorithm 3. The coverage probability of the MOE-powered confidence set is guaranteed by Theorem 5.

Algorithm 3 MOE-powered Inference on Quantile 1: Input: labelled and unlabelled predictor vectors fi , 1 ≤ i ≤ n and e fi , 1 ≤ i ≤ N , the observed

responses Yi , 1 ≤ i ≤ n, quantile level q ∈ (0, 1), bandwidth h, sigmoid function Sh (t) = 1/(1 + e−t/h ), initial set CαMOE = ∅; 2: for θ ∈ Θgrid do 3:

solve βbn = arg minβ Qn (θ, β);

P e⊤ b imputed gradient: geθ,βbn ← N −1 N i=1 Sh (θ − fi βn ); b b ← n−1 Pn Sh (θ − Yi ) − Sh (θ − f ⊤ βbn ) ; rectifier: ∆ i i=1 θ,βn 2 b n (θ, βbn ) ← n−1 Pn Sh (θ − Yi ) − Sh (θ − f ⊤ βbn ) 2 − ∆ b b ; sample variance: Q i i=1 θ,βn q b b − q ≤ zα/2 Q b n (θ, βbn )/n then if geθ,βbn + ∆ θ,βn

5: 6: 7:

update: CαMOE ← CαMOE ∪ {θ};

8: 9:

end if

10: end for MOE . 11: Output: confidence set Cα

Theorem 5 (quantile inference). Suppose that Assumption 3 and the conditions of Lemma 4 hold. The MOE-powered estimating function satisfies √

MOE

nm b βb

ep (θ∗ ) = Zn,f∗ + O

√

K + log n nh + √ + nh 2

K + log n + Nh

n log n , N

(21)

where n

1/2

n 1 X Sh (θ∗ − Yi ) − Sh (θ∗ − fi⊤ β∗ ) − E Sh (θ∗ − Yi ) − S(θ∗ − fi⊤ β∗ ) . Zn,f∗ := √ n i=1

Moreover, the coverage probability of CαMOE output by Algorithm 3 satisfies r r √ 2 K + log n K + log n n log n MOE √ P θ∗ ∈ Cα =1−α+O nh + + + . Nh N nh Remark 3. The optimality choice of the bandwidth in (18) is hn = n−1/3 , which implies r √ n log n MOE −1/6 e nm b βb (θ∗ ) = Zn,f∗ + Op n (K + log n) + , n N r n log n MOE −1/6 . P θ∗ ∈ Cα =1−α+O n (K + log n) + N

(22)

The confidence interval constructed in Algorithm 3 ignores the variance contributed by the unlabelled data, which also introduces the remainder term (n/N ) log n in (21) and (22). While this remainder term is negligible for large N , we can still eliminate it by slightly modifying Algorithm 3 and using a slightly wider confidence interval. Theorem 6 (quantile inference+). Suppose that Assumption 3 and conditions of Lemma 4 hold. The MOE-powered estimating function produced by Algorithm 3 satisfies √ 2 K + log n √ MOE e e , nm b βb (θ∗ ) = Zn,f∗ + ZN,f∗ + Op nh + √ n nh

(23)

where n

1/2

Zn,f∗ :=

n X

Sh (θ∗ − Yi ) − Sh (θ∗ − fi⊤ β∗ ) − E Sh (θ∗ − Y ) − S(θ∗ − f ⊤ β∗ ) .

i=1

eN,f := n1/2 N −1 n1/2 Z ∗

N X

Sh (θ∗ − e fi⊤ β∗ ) − ESh (θ∗ − f ⊤ β∗ ) .

i=1

⊤β ⊤β bn ) 2 be an estimate of Q(θ, e bn )2 − N −1 PN Sh (θ − e e β∗ ) = e n (θ, βbn ) := N f S (θ− f Let Q h i i i=1 i=1 Var Sh (θ − f ⊤ β∗ ) , then give the refined confidence set s e n (θ, βbn ) Qn (θ, βbn ) Q MOE MOE e Cα := θ ∈ Θgrid : m b βb (θ) ≤ zα/2 . + n n N PN −1

Then the coverage probability holds that √ 2 K + log n 1 P θ∗ ∈ CeαMOE = 1 − α + O nh + √ +√ nh N

3.3

! .

Linear regression and its inference

Given the pair of covariates and response (X, Y ) ∈ X × Y ⊂ Rd × R, the best linear model and its coefficients are defined by θ∗ = arg min E(Y − X ⊤ θ)2 = (EXX ⊤ )−1 E(XY ).

(24)

For simplicity, we assume θ∗ is unique (see Assumption 4). Motivated similarly as in Section 3.1, we consider the conventional, PPI, and our MOE-powered estimators defined as −1 θbcon := X⊤ X X⊤ y e ⊤X e −1 X e ⊤ f (X) e − X⊤ X −1 X⊤ (f (X) − y) θbfPPI := X e ⊤X e −1 X e ⊤F e βbn − X⊤ X −1 X⊤ (F⊤ βbn − y), θbβMOE := X b n

(25)

e = X e1 , · · · , X eN )⊤ ∈ respectively, where the feature matrices X = (X1 , · · · , Xn )⊤ ∈ Rn×d and X RN ×d are constructed using labelled and unlabelled data, respectively. The matrices F ∈ Rn×K and e ∈ RN ×K are constructed by applying experts to the feature matrices, similarly as in Section 3.1. F Here, f (X) ∈ Rn denotes the operation of applying the expert f (·) to the rows of X. The MOEpowered estimator is reliant on βbn , i.e., the estimation of the oracle mixture of experts., which we shall unfold shortly. The asymptotic properties of θbcon and θbfPPI are well-documented in the literature. It is clear that both estimators are biased, yet consistent and asymptotically normal. It was shown by White (1980) that √

d. n(θbcon − θ∗ ) −→ N 0, WY −X ⊤ θ∗ , as n → ∞, where WY −X ⊤ θ∗ := Σ−1 E (Y −X ⊤ θ∗ )2 XX ⊤ Σ−1 and Σ := EXX ⊤ . Note that we used the simple fact E(Y − X ⊤ θ∗ )X = 0 by the definition in (24). Similarly, Angelopoulos et al. (2023a) shows that √

d. n(θbfPPI − θ∗ ) −→ N 0, WY −f ,

as n, N → ∞ and n/N → 0,

2 where WY −f := Σ−1 E Y − f (X) − X ⊤ δ f XX ⊤ Σ−1 with the rectifier defined by δ f := Σ−1 E X Y − f (X) . Note that WY −f ̸= Σ−1 Cov X(Y − f (X)) Σ−1 in general. The recti fier serves as a debiasing treatment in the sense that E X(Y − f (X) − δ ⊤ f X) = 0. The PPI-based estimator enjoys a smaller variance than the conventional one when the expert f (·) is accurate, e.g., if (Y − f (X) − X ⊤ δ f )2 ≪ (Y − X ⊤ θ∗ )2 almost surely. Note that a finite-sample estimate of Cov(θbfPPI ) is given by c Y −f := n(X⊤ X)−1 · W

X n

bf Yi − f (Xi ) − Xi⊤ δ

Xi Xi⊤

· (X⊤ X)−1 ,

(26)

i=1

bf := (X⊤ X)−1 · where δ

i=1 Xi (Yi − f (Xi )) = (X

⊤ X)−1 X⊤ y − f (X) .

Our MOE-powered estimator (25) intends to minimize the (co)variance of θbMOE through a careb βn

fully chosen βbn . As discussed in Section 3.1, given a fixed mixture assignment β, the MOE-powered P estimator can be written as a PPI-based estimator via θbβMOE = θbFPPIβ , where Fβ := K k=1 βk fk . Moreover, Cov(θbβMOE ) is given by 2 WY −Fβ = Σ−1 · E Y − Fβ (X) − X ⊤ δ Fβ XX ⊤ · Σ−1 .

(27)

Since θbβMOE is multivariate, its variability can be viewed from several perspectives. For instance, for a fixed index s ∈ [d], one may consider the entrywise variance Var (θbβMOE )s = e⊤ s WY −Fβ es . Here, 20

es denotes the s-th canonical basis vector, whose dimension may vary at different appearances. Without loss of generality and for ease of exposition, we focus primarily on the total variance of θbMOE , defined as tr WY −Fβ . Define the oracle MOE weight vector by β∗ := arg minβ tr(WY −Fβ ), which aims to minimize the total variance of our MOE-powered estimator. Accordingly, define the oracle MOE by f∗ (X) = Fβ∗ (X). Its sample version is simply defined by c Y −F , βbn :=arg min tr W β

(28)

bF )X · (X⊤ X)−1 is as defined in (26). c Y −F = n(X⊤ X)−1 · X⊤ diag2 (y − Fβ (X) − Xδ where W β β Here, diag(v) represents a diagonal matrix constructed from the vector v. Note that Fβ (X) = Fβ with the matrix F ∈ Rn×K being defined as in Section 3.1. Hereafter, we shall write θbMOE for brevity whenever its dependence on βbn is clear from the context. Define Σ := EXX ⊤ and H := E (X ⊤ Σ−2 X)v(X)v⊤ (X) where v⊤ (X) = X ⊤ Σ−1 E(Xf ⊤ ) − f ⊤ , with f := f (X) = f1 (X), · · · , fK (X)

⊤

Assumption 4. There exist constants c0 , C0 > 0 such that |Y | ≤ C0 , ∥X∥ ≤ C0 , and |fk (X)| ≤ C0 almost surely for all k = 1, · · · , K, Σ ≽ c0 · Id , and H ≽ c0 · IK . The matrix Σ is invertible, so θ∗ is uniquely defined, whereas H is invertible to ensure the uniqueness of β∗ . Lemma 5. Suppose that Assumption 4 holds and n ≥ C1 K 2 d log n for a large absolute constant C1 > 0. The bias of the MOE-powered estimator by (25) is bounded by 2√ (K d + d) log n MOE Eθb − θ∗ = O . n Lemma 5 shows that the bias of MOE-powered estimator is negligible if sample size is large. The inference procedure is detailed in Algorithm 4.

Algorithm 4 MOE-powered Inference on Linear Regression 1:

e = (e Input: labelled and unlabelled predictor matrices F = (f1 , · · · , fn )⊤ and F f1 , · · · , e fN )⊤ , e = (X e1 , · · · , X eN )⊤ , the labelled and unlabelled covariate matrices X = (X1 , · · · , Xn )⊤ and X observed response vector y = (Y1 , · · · , Yn )⊤ , an index s ∈ [d].

b ← (X⊤ X)−1 X⊤ y; b ← n−1 X⊤ X, A b ← (X⊤ X)−1 X⊤ F, b 2: Σ b ← n−1 3: H 4: b r ← n−1

⊤ b −2 ⊤ b⊤ b⊤ i=1 (Xi Σ Xi )(fi − A Xi )(fi − A Xi ) ;

⊤ b −2 ⊤b b⊤ i=1 (Xi Σ Xi )(fi − A Xi )(Yi − Xi b);

b ≤ 0.01c0 then 5: if λmin (H) 6:

βbn ← 0;

7: else 8:

b −1b estimate MOE weight vector: βbn ← H r;

9: end if

e ⊤X e −1 X e ⊤F e βbn − X⊤ X −1 X⊤ (F⊤ βbn − y); X −1 bF = (X⊤ X)−1 X⊤ (y − bF )X · Σ b with δ c = n−1 Σ b −1 · X⊤ diag2 (y − Fβbn − Xδ 11: covariance W b b β β

bMOE ← 10: MOE point estimate θ

Fβbn ). MOE point estimate θbMOE and entrywise confidence interval Cα,s = q c s /n . We zα/2 e⊤ s

12: Return

θbsMOE ±

The non-asymptotic performance of the MOE-powered estimator and the coverage guarantee for its confidence intervals are given in Theorem 7. For simplicity, we focus on entrywise inference for fixed entries. Extending the results to simultaneous inference for all entries is possible, but it would substantially complicate the presentation of our method; therefore, we omit it. Theorem 7 (linear regression). Suppose Assumption 4 holds, n ≥ C0 dK 3 log2 n, and N ≥ C0 dn log n for a large constant C0 > 0. Let θbMOE be the MOE-powered estimator output by Algorithm 4. Then, for any fixed index s ∈ [d], √ 3 log n dn log n dK ep √ n θbs − θ∗,s = Zn,f∗ ,s + O + , (29) N n Pn ⊤ −1 · n−1/2 where Zn,f∗ ,s = e⊤ sΣ i=1 Yi − f∗ (Xi ) − Xi δ f∗ Xi . Moreover, the coverage probability √

MOE

of Cα,s output by Algorithm 4 satisfies

P θ∗,s ∈ Cα,s = 1 − α + O 22

dn log n + N

√ dK 3 log n √ . n

(30)

3.4

Logistic regression and its inference

Given the pair of covariates and response (X, Y ) ∈ X × Y ⊂ Rd × {0, 1} for binary classification, and target of inference for logistic regression is defined by θ∗ = arg min E − Y θ⊤ X + log 1 + exp(θ⊤ X) , θ∈Θ

where Θ ⊂ Rd . The population estimating function becomes m(θ) := E[ Y − S(θ⊤ X) X] with the sigmoid function S(t) = (1 + e−t )−1 . Therefore, θ∗ solves the equation m(θ) = 0. Let Θgrid be a fine-grained grid of Θ. e y, f , e e be constructed similarly as in Section 3.3. The conventional, PPILet X, X, f , F, and F based, and MOE-powered estimation functions are defined by b con (θ) := n−1 X⊤ y − S(Xθ) , m −1 e ⊤ e − S(Xθ) e b PPI X f (X) − n−1 X⊤ f (X) − y , m f (θ) := N −1 e ⊤ e e b MOE m X Fβ − S(Xθ) − n−1 X⊤ Fβ − y , β (θ) := N

(31)

respectively, where S(Xθ) = (S(X1⊤ θ), · · · , S(Xn⊤ θ))⊤ ∈ Rn . b con and m b PPI Note that both m f are both unbiased estimators of the population estimating func tion. Denote WY −X ⊤ θ := Cov X(Y − S(X ⊤ θ)) and WY −f := Cov X(Y − f (X)) . Then, under mild conditions, CLT dictates that d. b con (θ) − m(θ) → N 0, WY −X ⊤ θ , as n → ∞; n1/2 m d. b PPI n1/2 m f (θ) − m(θ) → N 0, WY −f , if n/N → 0 as n → ∞. For any given θ, the oracle MOE weight vector is defined as the minimizer of the total variance of b MOE m β (θ). Its empirical counterpart is defined in a similar fashion. More exactly, let β∗ := arg min tr WY −Fβ β

c Y −F , and βbn := arg min tr W β

(32)

P where the β-index mixture of experts Fβ (X) = f ⊤ β = k βk fk (X). The sample covariance matrix c Y −F = SCov X(f ⊤ β − Y )) with an explicit form as W β c Y −F = X⊤ diag(y − Fβ)Pn diag(y − Fβ)X, n·W β c where Pn = In − n−1 1n 1⊤ n . Note that these covariance matrices WY −f and WY −f are irrelevant to θ. 23

Assumption 5. There exist constant c0 , C0 such that ∥X∥ ≤ C0 , |Y | ≤ C0 and |fk (X)| ≤ C0 for all k ∈ [K] almost surely, and E f X ⊤ − Ef X ⊤ Xf ⊤ − EXf ⊤ ⪰ c0 Id . The last condition of Assumption 5 guarantees that the minimizer β∗ is unique. See the proof of Lemma 6 for the explicit form of β∗ . The following lemma provides the high-probability upper bound for βbn − β∗ . Lemma 6. Suppose that Assumption 5 holds, n ≥ C1 (d + K) log n, and N ≥ C1 n for a large constant C1 > 0. Then, (d + K) log n b e ∥βn − β∗ ∥ = Op n

and

b MOE (θ) − m(θ)) = O Em βb n

dK log n n

for any fixed θ ∈ Θ. The inference procedure is detailed in Algorithm 5. Note that the estimated MOE weight vector βbn is irrelevant of θ, whereas the imputed gradient must be recalculated for each θ ∈ Θgrid . Algorithm 5 outputs only the MOE-powered confidence set. If a point estimate of θ∗ is desired, we b b ∥, but its theoretical property is skipped. suggest θbMOE := arg minθ∈Θgrid ∥e gθ,βbn − ∆ βn Theorem 8 (Logistic regression). Suppose that Assumption 5 holds, n ≥ C1 (d + K)2 log2 n, and N ≥ C1 (d + K)n log n for a large constant C1 > 0. Then, for any fixed θ ∈ Θ, r √ (d + K) log n (d + K)n log n MOE e √ b βb (θ) − m(θ) = −Zn,f∗ + Op , (33) n m + n N n i P h where n1/2 Zn,f∗ := ni=1 Xi Fβ∗ (Xi )−Yi −EX Fβ∗ (X)−Y . Moreover, the coverage probability of CαMOE output by Algorithm 5 satisfies d(d + K) log n MOE √ P θ∗ ∈ Cα ≥1−α+O + n

d2 (d + K)n log n , N

(34)

given θ∗ ∈ Θgrid . The confidence interval constructed in Algorithm 5 ignores the variance contributed by the unlabelled data, which also introduces the remainder term (n/N ) log n in (33) and (34). While this remainder term is negligible for large N , we can still eliminate it by slightly modifying Algorithm 5 and using a slightly wider confidence interval. Denote e N,f := N −1/2 Z ∗

N h X i ei f∗ (X ei ) − S X ei⊤ θ − EX f∗ (X) − S(X ⊤ θ) , X i=1

Algorithm 5 MOE-powered Inference for Logistic Regression 1:

observed response vector y = (Y1 , · · · , Yn )⊤ , initial Set CαMOE = ∅; b ← n−1 Pn fi X ⊤ Xi f ⊤ − n−1 Pn fi X ⊤ n−1 Pn Xi f ⊤ ; 2: H i i i i i=1 i=1 i=1 P P P 3: b r ← n−1 ni=1 fi Xi⊤ Xi Yi − n−1 ni=1 fi Xi⊤ n−1 ni=1 Xi Yi ; b ≤ 0.01c0 then 4: if λmin (H) 5:

βbn ← 0;

6: else 7:

b −1b estimate MOE weight vector: βbn ← H r;

8: end if

fi⊤ βbn − Yi Xi ; b⊤ ; bb ∆ c Y −F ← n−1 Pn (f ⊤ βbn − Yi )2 Xi X ⊤ − ∆ 10: sample covariance: W i i=1 i b b βn β ← n−1 n

bb 9: rectifier: ∆ β

i=1

βn

11: for θ ∈ Θgrid do 12: 13:

P e⊤ b e⊤ e eθ,βbn ← N −1 N imputed gradient: g i=1 fi βn − S(Xi θ) Xi ; q b b ) ≤ zα/(2d) e⊤ W c Y −F es /n, ∀s ∈ [d] then gθ,βbn − ∆ if e⊤ s (e s b βn β n

14:

Cα

15:

end if

MOE

← Cα

MOE

∪ {θ};

16: end for MOE 17: Return confidence set Cα .

e N,f is asymptotically normal with covariance matrix where f∗ (·) = Fβ∗ (·). It is clear that Z ∗ ⊤ f f −X ⊤ θ = SCov X e f∗ (X) e − Wf∗ −X ⊤ θ = Cov X f∗ (X) − S(X θ) . Its empirical version W ∗ e ⊤ θ) is constructed solely from the unlabelled data. S(X Theorem 9 (Logistic regression+). Suppose that Assumption 5 holds, n ≥ C1 (d + K)2 log2 n, and b MOE N ≥ C1 n log2 n for a large constant C1 > 0. For any θ, the estimating function m b (θ) produced βn

by Algorithm 5 satisfies r √ MOE ne K) log n ep (d + √ b βb (θ) − m(θ)) = −Zn,f∗ + n m ZN,f∗ + O . n N n Construct the confidence set by √ ⊤ MOE ) ( b b (θ) n es m β n MOE ≤ zα/(2d) , ∀s ∈ [d] . Ceα := θ ∈ Θgrid : q f f −X ⊤ θ es /N c Y −F es /n + e⊤ W e⊤ W s s b ∗ β n

If θ∗ ∈ Θgrid , then d(d + K) log n d(d + K) log n MOE e √ √ P θ∗ ∈ C α ≥1−α+O + . n N

Numerical Experiments

This section empirically evaluates the performance of the MOE-powered inference framework through a comprehensive set of experiments. These experiments aim to demonstrate that: (i) MOE outperforms or at least matches PPI equipped with the best single predictor; (ii) MOE is robust under model misspecification; (iii) MOE-powered confidence sets achieve the desired coverage probabilities; and (iv) MOE requires the fewest labeled samples among all methods to achieve the same statistical efficiency.

4.1

Simulation studies

4.1.1

Data modeling

We consider two settings in which the labeled data are generated from a linear model and a nonlinear model, respectively. More specifically, (1) Linear model. The labeled data (X, Y ) satisfy Y = 10 + X ⊤ β + ε, where ε ∼ N (0, σ 2 ). The covariate X ∈ R20 has i.i.d. entries drawn from the standard normal distribution, and 26

β ∈ R20 has 10 nonzero entries in the first 10 coordinates, each sampled from U (0, 100). Here, we set σ = 10. The labeled data set L is drawn from fixed 50,000 samples. In addition, we generate another 50,000 samples to train the experts (predictors), which will be specified in Section 4.1.2. (2) Nonlinear model. The response variable Y contains a nonlinear component of X: Y = X ⊤ β +γg(X)+ε, where g(X) = 10 sin(πX11 X12 )+20(X13 −0.5)2 , and X and β are generated in the same way as above. The nonlinear function g(X) includes sinusoidal, polynomial, and interaction terms and corresponds to the nonlinear component of the Friedman1 model (Friedman, 1991). The parameter γ controls the strength of the nonlinearity. 4.1.2

Candidate expert predictors

We train six classical machine learning models to serve as the expert predictors {f1 , . . . , f6 }: (1) Deep neural network (DNN): 3-layer fully connected network (MLP) with ReLU activations, hidden dimensions [64, 32, 16], and trained with Adam optimizer for 100 epochs; (2) Linear model: the least squares estimator; (3) Quadratic regression: Degree-2 polynomial regression; (4) Random forest: 100 trees with max depth 10, minimum samples split 5; (5) XGBoost: Gradient boosting with 100 estimators, learning rate 0.1, max depth 6; (6) LightGBM: Light gradient boosting with 100 estimators, learning rate 0.1, max depth 6. All experts are trained independently on an independent auxiliary labeled training sample with size 50,000. Table 1 summarises the validation-set RMSE for each trained expert. It highlights Table 1: Validation RMSE of the K = 6 trained expert predictors Experts

Linear

Nonlinear

10.01 / 10.03

26.21 / 16.73

DNN

11.26

18.37

Random Forest

54.08

67.33

XGBoost

30.45

35.87

LightGBM

23.96

30.55

Linear / Polynomial

strikingly different performances across predictors, especially under nonlinear models. All the subsen o P K . quent simulation studies focus on linear mixture of experts F = Fβ (x) = K β f (x) : β ∈ R k k k=1 27

The empirical weight vector βb is learned by the variance-minimization principle detailed in previous sections. For each inferential task, we compare five methods: (1) Conventional: classical estimator based on labeled data; (2) PPI-worst: PPI with the single worst-performing expert; (3): PPI-mean: PPI with the average of K experts; (4) PPI-best: PPI with the single best-performing expert; (5): PPIMOE: our proposed MOE-powered inference method. Note that “best” and “worst” are realized in hindsight, consequently infeasible in practice, evaluated by the width of confidence interval. The results of point estimates are aggregated over 100 Monte Carlo runs. 4.1.3

Inference tasks

We study four inference tasks as summarized in Table 2, together with their population target (estimating function) and their empirical counterparts. For each task, we evaluate the variance of empirical inference target or estimating function by B = 1,000 bootstrap replicates. Table 2: Inference tasks, population targets, and their sample estimates. Inference Task

Population Target

Sample Estimate (MOE-powered)

Mean value

θ∗ = E[Y ]

−1 1⊤ Fβ bn − y eb N −1 1⊤ n N Fβn − n

Quantiles

θ∗ : inf {θ ∈ R : P(Y ≤ θ) ≥ q}

1 Pn i=1 n

Linear regression

−1 θ∗ := E(X ⊤ X) E(XY )

e ⊤F e βbn − X⊤ X −1 X⊤ (F⊤ βbn − y) e ⊤X e −1 X X

Logistic regression

θ∗ : E S(X ⊤ θ) − Y X = 0

e⊤ F e βbn − S(Xθ) e N −1 X − n−1 X⊤ Fβbn − y

4.1.4

1 PN e⊤ b Sh θ − Yi − Sh θ − fi⊤ βbn + N i=1 Sh θ − fi βn − q

Variance reduction

Figure 1 summarizes variance ratios across tasks and model regimes under fixed n = 200 and N = 20,000. Lower values indicate stronger efficiency gains over the conventional estimator. And black error bars denote the bootstrap standard deviation of the variance ratio estimates. Across mean and quantile inference, PPI methods consistently reduce variance, and PPI-MOE closely tracks—and occasionally exceeds—PPI-best. Under the linear (well-specified) setting, PPIbest always chooses Linear Model as its predictor since it is the oracle model. While in nonlinear (misspecified) setting, when Linear Model is no longer the oracle one, PPI-MOE can choose a better ensembled predictor and consequently exceeds PPI-best.

Figure 1: Variance reduction under Linear and Nonlinear regimes across inference tasks. For coefficient inference, performance is more sensitive to predictor quality: in linear settings weak single predictors can inflate variance, whereas under misspecification predictor diversity allows PPI-MOE and PPI-best to recover substantial gains. Overall, PPI-MOE is the most stable method across tasks because it adaptively combines complementary experts instead of relying on a single model.

Figure 2: Effect of the growth rates of n and N on mean-inference variance. 29

b Figure 2 decomposes the PPI-MOE variance into a rectifier variance term Var(Y − f ⊤ β)/n b and an imputed variance term Var(f ⊤ β)/N . As N increases with fixed n, total variance decreases quickly at first and then plateaus at the labeled-data floor. When n and N grow proportionally, both components decrease, confirming that PPI benefits from larger unlabeled samples but remains fundamentally constrained by labeled-sample information. 4.1.5

Coverage and interval width

For coverage probability guarantee, we sample constant number of labeled samples and repeat the procedure of mean inference for 500 times, and focus on whether the true value of mean falls in constructed confidence interval. Comparison of coverage rate and interval width will be shown in Figure 3 to the certainty will be improved when n increases and N/n is fixed to 10.

Figure 3: Coverage and confidence-interval width under linear and nonlinear settings. PPI-MOE delivers the narrowest confidence intervals in both regimes, with strong gains in precision relative to Conventional and robust performance against weak predictors. Coverage is near nominal in linear settings and remains competitive under misspecification, with a mild undercoverage trade-off at moderate n due to adaptive in-sample weight estimation. More extensive experiments and comparisons are summarized in Table 3, which show that PPI-MOE delivers a robust bias–variance trade-off across all inference tasks. 30

Table 3: Width comparison and coverage of PPI-MOE at n = 500.

Task-setting

MOE/Best

MOE/Mean

MOE/Worst

MOE/Conv

Coverage

Code

Mean (L)

0.98

0.82

0.56

0.20

0.940

Mean (N)

0.98

0.81

0.56

0.21

0.940

Quantile (L)

1.01

0.79

0.55

0.32

0.950

***

Quantile (N)

0.98

0.78

0.56

0.35

0.950

***

Linear Reg. (L)

0.99

0.51

0.30

0.99

0.930

Linear Reg. (N)

0.97

0.53

0.31

0.64

0.950

***

Logistic Reg. (L)

0.50

< 0.01

1.72

0.970

Logistic Reg. (N)

0.85

0.68

0.49

0.64

0.960

Notes. L/N in Task-setting column refers to data generation mode Linear/Nonlinear. Ratio columns report Width(PPI-MOE)/Width(comparator), with values below 1 indicating shorter intervals for PPI-MOE; values below 0.01 are shown as < 0.01. The last column gives a coverage-agreement code, where more stars mean closer p agreement of empirical coverage pb with the nominal target 0.95. With z = |b p − 0.95|/ 0.95(1 − 0.95)/R, we use ∗ ∗ ∗ for z ≤ 1.00, ∗∗ for 1.00 < z ≤ 1.64, ∗ for 1.64 < z ≤ 1.96, . for 1.96 < z ≤ 2.56, and leave the entry blank otherwise.

• For mean and quantile estimation, it achieves near-nominal coverage while yielding confidence intervals that are essentially as short as those of PPI-best and markedly shorter than those of the conventional estimator, PPI-mean, and PPI-worst. • For linear regression, PPI-MOE remains close to PPI-best in both settings. Its gain over the conventional estimator is limited in the well-specified linear setting (MOE/Conv = 0.99), where the conventional procedure is already near variance-optimal and leaves little room for further improvement; under nonlinear misspecification, however, PPI-MOE substantially shortens the interval width relative to the conventional estimator (MOE/Conv = 0.64). • A similar contrast appears in logistic regression. In the linear setting, PPI-MOE may offer little gain relative to the conventional estimator, and can even be slightly worse (MOE/Conv = 1.72), likely because the Bernoulli response introduces additional finite-sample variability that offsets the variance reduction from auxiliary predictions. Under nonlinear misspecification, by contrast, PPI-MOE again yields a clear improvement over the conventional estimator

(MOE/Conv = 0.64) while remaining far more stable than PPI-mean and PPI-worst. Overall, these results show that PPI-MOE stays close to the best-performing baseline while avoiding the severe instability of more aggressive alternatives, and that its gains are especially pronounced under misspecification. Overall, these results suggest that PPI-MOE preserves the efficiency gains of the strongest PPI baseline while being markedly more robust than simple averaging- or worst-case aggregation strategies. Full results are avaliable in Section B.1. 4.1.6

Sample size efficiency

A central practical question is how many labeled samples are needed to reach a target testing power when unlabeled data are abundant. We consider a mean-inference hypothesis test, H0 : µ = 0 versus H1 : µ ̸= 0, with true mean µ = 10. For each method, we increase the labeled sample size n (while fixing N for PPI methods) until the empirical rejection probability (power) reaches a target level, such as 80%. A method is considered more label-efficient if it requires fewer labeled samples to attain this target power.

Figure 4: Sample-size efficiency for mean inference under the linear setting. Figure 4 shows that all PPI variants improve power efficiency over Conventional, with the strongest gains from PPI-MOE and PPI-best. The leading term still scales as n−1 , but its coefficient is substantially reduced through rectification (from Var(Y ) toward Var(Y −Yb )), while the additional

N −1 imputation term is negligible when N ≫ n. In our experiments, PPI-MOE reaches target power with the fewest labeled samples and in some configurations slightly outperforms PPI-best.

4.2 4.2.1

Real data experiments Setup

We evaluate the methods on California Housing and Bike Sharing, using the same predictor pool and the same method definitions as in simulation. We report variance-based metrics for mean, median, and linear-coefficient inference, and then compare sample-size efficiency via power analysis under fixed unlabeled budgets. 4.2.2

Variance and power results

Figure 5: Real-data variance comparison across inference tasks. Figure 5 shows that on both datasets, PPI-MOE consistently attains the lowest variance (or is tied for the lowest) across tasks, while PPI-worst can substantially underperform due to poor predictor quality. This pattern reinforces the need for adaptive weighting rather than single-predictor selection.

Figure 6: Real-data power analysis and minimum labeled sample size for 80% power. Power curves in Figure 6 show that PPI-MOE and PPI-best require far fewer labeled samples than Conventional, with PPI-MOE frequently matching or exceeding PPI-best throughout the curve. These gains are operationally meaningful: the MOE strategy improves detection power while reducing annotation cost, without requiring oracle knowledge of the best single predictor. Note that PPI-worst performs considerably worse than the conventional approach under model misspecification. This suggests that blindly applying PPI may be detrimental, whereas the proposed MOE-based PPI provides a more reliable alternative.

4.3

Empirical summary

The numerical evidence strongly supports the central claims of this paper. 1. Near-oracle adaptation without oracle selection. Across tasks, PPI-MOE consistently achieves (and occasionally exceeds) PPI-best, whose predictor is oracle under well-specification, while avoiding infeasible oracle model selection. 2. Robustness under misspecification. When single predictors degrade under nonlinear misspecification, PPI-MOE remains stable and typically outperforms fixed single-model PPI baselines.

3. Reliable uncertainty reduction. For mean and quantile inference, PPI-MOE delivers substantial variance and interval-width reductions, with coverage close to nominal in most regimes. 4. Large practical gains in labeling efficiency. In the high-unlabeled regime, PPI-MOE attains target power with dramatically fewer labeled samples, yielding meaningful cost savings in data collection. Overall, these experiments show that inference-oriented expert aggregation is more effective than choosing a single predictor or averaging predictors uniformly. The MOE design provides a practical, adaptive, and robust route to high-quality prediction-powered inference.

Discussions

We introduce a general MOE-powered inference framework that exploits unlabeled data to achieve more efficient statistical inference in the presence of multiple predictors with unknown predictive performance. Our framework is flexible and can robustly incorporate additional predictors whenever they become available. Compared with the PPI framework, it leverages the collective predictive power of multiple predictors and enjoys a best-expert guarantee relative to standard PPI. We establish non-asymptotic theory for the normal approximation of MOE-powered estimators and for the coverage properties of MOE-powered confidence sets. For simplicity, we focus primarily on global weighting, under which the mixture-of-experts framework reduces to classical model averaging in the statistical literature. The full power of MOE is typically realized by assigning covariate-dependent weights to the experts; for example, P one may consider the linear mixture Fβ (x) = K k=1 βk (x)fk (x). This reflects the belief that each expert has a domain-specific advantage over the others. However, adaptive weights pose significant challenges for estimation and inference. One possible approach is to parameterize the weight vector {βk (x) : k ∈ [K]} through a softmax function and estimate it by maximizing variance reduction. Doing so would inevitably introduce additional bias and complicate the analysis of normal approximation and coverage properties. We leave these questions for future work.

References Josh Achiam, Steven Adler, Sandhini Agarwal, Lama Ahmad, Ilge Akkaya, Florencia Leoni Aleman, Diogo Almeida, Janko Altenschmidt, Sam Altman, Shyamal Anadkat, et al. Gpt-4 technical report. arXiv preprint arXiv:2303.08774, 2023. Zeyuan Allen-Zhu, Yuanzhi Li, and Zhao Song. A convergence theory for deep learning via overparameterization. In International conference on machine learning, pages 242–252. PMLR, 2019. Anastasios N Angelopoulos, Stephen Bates, Clara Fannjiang, Michael I Jordan, and Tijana Zrnic. Prediction-powered inference. Science, 382(6671):669–674, 2023a. Anastasios N Angelopoulos, John C Duchi, and Tijana Zrnic. Ppi++: Efficient prediction-powered inference. arXiv preprint arXiv:2311.01453, 2023b. Olivier Bousquet. A bennett concentration inequality and its application to suprema of empirical processes. Comptes Rendus Mathematique, 334(6):495–500, 2002. Abhishek Chakrabortty and Tianxi Cai. Efficient and adaptive linear regression in semi-supervised settings. The Annals of Statistics, 46(4):1541 – 1572, 2018. doi: 10.1214/17-AOS1594. URL https://doi.org/10.1214/17-AOS1594. Olivier Chapelle, Bernhard Scholkopf, and Alexander Zien. Semi-supervised learning (chapelle, o. et al., eds.; 2006)[book reviews]. IEEE Transactions on Neural Networks, 20(3):542–542, 2009. Gerda Claeskens and Nils Lid Hjort. Model selection and model averaging. Cambridge books, 2008. Damai Dai, Chengqi Deng, Chenggang Zhao, RX Xu, Huazuo Gao, Deli Chen, Jiashi Li, Wangding Zeng, Xingkai Yu, Yu Wu, et al. Deepseekmoe: Towards ultimate expert specialization in mixture-of-experts language models. In Proceedings of the 62nd Annual Meeting of the Association for Computational Linguistics (Volume 1: Long Papers), pages 1280–1297, 2024. Jerome H Friedman. Multivariate adaptive regression splines. The annals of statistics, 19(1):1–67, 1991. Sara A Geer. Empirical Processes in M-estimation, volume 6. Cambridge university press, 2000.

Robert A Jacobs, Michael I Jordan, Steven J Nowlan, and Geoffrey E Hinton. Adaptive mixtures of local experts. Neural computation, 3(1):79–87, 1991. Vladimir Koltchinskii. Oracle inequalities in empirical risk minimization and sparse recovery problems: Ecole D’Eté de Probabilités de Saint-Flour XXXVIII-2008, volume 2033. Springer, 2011. Vladimir Koltchinskii and Dong Xia. Perturbation of linear forms of singular vectors under gaussian noise. In High Dimensional Probability VII: The Cargèse Volume, pages 397–423. Springer, 2016. Xianhui Liu, Tong Wang, Yuanyuan Lin, and Zhanfeng Wang. Semi-supervised inference for the high-dimensional quantile regression. Science China Mathematics, pages 1–30, 2025. Fedor Nazarov. On the maximal perimeter of a convex set in with respect to a gaussian measure. In Geometric Aspects of Functional Analysis: Israel Seminar 2001-2002, pages 169–187. Springer, 2004. Valentin V Petrov. Sums of independent random variables. Springer Science & Business Media, 2012. Adrian E Raftery, David Madigan, and Jennifer A Hoeting. Bayesian model averaging for linear regression models. Journal of the American Statistical Association, 92(437):179–191, 1997. Martin Raič. A multivariate berry–esseen theorem with explicit constants. Bernoulli, 25(4A): 2824–2853, 2019. Gemini Team, Rohan Anil, Sebastian Borgeaud, Jean-Baptiste Alayrac, Jiahui Yu, Radu Soricut, Johan Schalkwyk, Andrew M Dai, Anja Hauth, Katie Millican, et al. Gemini: a family of highly capable multimodal models. arXiv preprint arXiv:2312.11805, 2023. T Tony Cai and Zijian Guo. Semisupervised inference for explained variance in high dimensional linear regression and its applications. Journal of the Royal Statistical Society Series B: Statistical Methodology, 82(2):391–419, 2020. Joel A Tropp. User-friendly tail bounds for sums of random matrices. Foundations of computational mathematics, 12(4):389–434, 2012.

Halbert White. Nonlinear regression on cross-section data. Econometrica: Journal of the econometric society, pages 721–746, 1980. Dong Xia and Ming Yuan. On polynomial time methods for exact low-rank tensor completion. Foundations of Computational Mathematics, 19(6):1265–1313, 2019. Zichun Xu, Daniela Witten, and Ali Shojaie. A unified framework for semiparametrically efficient semi-supervised learning. arXiv preprint arXiv:2502.17741, 2025. An Yang, Anfeng Li, Baosong Yang, Beichen Zhang, Binyuan Hui, Bo Zheng, Bowen Yu, Chang Gao, Chengen Huang, Chenxu Lv, et al.

Qwen3 technical report.

arXiv preprint

arXiv:2505.09388, 2025. Anru Zhang, Lawrence D. Brown, and T. Tony Cai. Semi-supervised inference: General theory and estimation of means. The Annals of Statistics, 47(5):2538 – 2566, 2019. doi: 10.1214/ 18-AOS1756. URL https://doi.org/10.1214/18-AOS1756. Yuqian Zhang and Jelena Bradic. High-dimensional semi-supervised learning: in search of optimal inference of the mean. Biometrika, 109(2):387–403, 2022. Xiaojin Jerry Zhu. Semi-supervised learning literature survey. 2005. Tijana Zrnic and Emmanuel J Candès. Cross-prediction-powered inference. Proceedings of the National Academy of Sciences, 121(15):e2322083121, 2024.

Proof of Lemmas and Theorems

A.1

Proof of Lemma 1

By definition, we have Qn (βbn ) ≤ Qn (β∗ ) which implies that Q(βbn ) − Q(β∗ ) ≤ Q(βbn ) − Qn (βbn ) − Q(β∗ ) − Qn (β∗ ) .

For notational brevity, denote m(β) = mθ (β, X, Y ), mi (β) = mθ (β, Xi , Yi ) for ∀i ∈ [n], m̄(β) = P n−1 ni=1 mθ (β, Xi , Yi ), and M(β) = Emθ (β, X, Y ). For any β ∈ B, we write n n 1 X ⊤ Qn (β) − Q(β) = · tr mi (β)m⊤ i (β) − Etr m(β)m (β) n−1 n | i=1 {z } J1 (β)

i n h + tr M(β)M⊤ (β) − tr m̄(β)m̄⊤ (β) n−1| {z } J2 (β)

i 1 h + Etr m(β)m⊤ (β) − tr M(β)M⊤ (β) , n−1| {z } =O(C12 )

where the last term is bounded by O(U 2 /n) because ∥m(β)∥ ≤ U for all β ∈ B and (X, Y ) ∈ X × Y under Assumption 1. Therefore, under the quadratic growth condition of Assumption 1, we get 2 n U b J2 (β∗ ) − J2 (βn ) + O . n−1 n Denote g(β) = tr m(β)m⊤ (β) and gi (β) = tr mi (β)m⊤ i (β) . Then, τ0 ∥βbn − β∗ ∥2 ≤ J1 (β∗ ) − J1 (βbn ) +

(35)

J1 (β∗ ) − J1 (βbn ) = Pn g(β∗ ) − g(βbn ) − P g(β∗ ) − g(βbn ) , P where P g(β∗ ) − g(βbn ) := E g(β∗ ) − g(βbn ) and Pn g(β∗ ) − g(βbn ) = n−1 ni=1 gi (β∗ ) − gi (βbn ) . For any δ > 0, define Pn − P g(β) − g(β∗ ) ,

γn (δ) := sup β∈B∗ (δ)

which represents the supremum of an empirical process indexed by β ∈ B∗ (δ) with B∗ (δ) := β ∈ B : ∥β − β∗ ∥ ≤ δ . Then, J1 (βbn ) − J1 (β∗ ) ≤ γn ∥βbn − β∗ ∥ . It remains to develop a uniform upper bound for γn (δ). Lemma 7. Suppose Assumption 1 holds. Then, there exist absolute constants C1 > 0 such that, for any δ ∈ n−1 , D , r γn (δ) ≤ C1 δτ1 U ·

q + n

which holds with probability at least 1 − n−10 . 39

log(Dn) log(Dn) + , n n

Denote the event in Lemma 7 by E1 . Therefore, on event E1 , r r q log(Dn) J1 (βbn ) − J1 (β∗ ) ≤ C1 ∥βbn − β∗ ∥τ1 U · , + n n assuming that n ≫ log(Dn). It remains to bound J2 (β∗ ) − J2 (βbn ) . By definition, 2

J2 (β) = −2 M(β), m̄(β) − M(β) − m̄(β) − M(β) , and as a result, J2 (βbn ) − J2 (β∗ ) ≤ 2 M(βbn ), m̄(βbn ) − M(βbn ) − M(β∗ ), m̄(β∗ ) − M(β∗ ) + ∥m̄(βbn ) − M(βbn )∥2 + ∥m̄(β∗ ) − M(β∗ )∥2 .

(36)

For any δ ∈ [n−1 , D], define (Pn − P m(β) − m(β∗ ) .

ψn (δ) := sup β∈B∗ (δ)

The proof of the following lemma is almost identical to that of Lemma 7 and hence skipped. Recall that m(β) ∈ Rp . Lemma 8. Suppose Assumption 1 holds. Then, there exist absolute constant C1 > 0 such that, for any δ ∈ [n−1 , D], r ψn (δ) ≤ C1 τ1 δ

pq + n

p log(Dn) p log(Dn) + , n n

which holds with probability at least 1 − n−10 . Denote the event in Lemma 8 by E2 . On event E2 , we have m̄(βbn ) − m̄(β∗ ) − M(βbn ) − M(β∗ ) ≤ C1 ∥βbn − β∗ ∥τ1 ·

pq + n

p log(Dn) , n

assuming n ≫ p log(Dn). Moreover, by Bernstein inequality, there exists an event E3 with P(E3 ) ≥ 1 − n−10 on which, r m̄(β∗ ) − M(β∗ ) ≤ C2 U

log n , n

for some constant C2 > 0. It implies, in the event E2 ∩ E3 , that r r r pq p log(Dn) log n m̄(βbn ) − M(βbn ) ≤ C1 ∥βbn − β∗ ∥τ1 · + . + C2 U n n n 40

Continuing from (36), under Assumption 1 and the bound 2 M(βbn ), m̄(βbn ) − M(βbn ) − M(β∗ ), m̄(β∗ ) − M(β∗ ) ≤2 M(βbn ) − M(β∗ ), m̄(β∗ ) − M(β∗ ) + 2 M(βbn ), (Pn − P)(m(βbn ) − m(β∗ )) r r r log n pq p log(Dn) b b ≤2 M(βn ) − M(β∗ ) · C2 U + 2C1 U ∥βn − β∗ ∥τ1 · + n n n r r r log n pq p log(Dn) , ≲τ1 U ∥βbn − β∗ ∥ + U ∥βbn − β∗ ∥τ1 · + n n n we get p log(Dn) 2 log n 2 2 pq b b J2 (βn ) − J2 (β∗ ) ≲ U + ∥βn − β∗ ∥ τ1 + n n n r r pq p log(Dn) b + τ1 U ∥βn − β∗ ∥ + . n n By (35) and the bounds of J1 and J2 , in the event E1 ∩ E2 ∩ E3 , we get −1 τ12 p(q + log(Dn)) 2 2 τ0 log n b , + 2· ∥βn − β∗ ∥ ≲ U n n τ0 assuming n ≥ C0 (τ12 /τ0 )p q + log(Dn) for a large C0 > 1, which concludes the proof.

A.2

Proof of Lemma 7

Let δ0 := n−1 and define δj = 2j−1 δ0 for j ∈ log(nD) . Following the discretization procedure as in (Xia and Yuan, 2019, Lemma 2), we derive a uniform upper bound for γn (δj ) holding for all j ∈ log(nD) , which can be easily extended for γn (δ) for all δ ∈ [n−1 , D]. Let us derive the upper bound for γn (δj ) with any fixed j. Note that g(β) − g(β∗ ) ≤

sup β∈B∗ (δj )

sup

m(β) − m(β∗ )

m(β) + m(β∗ ) ≤ 2τ1 U δj ,

β∈B∗ (δj )

where the last inequality is due to the Lipschitz and upper bound conditions from Assumption 1. Similarly, sup

Var g(β) − g(β∗ ) ≤

β∈B∗ (δj )

sup

E g(β) − g(β∗ )

≤ 4τ12 U 2 δj2 .

β∈B∗ (δj )

Applying Bousquet’s version of Talagrand’s concentration inequality (Bousquet, 2002), with probability at least 1 − e−t for all t > 0, r γn (δj ) ≤ 2Eγn (δj ) + 2τ1 U δj 41

t t + . n n

It suffices to upper bound Eγn (δj ). By the symmerization inequality (Koltchinskii, 2011), n 1X Eγn (δj ) ≤ 2E sup εi gi (β) − gi (β∗ ) , β∈B∗ (δj ) n i=1

where ε1 , · · · , εn are i.i.d. Rademacher random variables. Denote the function class G∗ (δj ) := g(β) − g(β∗ ) : β ∈ B∗ (δj ) . Conditioned on (X1 , Y1 ), · · · , (Xn , Yn ), and for h = g(β) − g(β∗ ) ∈ G∗ (δj ), we denote h(Xi , Yi ) = gi (β) − gi (β∗ ) and define the distance L2 (Pn ) in G∗ (δj ) by ∥h1 − h2 ∥2L2 (Pn ) :=

n 2 1X h1 (Xi , Yi ) − h2 (Xi , Yi ) , n

∀h1 , h2 ∈ G∗ (δj ).

i=1

By Dudley’s entropy bound (Koltchinskii, 2011, Theorem 3.10), we have Z √2σn q n 1X C1 Eγn (δj ) ≤ 2E sup εi h(Xi , Yi ) ≤ √ E log N G∗ (δj ); L2 (Pn ); ϵ dϵ, n 0 h∈G∗ (δj ) n i=1

where σn2 := suph∈G∗ (δj ) Pn h2 ≤ maxi∈[n] supβ∈B∗ (δj ) |gi (β) − gi (β∗ )|2 ≤ U 2 τ12 δj2 and C1 > 0 is an absolute constant. Note that N (F; d; ϵ) represents the ϵ-covering number of a set F under the distance d(·, ·). If h1 = g(β1 ) − g(β∗ ) and h2 = g(β2 ) − g(β∗ ), we have 1 ∥h1 − h2 ∥2L2 (Pn ) =

n X

gi (β1 ) − gi (β2 )

≤ U 2 τ12 ∥β1 − β2 ∥2 ,

i=1

2δj U τ1 q , where the last inequality implying that N G∗ (δj ); L2 (Pn ); ϵ ≤ N B∗ (δj ), ∥·∥, ϵ/(U τ1 ) ≤ ϵ is due to the fact that B∗ (δj ) ⊂ Rq has a diameter at most 2δj . Plugging it into the Dudley’s entropy bound, we get C1 Eγn (δj ) ≤ √ E n

Z 2U τ1 δj

√

r q

r 2δj U τ1 q log dϵ = O δj U τ1 · . ϵ n

Therefore, with probability at least 1 − e−t , r γn (δj ) ≤ C1 τ1 U δj

q + n

t t + . n n

By setting t = C2 log(Dn) for a large C2 > 0 and taking a union bound for all j ∈ [log(Dn)], we get with probability at least 1 − n−10 , r r q log(Dn) log(Dn) γn (δj ) ≤ C1 τ1 U δj + + , n n n

for all j = 1, 2, · · · , log(Dn).

By adjusting the constant C1 and chaining through all the intervals [Cj , Cj+1 ], we can extend the above bounds to all δ ∈ [n−1 , D], which concludes the proof. 42

A.3

Proof of Theorem 1

b b as eθ,βbn + ∆ By definition, we decompose the estimating function g θ,βn b b eθ,βbn + ∆ g θ,βn = N −1

−1

N X

ei , F b (X ei )) + n−1 gθ (X βn

i=1 N h X

n X

gθ (Xi , Yi ) − gθ (Xi , Fβbn (Xi ))

i=1 N h i i X −1 ei , F b (X e e ei )) − gθ (X ei , Fβ (X ei )) gθ (X gθ (Xi , Fβ∗ (Xi )) − Egθ (X, Fβ∗ (X)) + N ∗ βn i=1

i=1

e N −1/2 Z

}

Je1

N,θ,β∗

n h X i + n−1 gθ (Xi , Yi ) − gθ (Xi , Fβ∗ (Xi )) − Egθ (X, Y ) − Egθ (X, Fβ∗ (X)) i=1

| − n−1

}

n−1/2 Zn,θ,β∗ n h X

i gθ (Xi , Fβbn (Xi )) − gθ (Xi , Fβ∗ (Xi )) +Egθ (X, Y ).

i=1

} (37)

e N,θ,β , the existence of finite Cov gθ (X, Fβ (X) and the boundBy the definition of Zn,θ,β∗ and Z ∗ ∗ edness of mθ (β, X, Y ), the Central Limit Theorem indicates that d. d. e N,θ,β → (X)) Zn,θ,β∗ → N 0, Cov mθ (β∗ , X, Y ) and Z N 0, Cov g (X, F θ β ∗ ∗

(38)

b b . Note that eθ,βbn + ∆ Step 1: bias. First, we derive the bound of the bias of g θ,βn b b − Egθ (X, Y ) = EJe1 − EJ1 . eθ,βbn + ∆ E g θ,βn By Lemma 7, under the event in Lemma 1 with probability at least 1 − n−9 , we have Je1 − J1 ≤

Pn − P m(βbn ) − m(β∗ ) +

PN − P m(βbn ) − m(β∗ )

≤ ψn (δn ) + ψN (δn ) r r r r pq p log(Dn) p log(Dn) pq p log(DN ) p log(DN ) ≤ C1 τ1 δn + + + + + n n n N N N r p(q + log(Dn)) p2 (q + log(Dn))(q + log(DN )) ≲ + . (39) n nN q q τ0−1 log n p(q+log(Dn)) τ1 + τ0 · for a sufficiently large C2 > 0. where we set δn = C2 U n n 43

If n ≳ p q + log(Dn) and N ≳ n, then p(q + log(Dn)) EJe1 − EJ1 ≤ E∥Je1 − J1 ≲ . n

(40)

Take (40) back into the expectation form of (37), combining with (38), we finally have the bias as b b − Egθ (X, Y ) = O p(q + log(Dn)) . eθ,βbn + ∆ E g θ,βn n

(41)

Especially, if take θ = θ∗ so that Egθ∗ (X, Y ) = 0, we have b b eθ∗ ,βbn + ∆ E g θ∗ ,βn

p(q + log(Dn)) n

Step 2: normal approximation. By (37) and (39), we have r √ p q + log(Dn) ne b e √ eθ,βbn + ∆θ,βbn − Egθ (X, Y ) = Zn,θ,β∗ + ZN,θ,β∗ + Op n g . N n Bernstein Inequality implies there exists an event E1 with P(E1 ) ≥ 1 − n−10 , on which e N,θ,β ≤ O Z ∗

p p log n ,

implying that √ p(q + log(Dn)) b e √ eθ,βbn + ∆θ,βbn − Egθ (X, Y ) = Zn,θ,β∗ + Op n g + n

pn log(n) . N

ep (·), we get By Berry-Esseen bound and the high probability bound inherited from O r √ ⊤ e b b n es gθ∗ ,βbn + ∆ pn log(n) p(q + log(Dn)) θ ∗ ,β n q √ sup P + , ≤t =O N n t∈R e⊤ s Wθ∗ ,Y −Fβ∗ es where Wθ,Y −Fβ∗ = Cov gθ (X, Y ) − gθ (X, Fβ∗ (X)) . Therefore, for all s ∈ [p], √ P

b b eθ∗ ,βbn + ∆ n e⊤ s g θ ∗ ,β n q ⊤ es Wθ∗ ,Y −Fβ∗ es

≤ zα/(2p)

r α p(q + log(Dn)) pn log(n) √ ≥1− +O + . p N n

(42)

c θ ,Y −F − Wθ ,Y −F : It remains to bound W ∗ ∗ β∗ b β n

c θ ,Y −F − Wθ ,Y −F W ∗ ∗ β∗ b β n

= SCov mθ∗ (βbn , X, Y ) − Cov mθ∗ (β∗ , X, Y ) = n

−1

n X

mθ∗ (βbn , Xi , Yi )⊗2 − n−1

i=1

n−1

+ n−1

n X

mθ∗ (β∗ , Xi , Yi )⊗2

i=1 n X

n i⊗2 h i⊗2 X − n−1 mθ∗ (βbn , Xi , Yi ) mθ∗ (β∗ , Xi , Yi )

i=1 n X

i=1

mθ∗ (β∗ , Xi , Yi )⊗2 − Emθ∗ (β∗ , X, Y )⊗2

i=1 n X −1

i⊗2 h i⊗2 mθ∗ (β∗ , Xi , Yi ) − Emθ∗ (β∗ , X, Y ) ,

i=1

where x⊗2 = xx⊤ refers to the Kronecker product. By the Lipschitz condition of mθ (β, X, Y ) with respect to β, we have ( n n X X max n−1 mθ (βbn , Xi , Yi )⊗2 − n−1 mθ (β∗ , Xi , Yi )⊗2 ∗

h + n−1

∗

i=1 n X

mθ∗ (βbn , Xi , Yi )

i⊗2

− n

i=1 n X −1

i=1

i⊗2 mθ∗ (β∗ , Xi , Yi )

) (43)

i=1

≲ 2τ1 U βbn − β∗ = O

p log(Dn) n

under the event of Lemma 1. Since mθ (β, X, Y ) is bounded, by matrix Bernstein inequality, there exist an event E2 with P(E2 ) ≥ 1 − n−10 , on which ( n X max n−1 mθ∗ (β∗ , Xi , Yi )⊗2 − Emθ∗ (β∗ , X, Y )⊗2 i=1 n h i⊗2 h i⊗2 X −1 n mθ∗ (β∗ , Xi , Yi ) − Emθ∗ (β∗ , X, Y )

) =O

i=1

Combining (43) and (44), we conclude that c θ ,Y −F − Wθ ,Y −F W =O ∗ ∗ β∗ b β n

p log(Dn) . n

r p log n n

(44) .

For all s ∈ [p], √ ⊤ b b eθ∗ ,βbn + ∆ n es g θ∗ ,βn q c e⊤ s Wθ∗ ,Y −Fβbn es ! √ ⊤ ⊤ W c b e − W es eθ∗ ,βbn + ∆θ∗ ,βbn n es g θ∗ ,Y −Fβ∗ θ∗ ,Y −Fβb s n q 1+ q = q q ⊤ ⊤ ⊤ ⊤ es Wθ∗ ,Y −Fβ∗ es c c es Wθ∗ ,Y −Fβb es es Wθ∗ ,Y −Fβ∗ es + es Wθ∗ ,Y −Fβb es n n √ ⊤ b b log n eθ∗ ,βbn + ∆ n es g θ∗ ,βn ep p √ = q +O . n e⊤ W e θ∗ ,Y −Fβ∗ s s Plugging into (42), we have r √ ⊤ e b b n es gθ∗ ,βbn + ∆ p(q + log(Dn)) pn log(n) θ ∗ ,β n q √ , sup P ≤t =O + N n c t∈R e⊤ s Wθ∗ ,Y −Fβb es n

implying √ P

b b eθ∗ ,βbn + ∆ n e⊤ s g θ ∗ ,β n q c e⊤ s Wθ∗ ,Y −Fβb es

≤ zα/(2p)

r α p(q + log(Dn)) pn log(n) √ ≥1− +O + . p N n

The coverage probability satisfies that ( √ ⊤ )! p b b \ eθ∗ ,βbn + ∆ n es g θ ∗ ,β n MOE q P θ∗ ∈ Cα =P ≤ zα/(2p) ⊤W c e e s=1 θ∗ ,Y −Fβb s s n )! ( √ p ⊤ b b X eθ∗ ,βbn + ∆ n es g θ∗ ,βn q > zα/(2p) ≥1− P c θ ,Y −F es e⊤ W s=1 s ∗ bn β r !p p(q + log(Dn)) pn log(n) α √ +O + =1− p N n r 2 p (q + log(Dn)) p3 n log(n) √ =1−α+O + , N n if n ≳ p4 q 2 + log2 (Dn) and N ≳ p3 n log n, which concludes the proof.

A.4

Proof of Lemma 2

By (14) and the fact En−1 1⊤ n y = θ∗ , we have h i e − n−1 1⊤ F (F⊤ Pn F)−1 F⊤ Pn y EθbMOE − θ∗ =E N −1 1⊤ F N n h i e − n−1 1⊤ F (F⊤ Pn F)−1 F⊤ Pn y − Fβ∗ , =E N −1 1⊤ F N n 46

(45)

−1 ⊤ ⊤ e where the last equality holds since E[N −1 1⊤ N F] = E[n 1n F] = Ef1 . Then, ⊤ −1 ⊤ −1 ⊤ e EθbMOE − θ∗ ≤ E N −1 1⊤ N F − n 1n F · (F Pn F) F Pn y − Fβ∗

(46)

By the matrix Bernstein inequality (Tropp, 2012), there is an event E0 with P(E0 ) ≥ 1 − n−10 on which, r n 1X ⊤ 1 ⊤ K log n ⊤ max fi fi − Eff , F y − EY f , f̄ − Ef ≤ C1 , n n n i=1 √ and |Ȳ − EY | ≤ C1 / n for some absolute constant C1 > 0, where Ȳ = n−1 1⊤ n y. Therefore,

conditioned on event E0 , we have X n 1 ⊤ ⊤ F Pn F = n · fi fi − f̄ f̄ n ⊤

i=1

√ implying that F⊤ Pn F − n · Cov(f ) F = O K n log n . Similarly, conditioned on event E0 , we have n p 1X F Pn y − n · Cov(Y, f ) = n · Kn log n . Yi fi − Ȳ f̄ − Cov(Y, f ) = O n ⊤

i=1

As a result, conditioned on E0 and if n ≥ C2 K 2 log n for a large constant C2 > 0, we have (F⊤ Pn F)−1 = O(n−1 ) and F⊤ Pn (y − Fβ∗ ) ≤ F⊤ Pn y − n · Cov(Y, f ) + p ≤C1 K 3 n log n.

F⊤ Pn F − n · Cov(f ) β∗

By Bernstein inequality, there exist an event E1 with P(E1 ) ≥ 1 − n−10 on which, r −1 −1 ⊤ e − n 1 F ≤ C1 K log n , N 1N F n n we we assumed N ≫ n. Finally, continuing from (46), we get −1 ⊤ ⊤ −1 ⊤ e EθbMOE − θ∗ ≤E N −1 1⊤ N F − n 1n F · (F Pn F) F Pn y − Fβ∗

IE0 ∩E1

−1 ⊤ ⊤ −1 ⊤ e +E N −1 1⊤ N F − n 1n F · (F Pn F) F Pn y − Fβ∗

IE0c ∪E1c ≤ C1 ·

which concludes the proof.

K 2 log n , n

A.5

Proof of Theorem 3

−1 Recall that βbn = (F⊤ Pn F)−1 F⊤ Pn y and β∗ = Cov(f ) Cov(f , Y ). We first derive the upper bound for ∥βbn − β∗ ∥. Note that −1 ⊤ −1 −1 ⊤ βbn − β∗ = (F⊤ Pn F)−1 − nCov(f ) F Pn y + Cov(f ) n F Pn y − Cov(f , Y ) . As show in the proof of Lemma 2, there exists an event E0 with P(E0 ) ≥ 1 − n−10 , on which the following bounds hold p F⊤ Pn F − nCov(f ) F = O K n log n

and

F⊤ Pn y − nCov(f , Y ) = O

p Kn log n ,

and consequently ∥F⊤ Pn y∥ = O(n). Therefore, conditioned on event E0 , we get −1 nCov(f ) ∥F⊤ Pn y∥ r 3 K log n ⊤ F Pn y − nCov(f , Y ) = O , n

βbn − β∗ ≤∥(F⊤ Pn F)−1 ∥∥F⊤ Pn F − nCov(f )∥ +

−1 nCov(f )

(47)

where we assumed n ≫ K 2 log n. Denote the optimal MOE by f∗ (X) = ⟨β∗ , f ⟩ =

k=1 β∗,k fk (X).

We show that the difference

between θbMOE and θbfPPI (recall its definition in (10)) is negligible. By definition, ∗ = θbMOE − θbfPPI ∗

−1 ⊤ ∗ b e N −1 1⊤ N F − n 1n F (βn − β )

−1 ⊤ e ≤ N −1 1⊤ N F − n 1n F

βbn − β∗ .

By Bernstein inequality, there is an event E1 with P(E1 ) ≥ 1 − n−10 , on which the following bounds hold e N −1 1⊤ N F − Ef = O

p (K/N ) log n ,

and

n−1 1⊤ n F − Ef = O

Together with (47), conditioned on event E0 ∩ E1 , we get 2 K log n θbMOE − θbfPPI = O . ∗ n

p (K/n) log n .

(48)

Thus it suffices to focus on the classical PPI-based estimator θbfPPI equipped with the optimal MOE ∗ as the predictor. −1 ⊤ e Recall θbfPPI = N −1 1⊤ N Fβ∗ − n 1n Fβ∗ − y implying that ∗ √

√ n θbfPPI − θ∗ = n ∗

N n X 1 X ei ) − Ef∗ (X) − √1 f∗ (X f∗ (Xi ) − Yi − E f∗ (X) − Y (49) N n i=1

i=1

The classical central limit theorem dictates that n d 1 X − → N 0, Var f∗ (X) − Y Zn,f∗ := − √ f∗ (Xi ) − Yi − E f∗ (X) − Y n i=1

and the Chebyshev’s inequality gives that X r N √ n 1 e n f∗ (Xi ) − Ef∗ (X) = Op · Var f (X) . N N i=1

Moreover, by Bernstein inequality, we get √

ep n θbfPPI − θ∗ = Zn,f∗ + O ∗

n log n N

where Zn,f∗ →d N 0, Var(Y − f∗ (X)) as n → ∞. Bound (48) shows that 2 √ MOE K log n PPI b b e √ . n θ − θ f∗ = O p n Finally, we conclude that √

MOE

n θb

ep − θ∗ = Zn,f∗ + O

K 2 log n √ + n

n log n . N

ep (·), we get By the Berry-Esseen bound and the high probability bound inherited from O r √ bMOE 2 n(θ − θ∗ ) K log n n log n √ sup P ≤ t − Φ(t) = O , + σY −f∗ N n t∈R

(50)

(51)

where σY2 −f∗ := Var(Y − f∗ (X)). Therefore, r 2 h i σ σ K log n n log n Y −f Y −f ∗ ∗ MOE MOE √ P θ∗ ∈ θb − zα/2 √ , θb + zα/2 √ =1−α+O + . N n n n bY2 −f∗ − σY2 −f∗ . Define It suffices to bound σ εbi := Yi − ⟨fi , βbn ⟩

ε∗i := Yi − ⟨fi , β∗ ⟩,

and

∀i ∈ [n].

By definition, we write n

1X 2 σ bY2 −f∗ = εbi − n i=1

1X εbi n

!2 ,

i=1

and as a result # " n # n n X X X 1 1 1 σ bY2 −f∗ − σY2 −f∗ = εb2i − (ε∗i )2 + (ε∗i )2 − E(ε∗1 )2 n n n i=1 i=1 i=1   !2 !2   !2 n n n X X X 1 1 1 − εbi − ε∗i  −  ε∗i − (Eε∗1 )2  . n n n "

i=1

(52)

By (47), conditioned on event E0 , we get r

εbi − ε∗i

= ⟨fi , βbn − β∗ ⟩ = O

K 4 log n . n

By Bernstein inequality, there exists an event E2 with P(E2 ) ≥ 1 − n−10 , on which the following bounds hold n

1X ∗ 2 (εi ) − E(ε∗1 )2 = O n

i=1

K 4 log n n

1X ∗ εi − Eε∗1 = O n

and

i=1

K 2 log n , n

where we used the fact |ε∗i | = O(K) under Assumption 2. These bounds imply that σ bY2 −f∗ − σY2 −f∗

r ep =O

K 4 log n n

and

σ bY −f∗ − σY −f∗ ep =O σ bY −f∗

K 4 log n , n

where the second inequality holds assuming σY2 −f∗ > 0 and n ≫ K 4 log n. Observe that √

n(θbMOE − θ∗ ) = σ bY −f∗

√

bY −f∗ σY −f∗ − σ n(θbMOE − θ∗ ) · 1+ . σY −f∗ σ bY −f∗

(53)

Applying Bernstein inequality to Zn,f∗ in (50), we get that there exists an event E3 with P(E3 ) ≥ 1 − n−10 , on which the following bound holds √ bMOE p n(θ − θ∗ ) = O log n , σY −f∗ assuming that n ≫ K 4 log2 n and N ≫ n log n. Together with (53), we get √

n(θbMOE − θ∗ ) = σ bY −f∗

√

2 K log n n(θbMOE − θ∗ ) e √ + Op . σY −f∗ n

Plugging into (51), we conclude that r √ bMOE 2 n(θ − θ∗ ) K log n n log n √ sup P ≤ t − Φ(t) = O + , σ bY −f∗ N n t∈R implying that r 2 h i σ b σ b K log n n log n Y −f Y −f ∗ ∗ MOE MOE b b √ P θ∗ ∈ θ − zα/2 √ , θ + zα/2 √ =1−α+O + , N n n n which concludes the proof.

(54)

A.6

Proof of Theorem 4

From (48) and (49) in the Proof of Theorem 3, we get √

√ X n N 1 X n MOE b ei ) − Ef (X) e n θ − θ∗ = − √ f∗ (Xi ) − Yi − E f∗ (X) − Y + f (X N n i=1 i=1 K 2 log n ep √ +O . n

eN , we get Applying the Berry-Esseen bound to both Zn,f∗ and Z 1 1 and sup F ZeN (t) − Φ(t) = O √ , sup F Zn,f∗ (t) − Φ(t) = O √ √ n t t N σY −f n/N σ ∗ f∗ 2 where σf∗ := Var f∗ (X) and FY (·) denotes the c.d.f. of the random variable Y . Moreover, by the convergence rate for the sum of independent non-identically distributed random variables (Petrov, 2012), we get

√

n(θbMOE − θ∗ )

sup P q t σY2 −f∗ + (n/N )σf2∗

2 K log n 1 √ ≤ t − Φ(t) = O +√ . n N

The rest of the proof follows the same as the proof of Theorem 3.

A.7

Proof of Lemma 3

Recall that β∗ = β∗ (θ∗ ) and βbn = βbn (θ∗ ). The local strongly convexity and the separation condition (A2) in Assumption 3 ensures that β∗ is the unique minimizer of Q(θ∗ , β). Step 1: βbn is close to β∗ . By Bernstein inequality, for a given β ∈ B, there exists an event E0 with P(E0 ) ≥ 1 − n−10 , on which max

n 2 2 1 X , Sh (θ∗ − Yi ) − Sh (θ∗ − fi⊤ β) − E Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β) n i=1 ) r n X 1 log n Sh (θ∗ − Yi ) − Sh (θ∗ − fi⊤ β) − E Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β) ≤ C0 . n n i=1

Moreover, due to the Lipschitz continuity of Sh (θ∗ − f ⊤ β) with respect to β, we can apply Talagrand’s concentration inequality and Dudley’s entropy bound as in the proof of Lemma 1 to bound the following empirical process (see a sharper bound in Step 2) r K + log n sup |Qn (θ∗ , β) − Q(θ∗ , β)| ≤ C1 , n β∈B 51

which holds with probability at least 1 − n−10 . We denote this event by Ee0 . On the event Ee0 , we have inf

Qn (θ∗ , β) − Qn (θ∗ , β∗ )

inf

Q(θ∗ , β) − Q(θ∗ , β∗ ) − 2 sup |Qn (θ∗ , β) − Q(θ∗ , β)| ≥ η0 − C1

β∈B:∥β−β∗ ∥≥δ0

r ≥

β∈B:∥β−β∗ ∥≥δ0

β∈B

K + log n , n

which is positive as long as n ≳ η0−2 (K + log n). Since Qn (θ∗ , βbn ) ≤ Qn (θ∗ , β∗ ) by definition, we conclude that ∥βbn − β∗ ∥ < δ0 on event Ee0 , according to condition (A2) of Assumption 3. Step 2: upper bound for Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) . By definition, Qn (θ∗ , βbn ) ≤ Qn (θ∗ , β∗ ) implying that Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) ≤ Qn (θ∗ , βbn ) − Qn (θ∗ , β∗ ) − Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) ≤ 2 (Pn − P) Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) | {z } J1

+ (Pn − P) Sh2 (θ∗ − f ⊤ βbn ) − Sh2 (θ∗ − f ⊤ β∗ ) | {z } J2

−

2 2 Pn Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) − Pn Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) 2

P Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn )

2 + P Sh (θ∗ − Y ) − Sh (θ∗ − f β∗ ) ,

⊤

where we denote the third term on RHS by J3 . Note that the above LHS is lower bounded by 0.5h−1 ∥βbn − β∗ ∥2 · inf β∈B∗ (δ0 ) λmin h∇2 Q(θ∗ , β) , on the event Ee0 , under the condition (A2) of Assumption 3. All the three terms J1 , J2 , J3 can be bounded in a similar fashion. Take the term J1 for example. For any δ ∈ (n−1 , δ0 ), we study the supremum of an empirical processed indexed by β ∈ B∗ (δ): γn (δ) := sup

(Pn − P) Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) .

β∈B∗ (δ)

Then, |J1 | ≤ γn ∥βbn − β∗ ∥ . It suffices to develop an upper bound of γn (δ) uniformly for all δ ∈ (n−1 , δ0 ). Fix a δ ∈ (n−1 , δ0 ). Under Assumption 3, we have sup

Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) ≤ h−1 U sup ∥β − β∗ ∥ ≤ h−1 U δ.

β∈B∗ (δ)

Moreover, for any β ∈ B∗ (δ), we have 2 ESh2 (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) 2 ≤ E Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) Z 1 2 ⊤ −1 ′ θ∗ − f βt ⊤ =E h S f (β − β∗ )dt h 0 Z 1 ∥f ∥ ∥β − β∗ ∥ ′ θ∗ − f ⊤ βt 2 S dt ≤E h h 0 Z 1 ∥f ∥2 ∥β − β∗ ∥2 ′ θ∗ − f ⊤ βt 2 ≤E S dt h2 h 0 Z U 2 ∥β − β∗ ∥2 1 h ′ θ∗ − f ⊤ βt i2 E S dt = h2 h 0 Z ⊤ U 2 ∥β − β∗ ∥2 1 ′ θ∗ − f β t ≤ ES dt h2 h 0 where βt := β∗ + t(β − β∗ ) and S(t) := 1 + exp(−t)

(Sh (t) ≤ 1) (Fundamental theorem of calculus)

(Cauchy-Schwarz inequality) (Fubini’s theorem) (S ′ (t) < 1) −1

. It’s sufficient to control ES ′

θ−f ⊤ β h

for

all β ∈ B∗ (δ). By Assumption 3, Zβt := f ⊤ βt has a bounded density function pβ (z) for all β ∈ B∗ (δ). Thus, Z ⊤ ⊤ ′ θ∗ − f β t ′ θ ∗ − f βt = S pβt (z)dz ES h h Z Z ′ = h S (u)pβt (θ∗ − hu)du ≤ C0 h S ′ (u)du ≤ C2 h, where C0 is the upper bound of the density in Assumption 3. 2 Consequently, we have ESh2 (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) ≲ h−1 U 2 ∥β − β∗ ∥2 . This indicates that there exists an constant C1 > 0 such that C 1 U 2 Var Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) ≤ ∥β − β∗ ∥2 , h

∀β ∈ B∗ (δ).

(55)

Applying Bousquet’s version of Talagrand’s concentration inequality (Bousquet, 2002), with probability at least 1 − e−t for all t > 0, r γn (δ) ≤ 2Eγn (δ) + 2U δ

t t + . nh nh

Similarly, we need to bound Eγn (δ). By the symmerization inequality (Koltchinskii, 2011), n 1X εi Sh (θ∗ − Yi ) Sh (θ∗ − fi⊤ β) − Sh (θ∗ − fi⊤ β∗ ) , β∈B∗ (δ) n

Eγn (δ) ≤ 2E sup

i=1

where ε1 , .., εn are i.i.d. Rademacher random variables. Denote the function class G∗ (δ) := Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) : β ∈ B∗ (δ) . Conditioned on (X1 , Y1 ), . . . , (Xn , Yn ), and for g = Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) ∈ G∗ (δ), we denote g(Xi , Yi ) = Sh (θ∗ − Yi ) Sh (θ∗ − fi⊤ β) − Sh (θ∗ − fi⊤ β∗ ) and define the distance L2 (Pn ) in G∗ (δ) by n

∥g1 − g2 ∥2L2 :=

2 1X g1 (Xi , Yi ) − g2 (Xi , Yi ) , n

∀g1 , g2 ∈ G∗ (δ).

i=1

By Dudley’s entropy bound (Koltchinskii, 2011, Theorem 3.11), we have n

C2 1X εi g(Xi , Yi ) ≤ √ E Eγn (δ) ≤ 2E sup n g∈G∗ (δ) n i=1

Z 2σn q log N G∗ (δ); L2 (Pn ); ϵ dϵ, 0

where C2 > 0 is an absolute constant and σn2 := supg∈G∗ (δ) Pn g 2 ≲ h−2 U 2 δ 2 . Note that N (F; d; ϵ) represents the ϵ-covering number of a set F under the distance d(·, ·). If g1 = Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β1 ) − Sh (θ∗ − f ⊤ β∗ ) and g2 = Sh (θ∗ − Y ) Sh (θ∗ − f ⊤ β2 ) − Sh (θ∗ − f ⊤ β∗ ) , we have an constant C3 > 0 such that σn2 ≤ C3 h−2 U 2 δ 2 with probability at least 1 − n−11 , and 1 ∥g1 − g2 ∥2L2 (Pn ) =

n X i=1

2 C 3 U 2 Sh2 (θ∗ − Yi ) Sh (θ∗ − fi⊤ β1 ) − Sh (θ∗ − fi⊤ β2 ) ≤ ∥β1 − β2 ∥2 , h2

implying that √ hϵ 3 C3 δU K N G∗ (δ); L2 (Pn ); ϵ ≤ N B∗ (δ); ∥ · ∥, √ ≤ , hε C3 U where the last inequality is due to the fact that B∗ (δ) ⊂ RK has a diameter at most 2δ. Plugging them into the Dudley’s entropy bound, we get s √ Z 2√C3 ·U δ/h √ C2 3 C3 δU K log dϵ Eγn (δ) ≤ √ E hϵ n 0 √ Z r 2C2 KC3 U δ 1 3 √ = log dx 2x h n 0 r K Uδ ≲ . n h Therefore, there exists an constant C2 > 0, with probability at least 1 − e−t , r r C2 U δ K t t γn (δ) ≤ + + . h n n n

By setting t = C log(n) for a sufficiently large constant C > 0, we get with probability at least 1 − n−11 , CU δ γn (δ) ≤ h

K + log n , n

which holds for any fixed δ > 0. Following the discretization argument as in the proof of Lemma 7, we can extend it to all δ ∈ (n−1 , δ0 ) and with a union bound of probabilities. Consequently, there exists an event E1 with P(E1 ) ≥ 1 − n−10 , on which r K + log n b C U 2 J1 ≤ γn (∥βbn − β∗ ∥) ≤ ∥βn − β∗ ∥. h n

(56)

Similarly, we bound supβ∈B∗ (δ) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) . First, sup

Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) ≤ 2h−1 U sup ∥β − β∗ ∥ ≤ 2h−1 U δ,

β∈B∗ (δ)

(57)

β∈B∗ (δ)

and 2 sup E Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ )

sup Var Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) ≤ β∈B∗ (δ)

β∈B∗ (δ)

≤ h−1 U 2 sup ∥β − β∗ ∥2 ≲ h−1 U 2 δ 2 .

(58)

β∈B∗ (δ)

These bounds indicate that sup

(Pn − P) Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ )

β∈B∗ (δ)

≤ E sup

(Pn − P) Sh (θ∗ − f β) − Sh (θ∗ − f β∗ ) + 2U δ ⊤

⊤

β∈B∗ (δ)

Then, supβ∈B∗ (δ) Pn Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ )

t t . + nh nh

shares the same order with σn2 in the proof

of J1 , that is h−2 U 2 δ 2 . The same procedure concludes that there exist a constant C2 > 0, with probability that at least 1 − n−10 , C2 U sup (Pn − P) Sh (θ∗ − f β) − Sh (θ∗ − f β∗ ) ≤ h β∈B∗ (∥βbn −β∗ ∥) ⊤

⊤

K + log n b ∥βn − β∗ ∥. n

(59)

While J2 ≤ 2 (Pn − P) Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) , 2C2 U J2 ≤ h

K + log n b ∥βn − β∗ ∥. n

(60)

For J3 , we rewrite it as 2 2 Pn Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) − Pn Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) 2 2 − P Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) + P Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) = (Pn − P) Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) · (Pn + P) Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) − (Pn − P) Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) · (Pn + P) Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) ≤ (Pn − P) 2Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) Pn Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) + P 2Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) (Pn − P) Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) . Since 2Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ ) : β ∈ B∗ (δ0 ) is uniformly bounded and indexed by a finite-dimensional parameter, a standard empirical process bound yields r K + log n ⊤ ⊤ , sup (Pn − P) 2Sh (θ∗ − Y ) − Sh (θ∗ − f β) − Sh (θ∗ − f β∗ ) ≤ C4 n β∈B∗ (δ0 )

(61)

and P 2Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ )

≤ 2,

(62)

where C4 > 0 is a constant. Here, P(·) represents the expectation operator E(·). With the bound P Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) Z 1 θ − f ⊤β ∗ t ≤ sup E h−1 S ′ f ⊤ (β − β∗ )dt h b 0 β∈B∗ (∥βn −β∗ ∥) Z 1 θ − f ⊤β U ∥β − β∗ ∥ ∗ t ≤ sup ES ′ dt h h 0 β∈B∗ (∥βbn −β∗ ∥) ≤

U ∥β − β∗ ∥ · C0 h = U C0 ∥βbn − β∗ ∥, h β∈B∗ (∥βbn −β∗ ∥) sup

(63)

there exists a constant C5 > 0 such that Pn Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) ≤ (Pn − P) Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) + P Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ ) r C2 U K + log n b ≤ · ∥βn − β∗ ∥ + U C0 ∥βbn − β∗ ∥, h n 56

(64)

holding with probability at least 1 − n−11 . Combining (59), (61), (62), and (64), there exist a constant CJ3 > 0 and an event E3 with probability P(E3 ) ≥ 1 − 3n−11 ≥ 1 − n−10 , on which J3 has the bound r r r K + log n K + log n b C2 U K + log n b C U 2 J3 ≤ C4 · · ∥βn − β∗ ∥ + U C0 ∥βbn − β∗ ∥ + 2 · ∥βn − β∗ ∥ n h n h n r C U K + log n b ≤ 3 ∥βn − β∗ ∥, (65) h n for n ≳ K + log n. Combining (56), (60) and (65), conditioned on E1 ∩ E2 ∩ E3 with probability at least 1 − 3n−10 , there exist a constant C such that C Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) ≤ h

K + log n b ∥βn − β∗ ∥, n

since n ≫ log n. Step 3: upper bound for ∥βbn − β∗ ∥. Combining with the locally strong convexity of Q(θ∗ , β) in Assumption 3, we have 1 b 1 c0 ∥βbn − β∗ ∥2 ≤ ∥βn − β∗ ∥2 inf λmin h∇2 Q(θ∗ , β) 2h 2h β∈B(δ0 ) r C K + log n b ≤ Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) ≤ ∥βn − β∗ ∥. h n Overall, there exist a constant C and an event E with P(E) ≥ 1 − n−9 , on which r K + log n ∥βbn − β∗ ∥ ≤ C . n This completes the proof.

(66)

A.8

Proof of Lemma 4

Observe that (θ∗ ) Em b MOE βb n

1 n

n X i=1

N 1 X E Sh (θ∗ − Yi ) − Sh (θ∗ − fi⊤ βbn ) + ESh (θ∗ − e fi⊤ βbn ) − q N i=1

≤ ESh (θ∗ − Y ) − P(Y ≤ θ∗ ) {z } | J1

+ E |

1 n

n X

N 1 X Sh (θ∗ − fi⊤ βbn ) − Sh (θ∗ − fi⊤ β∗ ) − E Sh (θ∗ − e fi⊤ βbn ) − Sh (θ∗ − e fi⊤ β∗ ) . N i=1 i=1 {z }

For the first term, by the smoothness of the probability density fY in a neighborhood of θ∗ , we fix a small δ > 0 and write Z J1 = Sh (θ∗ − y) − I{y ≤ θ∗ } dFY (y) =: T1n + T2n + T3n , R

where Z

T1n :=

Sh (θ∗ − y) − 1 dFY (y),

(−∞,θ∗ −δ]

Sh (θ∗ − y) − I{y ≤ θ∗ } dFY (y),

T2n := (θ −δ,θ∗ +δ]

Z ∗

Sh (θ∗ − y)dFY (y).

T3n := (θ∗ +δ,∞)

For the local term T2n , using the change of variables y = θ∗ + hu, we obtain Z δ/h T2n = h ψ(u) fY (θ∗ + hu)du, ψ(u) := S(−u) − I{u ≤ 0}, −δ/h

−1 where S(u) = 1 + exp(−u) . Since ψ is odd, a Taylor expansion of fY (θ∗ + hu) around θ∗ yields Z Z Z 2 ′ 3 2 T2n = hfY (θ∗ ) ψ(u) du + h fY (θ∗ ) uψ(u) du + O h |u| |ψ(u)| du R R ZR = h2 fY′ (θ∗ ) u S(u) − I{u ≥ 0} du + o(h2 ) = O(h2 ). R

For the left tail term T1n , if y ≤ θ∗ − δ, then δ θ∗ − y ≥ . h h 58

Since S is a sigmoid function, sup

Sh (θ∗ − y) − 1 = 1 −

y≤θ∗ −δ

1 ≤ exp 1 + e−δ/h

−

δ , h

and consequently |T1n | ≤ q exp{−δ/h}. Similarly, for the right tail term T3n follows |T3n | ≤ (1 − q) exp{−δ/h}. Combining the above bounds, we obtain |J1 | ≲ h2 + e−δ/h = O(h2 ).

(67)

For J2 , we write |J2 | as n N i i 1 Xh 1 Xh ⊤b ⊤ E Sh (θ∗ − fi βn ) − Sh (θ∗ − fi β∗ ) − E Sh (θ∗ − e fi⊤ βbn ) − Sh (θ∗ − e fi⊤ β∗ ) n N

≤EL

i=1 n h X

1 n

i=1

N i i 1 Xh Sh (θ∗ − fi⊤ βbn ) − Sh (θ∗ − fi⊤ β∗ ) − EU Sh (θ∗ − e fi⊤ βbn ) − Sh (θ∗ − e fi⊤ β∗ ) N i=1

N i i 1 Xh Sh (θ∗ − fi⊤ β) − Sh (θ∗ − fi⊤ β∗ ) − EU Sh (θ∗ − e fi⊤ β) − Sh (θ∗ − e fi⊤ β∗ ) n N β∈B∗ (∥βbn −β∗ ∥) i=1 i=1 h i ≤E sup (Pn − P) Sh (θ∗ − fi⊤ β) − Sh (θ∗ − fi⊤ β∗ ) .

≤EL

n 1 Xh

sup

β∈B∗ (∥βbn −β∗ ∥)

By (59), we get r C2 U K + log n b |J2 | ≤ E∥βn − β∗ ∥ h n r r i C2 U K + log n h K + log n ≤ C P(E) + diam(B)P(E c ) h n n K + log n ≲ , nh

(68)

where the event E refers to the conclusion of Lemma 3 holds with P(E) ≥ 1 − n−10 . Together with (67) and (68), the bias satisfies K + log n MOE 2 , Em b βb (θ∗ ) = O h + n nh which concludes the proof.

A.9

Proof of Theorem 5

By Lemma 3, βbn is a good estimator of β∗ at θ = θ∗ , we then prove that m b MOE b (θ∗ ) is also a good βn

estimator of m b

MOE β∗

(θ∗ ). 59

m b MOE (θ∗ ) − m b MOE β∗ (θ∗ ) βb n

n N i i 1 Xh 1 Xh = Sh (θ∗ − fi⊤ β∗ ) − Sh (θ∗ − fi⊤ βbn ) − Sh (θ∗ − e fi⊤ β∗ ) − Sh (θ∗ − e fi⊤ βbn ) n N i=1 i=1 h i i h = Pn − P Sh (θ∗ − f ⊤ β∗ ) − Sh (θ∗ − f ⊤ βbn ) − PN − P Sh (θ∗ − f ⊤ β∗ ) − Sh (θ∗ − f ⊤ βbn ) .

Then, by (57), (58), (59) and Bernstein inequality, we have √

nm b MOE (θ∗ ) − m b MOE β∗ (θ∗ ) βb n

√ ≤ n

i h Pn − P Sh (θ∗ − f ⊤ β) − Sh (θ∗ − f ⊤ β∗ )

sup β∈B∗ (∥βbn −β∗ ∥)

√

i h n PN − P Sh (θ∗ − f ⊤ βbn ) − Sh (θ∗ − f ⊤ β∗ )

s r nU 2 ∥βb − β ∥2 √nU ∥βb − β ∥ √ C2 U K + log n n ∗ n ∗ ep ≤ n ∥βbn − β∗ ∥ + O + h n Nh Nh K + log n r K + log n ep √ ≤O + , Nh nh For a fixed β, let 2 b MOE mMOE β (θ∗ ) = E Sh (θ∗ − Y ) − q = O(h ), β (θ∗ ) = Em which is irrelevant to β. We then decompose

√

nm b MOE b (θ∗ ): βn

i i √ h √ h MOE √ √ MOE MOE MOE MOE (θ ) (θ ) + (θ ) + (θ ) − m nm b MOE (θ ) = n m n m b (θ ) − m b n m b ∗ ∗ ∗ ∗ ∗ ∗ β∗ β∗ β∗ β∗ βbn βbn r + log n i √ K + log n √ h MOE ep K √ =O( nh2 ) + O + + n m b β∗ (θ∗ ) − mMOE β∗ (θ∗ ) . Nh nh It suffices to study the last term. Define Zi = Sh (θ∗ − Yi ) − Sh (θ∗ − fi⊤ β∗ ) − ESh (θ∗ − Y ) − ESh (θ∗ − f ⊤ β∗ ) , ei = Sh (θ∗ − e Z fi⊤ β∗ ) − ESh (θ∗ − f ⊤ β∗ ), here, we suppress the dependence on β∗ for simplicity. Then MOE m b MOE β∗ (θ∗ ) − mβ∗ (θ∗ ) =

i=1

1 Xe 1X Zi + Zi , n N

(69)

ei }N are independent with mean zero. Since Sh (·) ∈ (0, 1), we have |Zi | ≤ 2 and {Zi }ni=1 and {Z i=1 and |Zei | ≤ 1. Moreover, Var Zi = Var Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ )

and

ei = Var Sh (θ∗ − e f ⊤ β∗ ) Var Z

are bounded as well. By central limit theorem, we get n

1 X d √ Zi → N 0, Var Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) . n i=1

By Bernstein inequality, we get ! r √ X N n n log n ei = O ep · Var Sh (θ∗ − f ⊤ β∗ ) . Z N N i=1

Therefore, we write √ X N √ MOE n MOE ei , Z n m b β∗ (θ∗ ) − mβ∗ (θ∗ ) = Zn,f∗ + N

(70)

i=1

where Zn,f∗ →d N

0, Var Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ )

as n → ∞.

Overall, combing (67), (69) and (70), we have √

! r K + log n r K + log n n log n ep ep nh2 + O ep √ nm b MOE (θ) = O + Zn,f∗ + O + βbn Nh N nh ! r r √ 2 K + log n K + log n n log n ep nh + √ + + . (71) = Zn,f∗ + O Nh N nh √

ep (·), we get By the Berry-Esseen bound and the high probability bound inherited from O ! r r √n m b MOE (θ∗ ) √ b K + log n K + log n n log n βn sup P p ≤ t − Φ(t) = O nh2 + √ + + , Nh N nh Q(θ∗ , β∗ ) t∈R where Q(θ∗ , β∗ ) = Var Sh (θ∗ − Y ) − Sh (θ∗ − f ⊤ β∗ ) . Therefore, √ ! ! r r nm b MOE (θ∗ ) √ 2 K + log n K + log n n log n βbn p P ≤ zα/2 = 1 − α + O nh + √ + + . Nh N nh Q(θ∗ , β∗ )

(72)

It remains to bound Qn (θ∗ , βbn ) − Q(θ∗ , β∗ ) . By (59), (66) and Lemma 3, Qn (θ∗ , βbn ) − Q(θ∗ , β∗ ) ≤

sup

Qn (θ∗ , β) − Q(θ∗ , β) + Q(θ∗ , βbn ) − Q(θ∗ , β∗ )

β∈B∗ (∥βbn −β∗ ∥)

≤ Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) +

sup

2 (Pn − P) Sh (θ∗ − f ⊤ β∗ ) − Sh (θ∗ − f ⊤ β)

β∈B∗ (∥βbn −β∗ ∥)

(Pn − P) Sh (θ∗ − f ⊤ β∗ ) − Sh (θ∗ − f ⊤ β) · (Pn + P) Sh (θ∗ − f ⊤ β∗ ) − Sh (θ∗ − f ⊤ β)

sup β∈B∗ (∥βbn −β∗ ∥)

≤ Q(θ∗ , βbn ) − Q(θ∗ , β∗ ) + 4

(Pn − P) Sh (θ∗ − f ⊤ β∗ ) − Sh (θ∗ − f ⊤ β)

sup β∈B∗ (∥βbn −β∗ ∥)

ep ≤O

K + log n nh

This bound implies q p Qn (θ∗ , βbn ) − Q(θ∗ , β∗ ) ep K + log n , q =O nh Qn (θ∗ , βbn ) where the second inequality holds assuming Q(θ∗ , β∗ ) > 0 and nh ≳ K + log n. Since

√ √ p p nm b MOE (θ∗ ) nm b MOE (θ∗ ) b Q(θ∗ , β∗ ) − Qn (θ∗ , β∗ ) βn βbn q p · 1+ . = p Q(θ , β ) Q (θ , β ) ∗ ∗ n ∗ ∗ b Qn (θ∗ , βn )

(73)

Applying Bernstein inequality to Zn,f∗ in (71), we get √ nm b MOE (θ∗ ) p βbn ep ( log n), p =O Q(θ∗ , β∗ ) assuming n → ∞, nh2 → ∞, nh4 → 0 and N ≳ n. Together with (73), we get √ √ √ nm b MOE (θ∗ ) nm b MOE (θ∗ ) log n (K + log n) log n βbn βbn e q = p + Op √ + . nh n Q(θ∗ , β∗ ) Qn (θ∗ , βbn ) Plugging into (72), we conclude that ! r r √n m b MOE (θ∗ ) √ 2 K + log n K + log n n log n βbn sup P q ≤ t − Φ(t) = O nh + √ + + , Nh N nh t∈R b Qn (θ∗ , βn ) implying that P(θ∗ ∈ CαMOE ) = P

√

nm b MOE (θ∗ ) βbn q ≤ zα/2 Qn (θ∗ , βbn )

! =1−α+O

which concludes the proof. 62

√

K + log n nh2 + √ + nh

K + log n + Nh

(74)

! n log n , N

A.10

Proof of Theorem 6

By combining (67), (69) and (70) in the proof of Theorem 5, we have √

eN,f + O ep = Zn,f∗ + Z nm b MOE ∗ βb n

√

nh2 +

K + log n √ . nh

Applying the Berry-Esseen bound to both Zn,f∗ and ZeN,f∗ , we get 1 1 and sup F , sup F √ Zn,f∗ (t) − Φ(t) = O √ (t) − Φ(t) = O √ e Z N,f∗ √ n t t N √ Q(θ∗ ,β∗ ) n/N

e ∗ ,β∗ ) Q(θ

FY (·) denotes the c.d.f. of the random variable Y . Moreover, by the convergence rate for the sum of independent non-identically distributed random variables (Petrov, 2012), we get √ nm b MOE (θ∗ ) √ 2 K + log n 1 βbn q √ √ sup P ≤ t − Φ(t) = O nh + + . nh t N e ∗ , β∗ ) Q(θ∗ , β∗ ) + (n/N )Q(θ The rest of the proof follows the same as the proof of Theorem 5.

A.11

Proof of Lemma 5

⊤ Let us derive the explicit form of β∗ and βbn , respectively. Recall f := f (X) = f1 (X), · · · , fK (X) so that Fβ (X) = f ⊤ β. By definition, we write 2 tr WY −Fβ = Σ−2 , E Y − Fβ (X) − X ⊤ δ Fβ XX ⊤ 2 = Σ−2 , E X ⊤ Σ−1 E(Xf ⊤ ) − f ⊤ β XX ⊤ +2E (Y − X ⊤ Σ−1 E(XY )) X ⊤ Σ−1 E(Xf ⊤ ) − f ⊤ βXX ⊤ + const. | {z }| {z } u(Y )

v⊤ (X)

b − H∥ and ∥b Step 1: boundedness of ∥H r − r∥. Denote h i H = E X ⊤ Σ−2 X v(X)v⊤ (X)

and

h i r = −E X ⊤ Σ−2 X u(Y )v(X) .

and their estimators b = 1 H n

n X b −2 Xi )b (Xi⊤ Σ vi (X)b vi⊤ (X)

and b r=−

i=1

1 X ⊤ b −2 (Xi Σ Xi )b ui (Y )b vi (X). n i=1

b = n−1 X⊤ X = n−1 Pn Xi X ⊤ . Similarly, we denote Σ e = N −1 X e ⊤ X. e Note that Σ i i=1 63

Suppose that Assumption 4 holds and n ≳ d log n. Bernstein inequality dictates that there exists an event E1 with P(E1 ) ≥ 1 − n−10 , under which the following bounds hold: r r d log n d log n −1 e ⊤ e −1 ⊤ N X X−Σ =O n X X−Σ =O and n N

(75)

Under event E1 , by Neumann series, we have b −1 − Σ−1 + Σ−1 (Σ b − Σ)Σ−1 − Σ−1 (Σ b − Σ)Σ−1 (Σ b − Σ)Σ−1 = O ∥Σ b − Σ∥3 . Σ b −1 IE c separately, we get b −1 IE and EΣ Therefore, by handling EΣ 1 1 d log n d log n −1 −1 −1 −1 b e EΣ − Σ =O and EΣ − Σ =O , n N

(76)

b = EΣ e = Σ. where we used the fact EΣ Similarly, there is an event E2 with P(E2 ) ≥ 1 − n−10 , under which the following bounds hold r r e ⊤ f∗ (X) e X⊤ f∗ (X) X X⊤ y log n log n −EY X + −EXf∗ (X) = O and −EXf∗ (X) = O . n n n N N b −1 − Σ−1 n−1 X⊤ F + X ⊤ Σ−1 n−1 X⊤ F − EXf ⊤ . Therebi⊤ (X) − vi⊤ (X) = Xi⊤ Σ Note that v i fore, on events E1 ∩ E2 , we have r d log n and ui (Y ) − u bi (Y ) = O n

r bi (X) = O vi (X) − v

Kd log n , n

∀i ∈ [n].

As a result, n 1 X ⊤ b −2 b bi (X)b H−H ≤ Xi Σ Xi v vi⊤ (X) − Xi⊤ Σ−2 Xi vi (X)vi⊤ (X) n i=1

n h i 1X + Xi⊤ Σ−2 Xi vi (X)vi⊤ (X) − E X ⊤ Σ−2 X v(X)v⊤ (X) , n i=1

where the second term is bounded by O

p K 2 log(n)/n on an event E3 satisfying P(E3 ) ≥ 1 − n−10 .

Therefore, on event E1 ∩ E2 ∩ E3 , we get r b − H∥ = O ∥H

dK 2 log n . n

We can bound ∥r − b r∥ in a similar fashion and conclude that there exists an event E4 with P(E4 ) ≥ 1 − n−10 such that the following bound holds r ∥b r − r∥ = O 64

Kd log n , n

under the event E1 ∩ E2 ∩ E3 . 2

b > 0.01c0 . Moreover, there exists an event E0 with P(E0 ) ≥ 1 − e−c1 n/(dK ) such that λmin (H) On the event E0c , we simply have βbn = 0 by the design of Algorithm 4. Step 2: closed-form solution for β∗ and βbn . Since H ⪰ c0 I, we have λmin (H) ≥ c0 . By Weyl’s inequality, b ≥ λmin (H) − ∥H b − H∥ ≥ c0 − ∥H b − H∥. λmin (H) b − H∥ < c0 , we have λmin (H) b > 0 and consequently Thus, for a sufficiently large n such that ∥H b is nonsingular. H By solving ∂tr WY −F β /∂β = 0, we get β∗ = H−1 r. Note that EXu(Y ) = 0 and EXv⊤ (X) = 0. Similarly, by (26), we get n E D −2 1 X c b Yi − Xi⊤ (X⊤ X)−1 X⊤ y + Xi⊤ (X⊤ X)−1 X⊤ F − fi⊤ β Xi Xi⊤ . tr WY −Fβ = Σ , n | {z } | {z } i=1

bi⊤ (X) v

u bi (Y )

b −1b Therefore, we have βbn = H r. b −1b Step 3: bounding βbn − β∗ . Note that βbn − β∗ = H r − H−1 r. To this end, define ui (Y ) = Yi − Xi⊤ Σ−1 E(XY )

and

vi⊤ (X) = Xi⊤ Σ−1 E(Xf ⊤ ) − fi⊤ ,

∀i ∈ [n].

It is then straightforward to show that 1/2

∥βbn − β∗ ∥2 = O

dK 3 log n , n

(77)

if n ≳ dK 2 log n. Note that we used the simple fact that ∥βbn ∥ is uniformly bounded. Indeed, on b > 0.01c0 , and ∥b the event E0 , we have λmin (H) r∥ is uniformly bounded under Assumption 4; on the event E0c , we have βbn = 0 by the design of Algorithm 4. Step 4: bias of PPI-based estimator. The bias of the oracle PPI-based estimator is EθbfPPI − θ∗ ≤ E(X⊤ X)−1 X⊤ y − θ∗ ∗ e ⊤ X) e −1 X e ⊤ f∗ (X) e − E(X⊤ X)−1 X⊤ f∗ (X) + E(X

(78)

Observe that −1 b − Σ−1 n−1 X⊤ y E(X⊤ X)−1 X⊤ y − θ∗ = E Σ −1 −1 b − Σ−1 (EY X) + E Σ b − Σ−1 (n−1 X⊤ y − EY X) ≤ E Σ d log n =O , n where the last inequality is due to (76) and the definition of events E1 ∩ E2 , assuming that N ≫ n ≳ d log n. The second term in (78) can be handled in a similar fashion. As a result, we conclude with d log n PPI b Eθf∗ − θ∗ = O . n

(79)

Step 5: bias of MOE-powered estimator. Recall that θbMOE = θbPPI b and f∗ = Fβ∗ . Therefore, βn

d log n EθbMOE −θ∗ = E θbFPPIb − θbFPPIβ∗ + O βn n h i d log n ⊤ e −1 e ⊤ e ⊤ −1 ⊤ b e = E (X X) X F − (X X) X F (βn − β∗ ) + O n 2 d log n 1/2 ⊤ −1 ⊤ ⊤ −1 ⊤ 1/2 2 e X) e e F e − (X X) X F E ∥βbn − β∗ ∥ + O ≤E (X X . n By the definition of events E1 and E2 , if n ≳ (d + K) log n, we have r K log n 2 1/2 ⊤ −1 ⊤ ⊤ −1 ⊤ e X) e e F e − (X X) X F = O E (X X . n

(80)

Together with (77), we conclude that 2√ (K d + d) log n MOE b Eθ − θ∗ = O , n which concludes the proof.

A.12

Proof of Theorem 7

b −1b r and β∗ = H−1 r. Recall from the proof of Lemma 5 that βbn = H Step 1: connecting θbMOE to θbfPPI . We first show that ∗ √

MOE

n(θb

− θ∗ ) =

√

ep n(θbfPPI − θ∗ ) + O ∗

√dK 2 log n √ n

(81)

where f∗ = Fβ∗ =

k=1 β∗,k fk .

As shown in the Step 3 of the proof of Lemma 5, we have ∥ = ∥θbFPPIb − θbFPPIβ∗ ∥ ≤ ∥θbMOE − θbfPPI ∗ βn

e ⊤ X) e −1 X e ⊤F e − (X⊤ X)−1 X⊤ F · ∥βbn − β∗ ∥ (X

Under the event E1 ∩ E2 ∩ E3 defined in the proof of Lemma 5, we have r r K log n dK 2 log n ⊤ e −1 e ⊤ e ⊤ −1 ⊤ b e (X X) X F − (X X) X F = O and ∥βn − β∗ ∥ = O , n n if n ≳ dK 2 log n, which immediately leads to (81). √

n(θbfPPI − θ∗ ) ∗ to a normal distribution. By decomposing Σ−1 E(XY ) = Σ−1 E(Xf∗ (X)) + Σ−1 E X(Y − f∗ (X)) Step 2: normal approximation of θbfPPI . It suffices to establish the convergence rate of ∗

and the definition of θbfPPI , we get ∗ e ⊤X e −1 X e ⊤ f∗ (X) e − Σ−1 E Xf∗ (X) − θ∗ = X θbfPPI ∗ | {z } J1

+ (X⊤ X)−1 X⊤ y − f∗ (X) − Σ−1 E X(Y − f∗ (X)) . | {z } J2

We first bound the term J1 . Note that J1 =

e e ⊤X N −1 X

−1

− Σ−1

e ⊤ f∗ (X) e − EXf∗ (X) . e ⊤ f∗ (X) e + Σ−1 N −1 X N −1 X

By the definition of the event E1 ∩ E2 ∩ E3 and assuming N ≳ d log n, we get r d log n e J1 = Op , N

(82)

where we used the fact ∥β∗ ∥ = O(1) and |f∗ (X)| = O(1) almost surely. This also implies that ∥δ f∗ ∥ ≲ E X(Y − f∗ (X)) = O(1) under Assumption 4. For J2 term, we write J2 = X⊤ X

−1

X⊤ y − f∗ (X) − Xδ f∗ ) =Σ−1 n−1 X⊤ y − f∗ (X) − Xδ f∗ ) + (X⊤ X/n)−1 − Σ−1 · n−1 X⊤ y − f∗ (X) − Xδ f∗ ) .

By Bernstein inequality, there exists an event E4 with P(E4 ) ≥ 1 − n−10 such that r log n 1 ⊤ X (y − f∗ (X) − Xδ f∗ ) = O , n n 67

(83)

where we used the fact E X(Y − f∗ (X) − X ⊤ δ f∗ ) = 0. Then, on event E1 ∩ E4 , we get s n X √ 1 d log2 n −1 ⊤ e nJ2 = Σ · √ Yi − f∗ (Xi ) − Xi δ f∗ Xi + Op n n i=1

Based on (82) and (83), we get √

− θ∗ ) = Σ−1 · n(θbfPPI ∗

s r n X 1 dn log n d log2 n ⊤ e √ Yi − f∗ (Xi ) − Xi δ f∗ Xi + Op . + N n n i=1

By the multivariate Berry-Esseen theorem (Raič, 2019), we get s r 2 d n log n d2 log2 n − θ ∈ U − P T ∈ U = O , (84) sup P n θbfPPI + ∗ f∗ ∗ N n U where Tf∗ ∼ N 0, Σ−1 E(Y − f∗ (X) − X ⊤ δ f∗ )2 XX ⊤ Σ−1 and U is taken supremum over all √

convex sets in Rd . Note that we also used the Nazarov-type inequality (Nazarov, 2004) concerning the Gaussian surface area of convex sets. Moreover, for a fixed index s ∈ [d], we can also get s √ bPPI r n θf∗ ,s − θ∗,s dn log n d log2 n sup P ≤ t − Φ(t) = O + , V∗,s N n t∈R 2 := e⊤ Σ−1 E(Y − f (X) − X ⊤ δ )2 XX ⊤ Σ−1 e . where V∗,s s ∗ f∗ s √ Step 3: Berry-Esseen bound for n θbMOE − θ∗ . Combining (81) and (85), we get √ 2 √ bMOE r n θs − θ∗,s dn log n dK log n √ , sup P ≤ t − Φ(t) = O + V∗,s N n t∈R

(85)

(86)

for any fixed index s ∈ [d]. Denote

X n 2 1 ⊤b ⊤ b −1 Yi − Fβbn (Xi ) − Xi δ Fβb Xi Xi Σ es , n n i=1 bF = Σ b −1 · n−1 Pn Yi − F b (Xi ) Xi . Under the event E1 ∩ E2 ∩ E3 , where Fβbn (Xi ) = fi⊤ βbn and δ i=1 b βn β 2 b −1 · = e⊤ Vb∗,s sΣ

we have Fβbn (Xi ) − f∗ (Xi ) =

fi⊤ (βbn − β∗ )

√ ≲

r K∥βbn − β∗ ∥ = O

dK 3 log n , n

and as a result n

bF − δ bf δ ∗ b β n

1X Fβbn (Xi ) − f∗ (Xi ) ∥Xi ∥ = O ≲ n i=1

dK 3 log n . n

∀i ∈ [n],

Therefore, on event E1 ∩ E2 ∩ E3 , r bf − δ f + O bF − δ f ≤ δ δ ∗ ∗ ∗ b β n

dK 3 log n n

n 1X · ≤ Σ Yi − f∗ (Xi ) Xi − E Y − f∗ (X) X + O n i=1 −1 b − Σ−1 · E Y − f∗ (X) X + Σ r dK 3 log n e =Op , n

b −1

dK 3 log n n

where the last bound holds due to the definition of the event E1 ∩ E2 ∩ E3 . This also implies that r 2 2 dK 3 log n ⊤b ⊤ =O , ∀i ∈ [n], Yi − Fβbn (Xi ) − Xi δ Fβb − Yi − f∗ (Xi ) − Xi δ f∗ n n on the same event. As a result, on the event E1 ∩ E2 ∩ E3 , we have X r n 2 1 dK 3 log n 2 ⊤ b −1 ⊤ ⊤ −1 b b V∗,s =es Σ · Yi − f∗ (Xi ) − Xi δ f∗ Xi Xi · Σ es + O . n n i=1

By Bernstein inequality (Tropp, 2012), there exists an event E5 with P(E5 ) ≥ 1 − n−10 on which, r n 1X d log n 2 ⊤ ⊤ ⊤ ⊤ Yi − f∗ (Xi ) − Xi δ f∗ Xi Xi − E Y − f∗ (X) − X δ f∗ XX =O . n n i=1

Therefore, on event ∩5k=1 Ek , we have r r 3 log n b∗,s V dK dK 3 log n 2 2 b V∗,s = V∗,s + O implying that =1+O . n V∗,s n Following the same argument as that in the proof of Theorem 3, we conclude that √ r √ bMOE n θs − θ∗,s dn log n dK 3 log n √ + . ≤ t − Φ(t) = O sup P N n t∈R Vb∗,s The rest of the proof is straightforward and hence omitted.

A.13

Proof of Lemma 6

Step 1: formulas for β∗ and βbn . Let us derive the explicit form and difference bound of β∗ and βbn . Note that tr WY −Fβ = tr E (f ⊤ β − Y )2 XX ⊤ − E (f ⊤ β − Y )X E (f ⊤ β − Y )X ⊤ = E (f ⊤ β − Y )2 X ⊤ X − E (f ⊤ β − Y )X ⊤ E (f ⊤ β − Y )X = β ⊤ Ef X ⊤ Xf ⊤ − (Ef X ⊤ )(EXf ⊤ ) β − 2 EY X ⊤ Xf ⊤ − (EY X ⊤ )(EXf ⊤ ) β + constant; {z } | {z } | r

n n n X X ⊤ ⊤ X ⊤ −1 2 ⊤ −2 ⊤ c tr WY −Fβ = tr n (fi β − Yi ) Xi Xi − n (fi β − Yi )Xi (fi β − Yi )Xi i=1

= n−1

n X

i=1

⊤ (fi β − Yi )2 Xi⊤ Xi − n−2

i=1

n X

⊤ (fi β − Yi )Xi⊤

i=1

i=1 n X

⊤ (fi β − Yi )Xi

i=1

n n n X X X ⊤ −1 ⊤ ⊤ −2 ⊤ ⊤ =β n fi Xi Xi fi − n fi Xi Xi fi β i=1

i=1

}

b H

n n n X X X −1 ⊤ ⊤ −2 ⊤ ⊤ −2 n Yi Xi Xi fi − n Yi Xi Xi fi β + constant. i=1

i=1

}

b r

Under Assumption 5, H is invertible so that β∗ = H−1 r exists and is unique. By Bernstein Inequality, there exist an event E1 with P(E1 ) ≥ 1 − n−11 such that max

n−1 n

−1

n X i=1 n X

fi Xi⊤ Xi fi⊤ − Ef X ⊤ Xf ⊤ , n−1 Yi Xi⊤ Xi fi⊤ − EY X ⊤ Xf ⊤

, n

n X

fi Xi⊤ − Ef X ⊤ ,

i=1 n X −1

i=1

Yi Xi⊤ − EY X ⊤ EXf ⊤

i=1

(d + K) log n . n

Conditioned on E1 , we have max

b −H , b H r−r

r =O

(d + K) log n . n

b − H ≤ c0 /2 on the event E1 . Then If n ≥ C1 (d + K) log n for a large enough C1 > 0, we have H Weyl’s inequality yields b ≥ λmin (H) − H b − H ≥ c0 , λmin (H) 2 b −1b r. implying that βbn = H b ≥ 0.01c0 . Moreover, there exists an event E0 with P(E0 ) ≥ 1 − ec1 n/(d+K) such that λmin (H) On the event E0 , ∥βbn ∥ is uniformly bounded under Assumption 5; on the event E0c , βbn = 0 by Algorithm 5. Therefore, we conclude that ∥βbn ∥ is uniformly bounded. Step 2: upper bound for ∥βbn − β∗ ∥. Recall that on event E1 , β∗ = H−1 r

b −1b and βbn = H r.

(87)

Then, b −1 b −1b r − H−1 r ≤ H βbn − β∗ = H

b −H H

r∥ + ∥H−1 ∥ · b r−r . H−1 · ∥b

Under Assumption 5, we have H−1 ≤ c−1 0 . Putting together the bounds above, we have conditional on E1 that r βbn − β∗ = O

(d + K) log n . n

(88)

b MOE Step 3: upper bound for Em b (θ) − m(θ)) . For any fixed θ ∈ Θ, by definition, we get βn

b MOE (θ) − EX(Y − S(X ⊤ θ)) Em βbn e⊤ F e βbn − S(Xθ) e = N −1 EX − n−1 EX⊤ Fβbn − y) h i e ⊤F e − n−1 X⊤ F (βbn − β∗ ) + N −1 EX e ⊤F e − n−1 EX⊤ F β∗ ≤ E N −1 X e X e ⊤ θ) − EXS(X ⊤ θ)) . + n−1 EX⊤ y − EXY + N −1 EXS(

(89)

(90)

e ⊤F e = n−1 EX⊤ F = EXf ⊤ , n−1 EX⊤ y = EXY and N −1 EXS( e X e ⊤ θ) = EXS(X ⊤ θ), Since N −1 EX the latter three terms of (90) vanish. It suffices to bound the first term. h i e ⊤F e − n−1 X⊤ F (βbn − β∗ ) E N −1 X e ⊤F e − n−1 X⊤ F (βbn − β∗ ) ≤ E N −1 X h i e ⊤F e − EXf · ∥βbn − β∗ ∥ ≤ E n−1 X⊤ F − EXf + N −1 X By Bernstein inequality, there exists an event E2 with P(E2 ) ≥ 1 − n−11 such that r r (d + K) log n (d + K) log n −1 e ⊤ e −1 ⊤ n X F − EXf = O and N X F − EXf = O . n N Combining with (88), we then have b MOE Em (θ) − EX(Y − S(θ⊤ X)) βbn h i e ⊤F e − EXf · ∥βbn − β∗ ∥IE ∩E ≤ E n−1 X⊤ F − EXf + N −1 X 1 2 h i e ⊤F e − EXf · ∥βbn − β∗ ∥IE c ∪E c + E n−1 X⊤ F − EXf + N −1 X 1 2 (d + K) log n =O , n where the last inequality holds due to the boundedness of X, f , and β. Note that we need an upper bound on E1/2 ∥βbn − β∗ ∥2 , which was established at the end of Step 1. 71

A.14

Proof of Theorem 8

b MOE Step 1: normal approximation of m b (θ). By definition, for any θ ∈ Θ, we have βn

√

√ √ MOE PPI PPI b MOE b b b (θ) + (θ) − m(θ) , (θ) − m(θ) = m (θ) − m m n m n n f f b b ∗ ∗ β β n

(91)

⊤ b PPI where f∗ = Fβ∗ and m(θ) = Em f∗ (θ) = EX(Y − S(X θ)).

Observe that e ⊤F e − n−1 X⊤ F βbn − β∗ N −1 X

b PPI b MOE (θ) − m m f∗ (θ) = βb n

≤

e ⊤F e − EXf + n−1 X⊤ F − EXf N −1 X

βbn − β∗ .

On the event E1 ∩ E2 defined in the proof of Lemma 6, we have √ (d + K) log n MOE PPI √ b βb (θ) − m b f∗ (θ) = O n m . n n

(92)

For second term in the RHS of (91) can be written as n i 1 Xh b (θ) − m(θ) = − √ Xi Fβ∗ (Xi ) − Yi − EX(Fβ∗ (X) − Y ) n m n i=1 √ X N h i n ei ) − S X e ⊤ θ − EX(Fβ (X) − S(X ⊤ θ)) ei Fβ (X + X i ∗ ∗ N i=1 r ne = −Zn,f∗ + ZN,f∗ , N h i e N,f := N −1/2 PN X ei Fβ (X ei ) − S X e ⊤ θ − EX Fβ (X) − S(X ⊤ θ) . where Z ∗ ∗ ∗ i i=1

√

PPI f∗

(93)

Since ∥X∥, |Y | and ∥f ∥ are bounded under Assumption 5, by Bernstein inequality, there exists an event E3 with P(E3 ) ≥ 1 − n−10 such that the following equality holds r (d + K)n log n √ PPI b β∗ (θ) − m(θ)) = −Zn,f∗ + O n m N on the event E1 ∩ E2 ∩ E3 . This proves the first claim that √

K) log n ep (d + √ b βb (θ) − m(θ)) = −Zn,f∗ + O n m + n n MOE

(d + K)n log n . N

(94)

Step 2: coverage probability. We continue from (94). Note that m(θ∗ ) = 0. By the Berry-Esseen ep (·), we get bound and the high probability bound inherited from O r MOE √n e⊤ m s b βbn (θ∗ ) (d + K) log n (d + K)n log n √ sup P q ≤ t − Φ(t) = O + , ∀s ∈ [d]. N n t∈R e⊤ W e s

Y −Fβ∗ s

(95)

Therefore, r MOE √n e⊤ m s b βbn (θ∗ ) (d + K)n log n (d + K) log n √ + , ∀s ∈ [d]. P q ≤ zα/(2d) = 1 − α/d + O N n e⊤ s WY −Fβ∗ es c Y −F from the c Y −F − WY −F . By the derivation of W Then it’s sufficient to bound W β β∗ b β n

proof of Lemma 6, we have c Y −F − WY −F W β∗ b β n

c Y −F − WY −F c Y −F − W c Y −F + W ≤ W β∗ β∗ β∗ bn β ⊤ ≤ βbn + β∗ SCov f ⊤ X, f ⊤ X · βbn − β∗ + SCov Y, f ⊤ X · βbn − β∗ + n

−1

n X

Fβ∗ (Xi ) − Yi Xi − E(Fβ∗ (X) − Y )X · n

i=1

+ n−1

−1

n X

Fβ∗ (Xi ) − Yi Xi + E(Fβ∗ (X) − Y )X

i=1

n X

2 fi⊤ β∗ − Yi Xi Xi⊤ − E(f ⊤ β∗ − Y )2 XX ⊤

i=1

By Bernstein inequality, there exists an event E4 with P(E4 ) ≥ 1 − n−10 such that ( n X −1 max n fi⊤ β∗ − Yi Xi − E(f ⊤ β∗ − Y )X , n

−1

i=1 n X

) 2 fi⊤ β∗ − Yi Xi Xi⊤ − E(f ⊤ β∗ − Y )2 XX ⊤

r d log n

i=1

Due to the bound (88) and the uniform boundedness of ∥X∥, |Y | and |fk (X)|, conditional on T4 i=1 Ei , we get r (d + K) log n c WY −Fβb − WY −Fβ∗ = O . n n Observe that √

b MOE n e⊤ sm b (θ∗ )

√

! c b MOE n e⊤ (θ∗ ) e⊤ sm s WY −Fβbn − WY −Fβ∗ es βbn q q q =q 1+ q ⊤ c c e⊤ e⊤ e⊤ e⊤ e⊤ s WY −Fβ∗ es s WY −Fβb es s WY −Fβ∗ es s WY −Fβb es + s WY −Fβ∗ es βn

By the normal approximation in (94) and the boundedness of covariance matrix, there exists an event E5 with P(E5 ) ≥ 1 − n−10 such that √ ⊤ MOE √ ⊤ MOE r b b (θ∗ ) b b (θ∗ ) n es m n es m d+K βn βn q =q +O log n n ⊤W c Y −F es e e e⊤ W Y −Fβ∗ s s s b β n

when n ≳ (d + K)2 log2 n and N ≳ (d + K)n log n. Together with (95), we conclude that r MOE √n e⊤ m s b βbn (θ∗ ) (d + K) log n (d + K)n log n √ , sup P q ≤ t − Φ(t) = O + N n c Y −F es t∈R e⊤ W s b β

(96)

implying that ( √

)! b MOE (θ∗ ) n e⊤ sm βbn q ≤ zα/(2d) c e⊤ s=1 s WY −Fβbn es √ ⊤ MOE ! d b b (θ∗ ) n es m Y βn ≥1− P q > zα/(2d) c e⊤ s=1 s WY −Fβbn es r d2 (d + K)n log n d(d + K) log n √ + , =1−α+O N n

P θ∗ ∈ CαMOE = P

d \

which concludes the proof.

A.15

Proof of Theorem 9

For any fixed θ ∈ Θ, by (92) and (93) in the Proof of Theorem 8, we get r √ ne K) log n MOE ep (d + √ b βb (θ) − m(θ) = −Zn,f∗ + n m ZN,f∗ + O . n N n e N,f , we get Applying the Berry-Esseen bound to both Zn,f∗ and Z ∗ r e⊤ d+K s Zn,f∗ sup P p ≤ t − Φ(t) = O n e⊤ t∈R s WY −f∗ es r e e⊤ d+K s ZN,f∗ ≤ t − Φ(t) = O , sup P ⊤ N es Wf∗ −X ⊤ θ es t∈R for any s ∈ [d]. By the convergence rate for the sum of independent non-identically distributed random variables, we get, for each s ∈ [d], that √ r b MOE n m (θ) (d + K) log n d+K s βbn q √ sup P ≤ t − Φ(t) = O + . N n ⊤ t e⊤ s WY −f∗ es + (n/N )es Wf∗ −X ⊤ θ es The rest of the proof follows the same as the proof of Theorem 8.

Experiments Result

B.1

Result of coverage and interval width

Table 4 reports the empirical coverage, interval width, width ratio, and coverage-agreement code for all methods across tasks and settings. Overall, for mean and quantile inference, the proposed PPIMOE method achieves coverage close to the nominal 95% level while producing intervals that are substantially shorter than those of the conventional estimator. This gain is especially pronounced in both linear and nonlinear settings, where the width ratio of PPI-MOE is often around 0.2–0.35 relative to the conventional baseline, indicating a large efficiency improvement without a noticeable loss in coverage. For linear regression inference with adequate samples n = 500, all methods behave reasonably well, and PPI-MOE remains competitive in terms of both coverage and interval width. However, when the samples size is limited to 200, coverage deteriorates substantially for all methods, including the conventional estimator, with several entries falling well below the nominal level and receiving weak or empty coverage-agreement codes. This suggests that the issue is not specific to one particular prediction-powered construction, but rather reflects a broader instability under small samples and high dimensions. For logistic regression inference, the results are more variable. In some settings, especially under the nonlinear design, PPI-MOE still attains near-nominal coverage with noticeably shorter intervals than the conventional method. In contrast, some alternative PPI procedures yield extremely large interval widths in the linear setting, suggesting numerical instability or near-separation phenomena in the underlying logistic regression fit. Relative to these unstable alternatives, PPI-MOE appears considerably more stable and practically usable. Taken together, the results indicate that PPI-MOE delivers the clearest gains for mean and quantile inference, where it consistently maintains coverage while substantially reducing interval width. For regression coefficient inference, the benefit of MOE-powered prediction-powered methods depends more strongly on whether the inferential model is well specified. Nevertheless, PPI-MOE still adaptively chooses a suitable weight, moving toward the variance-minimizing choice when prediction-powered correction is beneficial and degenerating to the conventional estimator when it is not.

Table 4: Full Coverage and Width results across tasks. The width ratio is computed relative to the Conventional estimator within the same task, setting, and sample size. Column Code is defined in Table 3.

Task

Data Mode

Method

Cov.

Width

Ratio

Code

Mean inference Mean

Linear

200

Conventional

0.934

45.9838

1.000

Mean

Linear

200

PPI-best

0.956

9.6054

0.209

***

Mean

Linear

200

PPI-mean

0.953

11.5188

0.250

***

Mean

Linear

200

PPI-worst

0.954

16.7861

0.365

***

Mean

Linear

200

PPI-MOE

0.952

9.4184

0.205

***

Mean

Linear

500

Conventional

0.960

29.1178

1.000

Mean

Linear

500

PPI-best

0.946

6.0713

0.209

***

Mean

Linear

500

PPI-mean

0.944

7.2855

0.250

Mean

Linear

500

PPI-worst

0.944

10.6361

0.365

Mean

Linear

500

PPI-MOE

0.940

5.9590

0.205

Mean

Nonlinear

200

Conventional

0.940

55.1793

1.000

Mean

Nonlinear

200

PPI-best

0.954

11.9046

0.216

***

Mean

Nonlinear

200

PPI-mean

0.949

14.3381

0.260

***

Mean

Nonlinear

200

PPI-worst

0.954

20.8877

0.379

***

Mean

Nonlinear

200

PPI-MOE

0.944

11.5958

0.210

Mean

Nonlinear

500

Conventional

0.934

34.8191

1.000

Mean

Nonlinear

500

PPI-best

0.958

7.5351

0.216

***

Mean

Nonlinear

500

PPI-mean

0.948

9.0691

0.260

***

Mean

Nonlinear

500

PPI-worst

0.934

13.1799

0.379

Mean

Nonlinear

500

PPI-MOE

0.940

7.3539

0.211

Quantile inference (q = 0.5) Quantile

Linear

200

Conventional

0.954

61.5012

1.000

***

Quantile

Linear

200

PPI-best

0.954

18.4687

0.300

***

Quantile

Linear

200

PPI-mean

0.955

23.9068

0.389

***

Quantile

Linear

200

PPI-worst

0.954

34.3709

0.559

***

Continued on next page

Table 4: Full Coverage and Width results across tasks (continued).

Task

Data Mode

Method

Cov.

Width

Ratio

Code

Quantile

Linear

200

PPI-MOE

0.954

18.7296

0.305

***

Quantile

Linear

500

Conventional

0.936

37.4121

1.000

Quantile

Linear

500

PPI-best

0.962

11.8119

0.316

Quantile

Linear

500

PPI-mean

0.949

15.1306

0.404

***

Quantile

Linear

500

PPI-worst

0.950

21.6758

0.579

***

Quantile

Linear

500

PPI-MOE

0.950

11.8886

0.318

***

Quantile

Nonlinear

200

Conventional

0.938

70.3177

1.000

Quantile

Nonlinear

200

PPI-best

0.950

24.2521

0.345

***

Quantile

Nonlinear

200

PPI-mean

0.940

30.2135

0.430

Quantile

Nonlinear

200

PPI-worst

0.926

42.4006

0.603

Quantile

Nonlinear

200

PPI-MOE

0.942

23.7573

0.338

***

Quantile

Nonlinear

500

Conventional

0.934

43.1535

1.000

Quantile

Nonlinear

500

PPI-best

0.952

15.2533

0.353

***

Quantile

Nonlinear

500

PPI-mean

0.938

19.1010

0.443

Quantile

Nonlinear

500

PPI-worst

0.922

26.7733

0.620

Quantile

Nonlinear

500

PPI-MOE

0.950

14.9006

0.345

***

Linear regression coefficient (βb1 ) Lin. Reg.

Linear

200

Conventional

0.920

2.7534

1.000

Lin. Reg.

Linear

200

PPI-best

0.920

2.7535

1.000

Lin. Reg.

Linear

200

PPI-mean

0.918

5.3552

1.945

Lin. Reg.

Linear

200

PPI-worst

0.910

9.0615

3.290

Lin. Reg.

Linear

200

PPI-MOE

0.910

2.7072

0.983

Lin. Reg.

Linear

500

Conventional

0.945

1.7679

1.000

***

Lin. Reg.

Linear

500

PPI-best

0.945

1.7679

1.000

***

Lin. Reg.

Linear

500

PPI-mean

0.933

3.4479

1.950

Lin. Reg.

Linear

500

PPI-worst

0.940

5.8670

3.318

Lin. Reg.

Linear

500

PPI-MOE

0.930

1.7563

0.993

Lin. Reg.

Nonlinear

200

Conventional

0.905

7.1322

1.000

Lin. Reg.

Nonlinear

200

PPI-best

0.920

4.7344

0.664

Continued on next page

Table 4: Full Coverage and Width results across tasks (continued).

Task

Data Mode

Method

Cov.

Width

Ratio

Code

Lin. Reg.

Nonlinear

200

PPI-mean

0.914

8.4522

1.185

Lin. Reg.

Nonlinear

200

PPI-worst

0.885

14.4147

2.021

Lin. Reg.

Nonlinear

200

PPI-MOE

0.910

4.5400

0.636

Lin. Reg.

Nonlinear

500

Conventional

0.970

4.5206

1.000

Lin. Reg.

Nonlinear

500

PPI-best

0.945

2.9830

0.660

***

Lin. Reg.

Nonlinear

500

PPI-mean

0.948

5.4124

1.197

***

Lin. Reg.

Nonlinear

500

PPI-worst

0.940

9.2308

2.042

Lin. Reg.

Nonlinear

500

PPI-MOE

0.950

2.8814

0.637

***

Logistic regression coefficient (βb1 ) Log. Reg.

Linear

500

Conventional

0.920

3.8692

1.000

Log. Reg.

Linear

500

PPI-best

0.985

13.3119

3.441

Log. Reg.

Linear

500

PPI-mean

0.990

102660.4034

26533.991

Log. Reg.

Linear

500

PPI-worst

0.990

325903.3788

84231.167

Log. Reg.

Linear

500

PPI-MOE

0.970

6.6481

1.718

Log. Reg.

Linear

1000

Conventional

0.885

2.1308

1.000

Log. Reg.

Linear

1000

PPI-best

0.930

2.6450

1.241

Log. Reg.

Linear

1000

PPI-mean

0.952

7588.0250

3561.092

***

Log. Reg.

Linear

1000

PPI-worst

0.975

45510.8492

21357.191

Log. Reg.

Linear

1000

PPI-MOE

0.915

2.9004

1.361

Log. Reg.

Nonlinear

500

Conventional

0.945

0.4654

1.000

***

Log. Reg.

Nonlinear

500

PPI-best

0.940

0.3535

0.760

***

Log. Reg.

Nonlinear

500

PPI-mean

0.951

0.4392

0.944

***

Log. Reg.

Nonlinear

500

PPI-worst

0.975

0.6089

1.308

Log. Reg.

Nonlinear

500

PPI-MOE

0.960

0.2989

0.642

***

Log. Reg.

Nonlinear

1000

Conventional

0.935

0.3218

1.000

Log. Reg.

Nonlinear

1000

PPI-best

0.950

0.2401

0.746

***

Log. Reg.

Nonlinear

1000

PPI-mean

0.944

0.2973

0.924

***

Log. Reg.

Nonlinear

1000

PPI-worst

0.925

0.4074

1.266

Log. Reg.

Nonlinear

1000

PPI-MOE

0.940

0.2068

0.643

***

Record · ID 146949 · SHA-256 5c496da12192c422
Retrieved via Conceptio — every document is proof-bundled with source, license, and retrieval metadata.