Optimal Model Averaging for Semiparametric Partially Linear Models with Censored Data

Abstract

In the past few decades, model averaging has received extensive attention, and has been regarded as a feasible alternative to model selection. However, this work is mainly based on parametric model framework and complete dataset. This paper develops a frequentist model-averaging estimation for semiparametric partially linear models with censored responses. The nonparametric function is approximated by B-spline, and the weights in model-averaging estimator are picked up via minimizing a leave-one-out cross-validation criterion. The resulting model-averaging estimator is proved to be asymptotically optimal in the sense of achieving the lowest possible squared error. A simulation study demonstrates that the method in this paper is superior to traditional model-selection and model-averaging methods. Finally, as an illustration, the proposed procedure is further applied to analyze two real datasets.

1. Introduction

The semiparametric partially linear model (PLM), which was proposed by [1], has attracted extensive attention in statistics because it combines the interpretability of the linear model with the flexibility of the nonparametric model. A large collection of literature has explored the estimation methods of this model, including the parametric part and the nonparametric function, such as [2,3,4,5], and so on. The premise of using methods in the aforementioned literature is that the model is correctly specified. However, in real data analysis, researchers are always able to collect a variety of variables and are not sure which variables are to be included in the true model. This kind of uncertainty is generally referred to as model uncertainty, which brings great trouble to statistical analysis.

Refs. [6,7,8] pointed out that model selection and model averaging are two mainstream methods to deal with model uncertainty. Model selection, which has a long history, selects a model from a series of candidate models through selection criteria, for example, Akaike information criterion (AIC [9]), Bayesian information criterion (BIC [10]), focused information criterion (FIC [11]), and so on. In addition, shrinkage-estimation-based variable selection was also applied to determine which variables are neededto build a PLM (see, e.g., [12,13,14], among others). These model-selection methods can be viewed as combining a series of candidate models, and assigning a weight of 1 to the selected model and 0 to other models.

As an important alternative to model selection, model averaging incorporates the model uncertainty in statistical analysis by assigning nonzero weight vector to a set of candidate models, which frequently leads to more effective results (see [15]). Bayesian model averaging, which is an important branch of model averaging, has been fully developed in the past decades; see [16] for details. In the current paper, we focus on model-averaging approach in PLM from a frequentist perspective. Since ref. [17] pioneered the use of Mallows criterion for weight choice in model averaging, there is rapidly growing research on asymptotically optimal model averaging. Different kinds of optimal model-averaging methods were proposed, including jackknife model averaging (JMA [7]), Kullback–Leibler model averaging [18], generalized least-squares model averaging [19], leave-subject-out cross-validation [20], K-fold cross-validation [21], and so on. In addition, the optimal model-averaging methods were extended to quantile regression [22], semiparametric models [23,24], missing data [25,26], functional data [27], measurement error data [28], and high-dimensional data [29,30].

Censored data are ubiquitous in the fields of biomedicine, industry, econometrics, etc. For example, in biomedicine, when some sampled individuals are lost to follow-up before the end of the study or drop out during the study, the survival time will be subject to censoring. Compared with the complete data, censoring makes some data unable to be observed completely, which increases the difficulty of statistical analysis. Although an extensive body of literature focusing on the estimation methods in the presence of censored data has been developed, such as [31,32], there is little work on the development of model-averaging approaches with censored data. Based on FIC advocated in [11,33,34,35], theydeveloped model-averaging methods for different regression models with censored data under the local misspecification framework, where the weights of the model-averaging estimators were constructed by the information criterion values, rather than being selected in a data-driven fashion. Moreover, the framework mentioned above requires that the distance between the candidate model and the true model is $O (1/\sqrt{n})$. This means the candidate model is close to the true model when the sample size is large, which is unrealistic.

Without local misspecification framework, ref. [36] constructed the optimal model-averaging estimator for a high-dimensional linear model with censored data by adapting a leave-one-out cross-validation criterion. Ref. [37] studied Mallows model-averaging method for linear models with censored responses, and the resulting model-averaging estimator was proved to be asymptotically optimal in terms of minimizing the squared error loss. The optimal model-averaging methods for censored data mentioned above are based on classical linear models. The primary object of the current paper is to construct an optimal model-averaging estimator for semiparametric PLM with censored responses, in which the weight vector is selected by minimizing a leave-one-out cross-validation criterion. Compared with [36,37], we confront two major challenges. Firstly, the nonparametric function in PLM significantly complicates the construction of the model-averaging estimator and the development of the weight choice criterion. Secondly, our proof of optimality cannot follow the approach of the linear model, since their proof techniques cannot be directly applied when a nonparametric part is present.

The plan of this article is as follows. In Section 2, we describe the model set and introduce the parametric estimation method in the candidate PLM. Section 3 constructs the model-averaging estimator and proposes a weight choice criterion. Section 4 establishes the asymptotic optimality of the model-averaging estimator. Section 5 explores the finite sample performance of our method by the simulation study. Section 6 applies the proposed method to a real dataset. Section 7 gives some conclusions. Proofs are listed in the Appendix A.

2. Model Setup and Parametric Estimation

To facilitate presentation, we first list the basic notations used in this paper in Table 1. Then, we consider the following PLM:

(1)$Y_i={\mu }_i+{ϵ}_i=\sum _{j=1}^{\infty }{x}_{ij}{\beta }_j+g\left(U_i\right)+{ϵ}_i,i=1,\dots ,n,$

where $Y_i$ is a response variable with a continuous distribution function $F(·)$, a countably infinite vector $X_i={(x_{i1},x_{i2},\dots )}^{\top }$ is linearly related to $Y_i$, $U_i$ is a covariate nonlinearly related to $Y_i$, $g(·)$ is an unknown smooth function, and ${ϵ}_i$ is the model error with $E\left({ϵ}_i\right|X_i,U_i)=0$ and $E\left({ϵ}_i^2\right|X_i,U_i)={\sigma }_i^2$. The covariate $U_i$ is distributed on a compact interval $[a^{★},b^{★}]$. Without loss of generality, we take $[a^{★},b^{★}]=[0,1]$. The conditional expectation of the response is denoted as $E\left(Y_i\right|X_i,U_i)={\mu }_i={\sum }_{j=1}^{\infty }{x}_{ij}{\beta }_j+g\left(U_i\right)$.

In survival analysis, we assume $Y_i$ to be a known monotonic transformation of the survival time $T_i$, for example, the commonly used logarithm $Y_i=log{T}_i$. $Y_i$ may be censored by a censoring time $C_i$ and hence cannot be observed completely. We consider a sample of independent observations $(Z_i,{\delta }_i,X_i,U_i)$ for $i=1,\dots ,n$, where $Z_i=min(Y_i,C_i)$ and ${\delta }_i=I(Y_i\le C_i)$ is the censoring indicator.

Let $G(·)$ be the cumulative distribution function of the censoring time, and $Z_{G,i}=Z_i{\delta }_i/\{1-G\left(Z_i\right)\}$ be a synthetic response. Then, from [38], it is not difficult to verify that $E\left(Z_{G,i}\right|X_i,U_i)=E\left(Y_i\right|X_i,U_i)={\mu }_i$. Therefore, under model (1), we obtain

(2)$Z_{G,i}={\mu }_i+e_{G,i}=\sum _{j=1}^{\infty }{x}_{ij}{\beta }_j+g\left(U_i\right)+e_{G,i},i=1,\dots ,n,$

where $E\left(e_{G,i}\right|X_i,U_i)=0$ and $E\left(e_{G,i}^2\right|X_i,U_i)={\sigma }_{G,i}^2$. Model (2) can be expressed in matrix form as

(3)$Z_G=\mu +e_G=X\beta +M+e_G,$

where $Z_G={(Z_{G,1},\dots ,Z_{G,n})}^{\top }$ is an n-dimensional synthetic response vector, conditional mean $\mu ={({\mu }_1,\dots ,{\mu }_n)}^{\top }$ is an n-dimensional vector, $M={\{g\left(U_1\right),\dots ,g\left(U_n\right)\}}^{\top }$, $X={(X_1,\dots ,X_n)}^{\top }$ is a linear covariate matrix, $U={(U_1,\dots ,U_n)}^{\top }$, and $e_G={(e_{G,1},\dots ,e_{G,n})}^{\top }$ is an n-dimensional error vector satisfying $E\left(e_G\right|X,U)=0$ and $E\left(e_G{e}_G^{\top }\right|X,U)={\Omega }_G=\mathrm{diag}({\sigma }_{G,1}^2,\dots ,{\sigma }_{G,n}^2)$.

Assume that we have a total of S candidate PLMs to approximate the true data generating process, where S is allowed to go to infinity. Suppose the sth candidate model is

(4)$Z_G={\mu }_{\left(s\right)}+e_{G,\left(s\right)}=X_{\left(s\right)}{\beta }_{\left(s\right)}+M_{\left(s\right)}+e_{G,\left(s\right)},$

where $X_{\left(s\right)}={(X_{\left(s\right),1},\dots ,X_{\left(s\right),n})}^{\top }$ is an $n\times p_s$ covariate matrix which includes $p_s$ columns of X with full column rank, ${\beta }_{\left(s\right)}={({\beta }_{\left(s\right)1},\dots ,{\beta }_{\left(s\right)p_s})}^{\top }$ is the corresponding $p_s\times 1$ unknown linear regression coefficient vector, $M_{\left(s\right)}={\{g_{\left(s\right)}\left(U_1\right),\dots ,g_{\left(s\right)}\left(U_n\right)\}}^{\top }$ is an $n\times 1$ unknown nonparametric function vector, and $e_{G,\left(s\right)}$ is the model error.

To obtain the estimator of $\mu$ under the sth candidate model, we should estimate the coefficient vector ${\beta }_{\left(s\right)}$ and nonparametric function vector $M_{\left(s\right)}$ firstly. There are many estimation methods for model (4), including kernel smoothing, polynomial spline smoothing, and so on. Recently, ref. [39] pointed out that using B-splines to approximate nonparametric functions in the area of model averaging has great advantages; therefore, in this paper, we adopt the spline technique to estimate the unknowns in model (4).

Denote ${\psi }_n=\{0=u_0<u_1<\dots <u_{J_n}<u_{J_n+1}=1\}$ as a partition of $[0,1]$, where $J_n$ is the number of interior knots. Let ${\zeta }_n$ be the polynomial spline space on interval $[0,1]$ of degree r. From [40], the nonparametric function in PLM can be well estimated by a B-spline expansion. Then, for $g_s\left(u\right)$, one can write

(5)$g_{\left(s\right)}\left(u\right)\approx B_{\left(s\right)}^{\top }\left(u\right){\alpha }_{\left(s\right)},$

where $B_{\left(s\right)}(·)={(B_{\left(s\right)1}(·),\dots ,B_{\left(s\right)k_n}(·))}^{\top }$ is the normalized B-spline basis function vector in the sth candidate model, $k_n=J_n+r+1$, and ${\alpha }_{\left(s\right)}={({\alpha }_{\left(s\right)1},\dots ,{\alpha }_{\left(s\right)k_n})}^{\top }$ is the vector of spline coefficient. Define the $n\times k_n$ matrix $B_{\left(s\right)}={(B_{\left(s\right)}\left(U_1\right),\dots ,B_{\left(s\right)}\left(U_n\right))}^{\top }$. Therefore, there exists a design matrix $X_{\left(s\right)}^{*}=(X_{\left(s\right)},B_{\left(s\right)})$ and the corresponding unknown parameter vector $\gamma ={({\beta }^{\top },{\alpha }^{\top })}^{\top }$ such that

(6)${\mu }_{\left(s\right)}\approx X_{\left(s\right)}^{*}\gamma ,$

where $X_{\left(s\right)}^{*}$ is supposed to be of full-column rank. By regressing $Z_G$ on $X_{\left(s\right)}^{*}$, we can obtain the least-squares estimators of ${\beta }_{\left(s\right)}$ and ${\alpha }_{\left(s\right)}$:

(7)${\widehat{\beta }}_{G,\left(s\right)}={\left\{X_{\left(s\right)}^{\top }(I-Q_{\left(s\right)})X_{\left(s\right)}\right\}}^{-1}{X}_{\left(s\right)}^{\top }(I-Q_{\left(s\right)})Z_G,$

and

(8)${\widehat{\alpha }}_{G,\left(s\right)}={\left(B_{\left(s\right)}^{\top }{B}_{\left(s\right)}\right)}^{-1}{B}_{\left(s\right)}^{\top }(Z_G-X_{\left(s\right)}{\widehat{\beta }}_{\left(s\right)}),$

where $Q_{\left(s\right)}=B_{\left(s\right)}{\left(B_{\left(s\right)}^{\top }{B}_{\left(s\right)}\right)}^{-1}{B}_{\left(s\right)}^{\top }$. Therefore, the estimator of $\mu$ under the sth candidate model is given by

(9)${\widehat{\mu }}_{G,\left(s\right)}=X_{\left(s\right)}{\widehat{\beta }}_{G,\left(s\right)}+B_{\left(s\right)}{\widehat{\alpha }}_{G,\left(s\right)}=\{Q_{\left(s\right)}+{\tilde{X}}_{\left(s\right)}{\left({\tilde{X}}_{\left(s\right)}^{\top }{\tilde{X}}_{\left(s\right)}\right)}^{-1}{\tilde{X}}_{\left(s\right)}^{\top }\}{Z}_G=P_{\left(s\right)}{Z}_G,$

where ${\tilde{X}}_{\left(s\right)}=(I-Q_{\left(s\right)})X_{\left(s\right)}$, and $P_{\left(s\right)}=Q_{\left(s\right)}+{\tilde{X}}_{\left(s\right)}{\left({\tilde{X}}_{\left(s\right)}^{\top }{\tilde{X}}_{\left(s\right)}\right)}^{-1}{\tilde{X}}_{\left(s\right)}^{\top }$. From Equation (9), we find that ${\widehat{\mu }}_{G,\left(s\right)}$ is linearly dependent on $Z_G$.

3. Model-Averaging Estimator and Weight Choice Criterion

Let $\omega ={({\omega }_1,\dots ,{\omega }_S)}^{\top }$ be a weight vector belonging to the set $W=\{\omega \in {[0,1]}^S:{\sum }_{s=1}^S{\omega }_s=1\}$, then the model-averaging estimator of $\mu$ can be formulated as

(10)${\widehat{\mu }}_G\left(\omega \right)=\sum _{s=1}^S{\omega }_s{\widehat{\mu }}_{G,\left(s\right)}=P\left(\omega \right)Z_G,$

where $P\left(\omega \right)={\sum }_{s=1}^S{\omega }_s{P}_{\left(s\right)}$.

Motivated by [7], we propose a leave-one-out cross-validation criterion to select the weight vector $\omega$ for Equation (10) in the PLM framework. Let $X_{\left(s\right)}^{[-i]}$, $B_{\left(s\right)}^{[-i]}$ and $Z_G^{[-i]}$ be the matrices/vectors $X_{\left(s\right)}$, $B_{\left(s\right)}$ and $Z_G$ with the ith row deleted. The leave-one-out estimator of ${\mu }_i$ in the sth candidate model is given by

(11)${\tilde{\mu }}_{G,\left(s\right),i}=X_{\left(s\right),i}^{\top }{\widehat{\beta }}_{G,\left(s\right)}^{[-i]}+B_{\left(s\right)}^{\top }\left(U_i\right){\widehat{\alpha }}_{G,\left(s\right)}^{[-i]},$

where

${\widehat{\beta }}_{G,\left(s\right)}^{[-i]}={\left\{X_{\left(s\right)}^{[-i]\top }(I-{\tilde{Q}}_{\left(s\right)})X_{\left(s\right)}^{[-i]}\right\}}^{-1}{X}_{\left(s\right)}^{[-i]\top }(I-{\tilde{Q}}_{\left(s\right)})Z_G^{[-i]},$

${\widehat{\alpha }}_{G,\left(s\right)}^{[-i]}={\left(B_{\left(s\right)}^{[-i]\top }{B}_{\left(s\right)}^{[-i]}\right)}^{-1}{B}_{\left(s\right)}^{[-i]\top }\left(Z_G^{[-i]}-X_{\left(s\right)}^{[-i]}{\widehat{\beta }}_{G,\left(s\right)}^{[-i]}\right),$

and ${\tilde{Q}}_{\left(s\right)}=B_{\left(s\right)}^{[-i]}{\left(B_{\left(s\right)}^{[-i]\top }{B}_{\left(s\right)}^{[-i]}\right)}^{-1}{B}_{\left(s\right)}^{[-i]\top }$. Denote the sth jackknife estimator and the corresponding jackknife version of the averaging estimator as ${\tilde{\mu }}_{G,\left(s\right)}={({\tilde{\mu }}_{G,\left(s\right),1},\dots ,{\tilde{\mu }}_{G,\left(s\right),n})}^{\top }$ and ${\tilde{\mu }}_G\left(\omega \right)={\sum }_{s=1}^S{\omega }_s{\tilde{\mu }}_{G,\left(s\right)}$. The leave-one-out cross-validation weight choice criterion is

(12)$C{V}_G\left(\omega \right)={\parallel Z_G-{\tilde{\mu }}_G\left(\omega \right)\parallel }^2,$

then minimizing $C{V}_G\left(\omega \right)$ over the space W yields the optimal weight vector. However, in practice, such a minimization process is computationally infeasible because the cumulative distribution function $G(·)$ in Equation (12) is unknown and needs to be estimated. Similar to [41], we can estimate $G(·)$ by the commonly used Kaplan–Meier estimator

(13)${\widehat{G}}_n\left(z\right)=1-\prod _{i=1}^n{\left[{\displaystyle \frac{n-i}{n-i+1}}\right]}^{I[Z_{\left(i\right)}\le z,{\delta }_{\left(i\right)}=0]},$

where $Z_{\left(1\right)}\le Z_{\left(2\right)}\le \dots \le Z_{\left(n\right)}$ denote the order statistics of $Z_1,Z_2,\dots ,Z_n$, and ${\delta }_{\left(i\right)}$ is the indicator corresponding to $Z_{\left(i\right)}$. In what follows, a letter subscripted by ${\widehat{G}}_n$ denotes that it is obtained by replacing G in its corresponding estimator with ${\widehat{G}}_n$. For instance, $Z_{{\widehat{G}}_n}$ is obtained by replacing G with its estimator ${\widehat{G}}_n$ in $Z_G$. Then a feasible counterpart of $C{V}_G\left(\omega \right)$ is given by

(14)$C{V}_{{\widehat{G}}_n}\left(\omega \right)={\parallel Z_{{\widehat{G}}_n}-{\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)\parallel }^2,$

where ${\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)={\sum }_{s=1}^S{\omega }_s{\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$, and ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}={({\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right),1},\dots ,{\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right),n})}^{\top }$. Minimizing $C{V}_{{\widehat{G}}_n}\left(\omega \right)$ with respect to $\omega$ over the set W leads to the jackknife choice of weight vector

(15)$\widehat{\omega }=arg\underset{\omega \in W}{min}C{V}_{{\widehat{G}}_n}\left(\omega \right).$

Plugging ${\widehat{G}}_n$ and $\widehat{\omega }$ into Equation (10) yields the model-averaging estimator of $\mu$, written as ${\widehat{\mu }}_{{\widehat{G}}_n}\left(\widehat{\omega }\right)$, which is named the censored partially linear model averaging (CPLMA) estimator hereafter.

However, minimizing the weight choice criterion (14) is not easy because the computation of ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$ requires n separate regressions, which is especially cumbersome when the number of candidate models and the sample sizes are large. Motivated by the computationally efficient cross-validation criterion introduced by [7] for linear regression model, we express ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$ in a simple form which yields an enormous reduction in calculation time. Let ${\varphi }_{ii}^{\left(s\right)}$ be the ith diagonal entry of $P_{\left(s\right)}$. From [20,42], ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$ can be conveniently written as

(16)${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}=\{P_{\left(s\right)}-D_{\left(s\right)}{A}_{\left(s\right)}\}{Z}_{{\widehat{G}}_n}={\tilde{P}}_{\left(s\right)}{Z}_{{\widehat{G}}_n},$

where $D_{\left(s\right)}=\mathrm{diag}(D_{11}^{\left(s\right)},\dots ,D_{nn}^{\left(s\right)})$, $D_{ii}^{\left(s\right)}={\varphi }_{ii}^{\left(s\right)}+{\left({\varphi }_{ii}^{\left(s\right)}\right)}^2+{\left({\varphi }_{ii}^{\left(s\right)}\right)}^3+\dots ={\varphi }_{ii}^{\left(s\right)}/(1-{\varphi }_{ii}^{\left(s\right)})$, and $A_{\left(s\right)}=I-P_{\left(s\right)}$. The shortcut formula of ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$ given in (16) indicates that all elements in ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$ can be simultaneously calculated based on all observations, which is much more convenient and time-saving than the standard method based on Equation (11). Let ${{\Lambda }}_{\left(s\right)}=D_{\left(s\right)}{A}_{\left(s\right)}$, ${\Lambda }\left(\omega \right)={\sum }_{s=1}^S{\omega }_s{{\Lambda }}_{\left(s\right)}$, and $A\left(\omega \right)={\sum }_{s=1}^S{\omega }_s{A}_{\left(s\right)}=I-P\left(\omega \right)$. The corresponding computational shortcut formula for the feasible jackknife criterion (14) then follows as

(17)$\begin{array}{cc}\hfill C{V}_{{\widehat{G}}_n}\left(\omega \right)& =\parallel Z_{{\widehat{G}}_n}-\{P\left(\omega \right)-{\Lambda }\left(\omega \right)\}{Z}_{{\widehat{G}}_n}{\parallel }^2\hfill \\ \hfill & =\parallel \{A\left(\omega \right)+{\Lambda }\left(\omega \right)\}{Z}_{{\widehat{G}}_n}{\parallel }^2\hfill \\ \hfill & =Z_{{\widehat{G}}_n}^{\top }{\{A\left(\omega \right)+{\Lambda }\left(\omega \right)\}}^{\top }\{A\left(\omega \right)+{\Lambda }\left(\omega \right)\}{Z}_{{\widehat{G}}_n}\hfill \\ \hfill & ={\omega }^{\top }{H}_{{\widehat{G}}_n}^{\top }{H}_{{\widehat{G}}_n}\omega ,\hfill \end{array}$

where $H_{{\widehat{G}}_n}=\{(A_{\left(1\right)}+{{\Lambda }}_{\left(1\right)})Z_{{\widehat{G}}_n},\dots ,(A_{\left(S\right)}+{{\Lambda }}_{\left(S\right)})Z_{{\widehat{G}}_n}\}$ is an $n\times S$ matrix. From Equation (17), we observe that the minimization of $C{V}_{{\widehat{G}}_n}\left(\omega \right)$ is a standard quadratic programming problem, which can be performed by various existing software packages, for example, the quadprog package in R [43].

4. Asymptotic Optimality

In this section, we demonstrate that the resulting weight vector, which is obtained by minimizing the weight choice criterion $C{V}_{{\widehat{G}}_n}\left(\omega \right)$, is asymptotically optimal under some mild conditions.

Define the squared loss as $L_G\left(\omega \right)={\parallel {\widehat{\mu }}_G\left(\omega \right)-\mu \parallel }^2$, and the corresponding risk function as $R_G\left(\omega \right)=E\left(L_G\left(\omega \right)\right|X,U)={\parallel A\left(\omega \right)\mu \parallel }^2+\mathrm{tr}\left\{P\left(\omega \right){\Omega }_{G}P\left(\omega \right)\right\}$. Let $\overline{p}={max}_{1\le s\le S}{p}_s$, $\overline{q}={max}_{1\le s\le S}\mathrm{rank}\left(B_{\left(s\right)}\right)$, ${\xi }_G={inf}_{\omega \in W}{R}_G\left(\omega \right)$, and ${\omega }_s^0$ be a $S\times 1$ vector with the sth entry taking on 1 and the others taking on 0. To prove the asymptotic optimality of the model-averaging estimator ${\widehat{\mu }}_{{\widehat{G}}_n}\left(\widehat{\omega }\right)$, we list the following regularity conditions, where all limiting processes correspond to $n\to \infty$.

• (Condition (C.1)) ${\tau }_F<{\tau }_G$, where ${\tau }_L=inf\{t:L\left(t\right)=1\}$ for any distribution function L.

• (Condition (C.2)) $\overline{\lambda }\left({\Omega }_G\right)\le C_G$, where $\overline{\lambda }(·)$ denotes the maximum singular value of a matrix, and $C_G$ is a constant.

• (Condition (C.3)) ${\xi }_G^{-2}{\sum }_{s=1}^S{R}_G\left({\omega }_s^0\right)\to 0$, a.s.

• (Condition (C.4)) $\overline{p}{\xi }_G^{-1}\to 0$ and $\overline{q}{\xi }_G^{-1}\to 0$, a.s.

• (Condition (C.5)) ${\parallel \mu \parallel }^2=O\left(n\right)$, a.s.

• (Condition (C.6)) ${\varphi }_{ii}^{\left(s\right)}\le C_s{n}^{-1}\mathrm{tr}\left(P_{\left(s\right)}\right)$, a.s., where $C_s$ is a constant.

• (Condition (C.7)) The function g belongs to a class of functions $\mathcal{A}$, whose rth derivative $g_{★}^{\left(r\right)}$ exsits and is Lipschitz of order ${\alpha }_0$. That is,

  $\mathcal{A}=\{g_{★}(·):|g_{★}^{\left(r\right)}\left(t_1\right)-g_{★}^{\left(r\right)}\left(t_0\right)|\le M^{★}|t_1-t_0{|}^{{\alpha }_0}for{t}_1,t_0\in \mathcal{U}\},$

  for some positive constant $M^{★}$, where $\mathcal{U}$ is the support of U, r is a nonnegative integer and ${\alpha }_0\in (0,1]$ such that $r+{\alpha }_0>0.5$.

Condition (C.1), which is the same as condition (C5) in [35], is widely used to ensure the uniform convergence of the Kaplan–Meier estimator ${\widehat{G}}_n(·)$ in studies of the censored data. Condition (C.2) imposes a mild restriction on the maximum singular value of the covariance matrix ${\Omega }_G$, which is also used by [44]. Condition (C.3), which is from condition (21) of [45], is less restrictive than the commonly used condition $S{\xi }_G^{-2N}{\sum }_{s=1}^S{\left\{R_G\left({\omega }_s^0\right)\right\}}^N\to 0,a.s.$, for some constant $N\ge 1$ in model-averaging references. Condition (C.4) places constraints on the growth rates of $\overline{p}$ and $\overline{q}$, which is similar to condition (22) in [45]. Condition (C.5) discusses the sum of ${\mu }_i^2$, and is frequently used in the model-averaging literature, such as [23,24], and so on. Condition (C.6) is a common assumption utilized to guarantee the asymptotic optimality of cross-validation; see [7,24], for instance. Conditions (C.3)–(C.6) require almost sure convergence, which ensure that the result (18) holds whether the covariates X and U are random or not. Specifically, when X and U are nonstochastic, we only need to assume convergence in probability in conditions (C.3)–(C.6); see [46]. Otherwise, we should impose almost sure convergence to guarantee that the proof method, which is used in the case of nonstochastic, is still effective. Condition (C.7) is required for the B-spline approximation in PLM; see [39,47].

Theorem 1 indicates that the CPLMA estimator proposed in this paper is asymptotically optimal in the sense that its squared error loss is asymptotically equivalent to that of the infeasible best possible model-averaging estimator in PLM framework. The proof of Theorem 1 is shown in the Appendix A.

Theorem 1.

Under Conditions (C.1)–(C.7), we have

(18) ${\displaystyle \frac{L_{{\widehat{G}}_n}\left(\widehat{\omega }\right)}{{inf}_{\omega \in W}{L}_{{\widehat{G}}_n}\left(\omega \right)}}\to 1$

in probability as $n\to \infty$.

5. A Simulation Study

In this section, a simulation experiment is conducted to investigate the finite sample performance of the CPLMA estimator, which arises from the proposed leave-one-out cross-validation weight choice approach, in PLM with censored responses. We compare it with several popular information-criterion-based model-selection methods as well as other model-averaging procedures.

5.1. The Design of Simulation

The data-generating process in this part is similar to the infinite-order regression model proposed by [17], except that responses are subject to censoring and a nonparametric function is included in addition to the linear part. Specifically, the data are generated by the following regression model:

(19)$Y_i={\mu }_i+{ϵ}_i=\sum _{j=1}^{200}{x}_{ij}{\beta }_j+g\left(U_i\right)+{ϵ}_i,$

where $X_i={(x_{i1},\dots ,x_{i200})}^{\top }$, the covariates in linear component, follows a multivariate normal distribution with mean 0 and covariance ${0.5}^{\left|j_1-j_2\right|}$ between $x_{i{j}_1}$ and $x_{i{j}_2}$. The coefficients of the linear part are set as ${\beta }_j=1/j^{\alpha }$, and the parameter $\alpha$ varies between 2 and $0.5$. The larger $\alpha$ implies that the coefficient decays more quickly as j increases. The nonparametric function is $g\left(U_i\right)=sin\left(2\pi U_i^2\right)$, where $U_i$ is generated from the uniform distribution on $[0,1]$. The model error ${ϵ}_i$ follows a normal distribution $N(0,{\eta }^2(x_{i2}^2+0.01))$. We choose the value of $\eta$ so that $R^2=var\left({\mu }_i\right)/var\left(Y_i\right)$ varies from $0.1$ to $0.9$, where $var(·)$ denotes the sample variance. In addition, the censoring variable $C_i$ is generated from a uniform distribution on interval $[a_1,a_2]$, where different values of $a_1$ and $a_2$ are selected to yield a censoring rate (CR) of either $20\%$ or $40\%$. In order to evaluate the performance of the methods as comprehensively as possible, we consider two designs and set the sample size $n=50$, 75, 100, 200, 300, and 400.

Design 1

(non-nested setting).The linear parts of all candidate models are a subset of$\{x_{i1},\dots ,x_{i5}\}$(with the rest of$X_i$being ignored), so the number of candidate models is$2^5=32$.

Design 2

(nested setting).The sth candidate model includes the first s linear variables. The number of the candidate models is determined by$S=\mathrm{INT}\left(3n^{1/3}\right)$, where$\mathrm{INT}\left(b\right)$denotes the nearest integer of b. Therefore,$S=11$, 13, 14, 18, 20, and 22, for $n=50$, 75, 100, 200, 300, and 400, respectively.

5.2. Estimation and Comparison

Suggested by [35], the cubic B-spline is used to approximate the nonparametric function, and the spline basis matrix is generated by $bs(·)$ in splines package with R project [48]. To select the number of knots, we set $\eta =1$, $\alpha =2$, and investigate the impact of the number of knots on the risk of CPLMA estimator under different scenarios. Figure 1 shows the variation for the mean of risks when the number of knots varies, over 500 replications for four combinations of designs and CRs considered. From Figure 1, we see that in almost all cases, 1 knot yields the smallest mean of risk, except that in the left lower panel, 1 knot leads to a mean of risk that is second to 2 knots but best among the remaining number of knots. In addition, in all cases, the mean of risk increases with the number of knots when the number of knots exceeds 2. This observation coincides with the findings in [39] that the larger number of knots results in a more serious overfitting effect. Therefore, the number of knots is set to be 1 in the simulation studies.

We compare the performance of the CPLMA method with two traditional model-selection methods (AIC and BIC) and two model-averaging methods based on scores of information criteria (SAIC and SBIC). For the sth model, we calculate AIC and BIC scores by

$AI{C}_s=log\left({\widehat{\sigma }}_{{\widehat{G}}_n,\left(s\right)}^2\right)+2n^{-1}tr\left(P_{\left(s\right)}\right),$

and

$BI{C}_s=log\left({\widehat{\sigma }}_{{\widehat{G}}_n,\left(s\right)}^2\right)+n^{-1}tr\left(P_{\left(s\right)}\right)log\left(n\right),$

respectively, where ${\widehat{\sigma }}_{{\widehat{G}}_n,\left(s\right)}^2=n^{-1}{\parallel Z_{{\widehat{G}}_n}-{\widehat{\mu }}_{{\widehat{G}}_n,\left(s\right)}\parallel }^2$, and ${\widehat{\mu }}_{{\widehat{G}}_n,\left(s\right)}$ is obtained by replacing G in Equation (9) with ${\widehat{G}}_n$. Both methods pick the model with the smallest information criterion score.

For SAIC and SBIC, the weights of sth model are defined as

${\omega }_{AI{C}_s}=exp(-AI{C}_s/2)/\sum _{s^{\prime }=1}^{S}exp(-AI{C}_{s^{\prime }}/2),$

and

${\omega }_{BI{C}_s}=exp(-BI{C}_s/2)/\sum _{s^{\prime }=1}^{S}exp(-BI{C}_{s^{\prime }}/2),$

respectively. To evaluate these five methods, we draw 500 independent samples of size n, and compute the risk of the estimator of $\mu$. For comparison convenience, the risks of all estimators are normalized by the risk produced by the AIC method.

5.3. Results

The simulation results for Design 1 are presented in Figure 2 and Figure 3 for the censoring rate of $20\%$ and in Figure 4 and Figure 5 for the censoring rate of $40\%$. These four figures show that our method CPLMA leads to the smallest risk in most cases, except that both SAIC and SBIC sometimes have marginal advantages over ours when $R^2$ is large, and the advantage of SBIC is more obvious when n is small. In particular, comparing the simulation results between $\alpha =0.5$ and $\alpha =2$, we find that the performance of CPLMA is better when $\alpha$ is small. As expected, SAIC and SBIC invariably produce vastly more accurate outcomes than their respective model-selection counterparts.

The simulation results for Design 2 are depicted in Figure 6, Figure 7, Figure 8 and Figure 9 with $CR\approx 20\%$ and $40\%$, from which we see that in most cases our proposed CPLMA method still outperforms its rivals in terms of risk. The superiority of the CPLMA method over the other methods is more apparent than in Design 1. Additionally, we find that BIC-based model selection and averaging estimators have much worse risk performance than the other three estimators when $R^2$ is small, which is different from the simulation results in Design 1. We also note that SAIC and our method CPLMA almost perform equally well when $R^2$ is very large.

In a word, no matter whether the candidate models are nested or not, our proposal, CPLMA, is superior to the traditional model-selection and model-averaging methods for all the combinations of censoring rates and sample sizes considered.

6. Real Data Analysis

In this section, we apply the proposed CPLMA method to analyze two real datasets by R software. The first real dataset can be found in the R package “survival” [49], and the other is available at http://llmpp.nih.gov/MCL (accessed on 18 January 2023).

6.1. Primary Biliary Cirrhosis Dataset Study

The primary biliary cirrhosis (PBC) dataset includes information of 424 patients, which were collected at Mayo Clinic from January 1974 to May 1984, and has been extensively explored by [34,35,37,50,51]. Following the related literature, we restrict our attention to the $n=276$ patients without missing observations, each of whom contains the data of 17 covariates. There are 111 deaths among 276 patients, which leads to $60\%$ of censoring.

In this dataset, the dependent variable is the log number of days between registration and the earlier of death or study analysis time in 1986. The 17 covariates include age (in years), albumin (serum albumin in g/dL), alk.phos (alkaline phosphotase in U/L), bili (serum bilirunbin in mg/dL), chol (serum cholesterol in mg/dL), copper (urine copper in ug/day), platelet (platelet count), protime (standardized blood clotting time in seconds), ast (aspartate aminotransferase, once called SGOT in U/mL), trig (triglycerides in mg/dL), ascites (presence of ascites, $0=$ no, $1=$ yes), edema ($0=$ no; $0.5=$ yes, but responded to diuretic treatment; $1=$ yes, did not respond to treatment), hepato (presence of hepatomegaly, $0=$ no, $1=$ yes), sex ($0=$ male, $1=$ female), spiders (presence of spiders, $0=$ no, $1=$ yes), stage (histologic state of disease, graded 1, 2, 3, or 4), and trt (treatment code, $1=$ D-pencillamine, $2=$ placebo). The first 10 variables are continuous, which are standardized to have mean of 0 and variance of 1 in the analysis.

A total of 17 covariates lead to a huge number of candidate models, which brings heavy calculation burden. Ref. [50] pointed out that only eight covariates, that is, ge, edema, bili, albumin, copper, ast, protime, and stage, have significant impact on the response variable, and ref. [35] found that albumin has a functional impact on the response variable. Thus, we only consider eight significant covariates. Specifically, we assign albumin in the nonparametric part while including others in the linear part of PLM, and we run model selection and averaging on the covariates in the linear part. Accordingly, there are $2^7=128$ prepared models as candidates. Similar to [35], we also use the cubic B-spline with two knots to approximate the nonparametric component.

To evaluate the prediction effect of two model-selection methods (AIC and BIC) and three model-averaging methods (SAIC, SBIC, and CPLMA), we randomly separate the data into a training sample and a test sample. Let $n_0$ be the size of the training sample, and $n_1=n-n_0$ be the size of the test sample. We set $n_0$ to $140,160,180,200,220$, and 240. The mean-squared prediction error (MSPE) is used to describe the out-of-sample prediction performance of the proposed CPLMA and its competitors. We further calculate the mean and the median of the MSPE for each method based on 1000 replications. Specifically,

(20)${\mathrm{MSPE}}_{\mathrm{mean}}=\frac{1}{1000}\sum _{d=1}^{1000}{\mathrm{MSPE}}^{\left(\mathrm{d}\right)},$

and

(21)${\mathrm{MSPE}}_{\mathrm{median}}={\mathrm{median}}_{\mathrm{d}=1,2,\dots ,1000}{\mathrm{MSPE}}^{\left(\mathrm{d}\right)},$

where

(22)${\mathrm{MSPE}}^{\left(\mathrm{d}\right)}={\displaystyle \frac{1}{n_1}}\sum _{i=n_0+1}^n{(Z_{{\widehat{G}}_n,i}^{\left(d\right)}-{\widehat{\mu }}_i^{\left(d\right)})}^2,$

and ${\widehat{\mu }}_i^{\left(d\right)}$ is the predicted value of ${\mu }_i$ in the dth repetition.

To facilitate comparison, we calculate the ratio of MSPE for a given method to the MSPE produced by AIC, which is referred to relative MSPE (RMSPE). Table 2 reports the mean and median of RMSPE across 1000 repetitions. We see that our proposed CPLMA always yields the lowest mean and median of RMSPE for all considered training sample sizes. In all cases, the values of RMSPE for BIC are bigger than 1 and the values of RMSPE for three model-averaging methods are smaller than 1, which indicates that in terms of prediction performance, AIC is the overwhelming favorite of two model-selection methods, and model-averaging methods outperform model-selection methods.

Table 3 presents the Diebold and Mariano (DM) [52] test results for the differences in MSPE, where a positive DM statistic implies that the method in the numerator yields a larger MSPE than the method in the denominator. The results in columns 6, 9, 11, and 12 indicate that the differences between CPLMA and its competitors are statistically significant, and our method always produces smaller MSPE than other four methods, which again demonstrates the superiority of our proposal. The results in column 3 show that AIC is significantly better than BIC, which coincides with the finding in Table 2. Columns 4 and 8 indicate that SAIC and SBIC are significantly different from their respective model-selection counterparts.

6.2. Mantle Cell Lymphoma Data Analysis

The mantle cell lymphoma (MCL) dataset contains 92 patients who were classified as having MCL based on established morphologic and immunophenotype criteria. Since 2003, this dataset has been widely studied by [53,54]. The response variable of interest is time (time of follow-up in year). The variable status denotes patient status at follow-up ($1=$ death, $0=$ censored). The six covariates include indicator of INK/ARF deletion ($1=$ yes, $0=$ no), indicator of ATM deletion ($1=$ yes, $0=$ no), indicator of P-53 deletion ($1=$ yes, $0=$ no), cyclinD-1 taqman result, BMI expression, and proliferation signature averages. After removing seven records with missing covariates, we focus on the information of the remaining 85 patients, and the censoring rate is $29.41\%$.

Ref. [54] found that BMI expression has a functional impact on the response variable; therefore, we establish a full PLM model with BMI expression as nonparametric variable and other covariates as linear variables. Then, we can obtain $2^5=32$ candidate models to conduct model selection and model averaging. Let the size of training sample be $n_0=55$ or 65, and the mean and median of the RMSPE across 1000 repetitions are shown in Table 4. It can be seen from Table 4 that, both in terms of mean and median, our method CPLMA is significantly superior to other competitive methods. Figure 10 shows that the variation of MSPE for CPLMA is minor relative to that of other methods, regardless of if $n_0$ is 55 or 65.

7. Conclusions

In the context of the semiparametric partially linear models with censored responses, we develop a jackknife model-averaging method which selects the weights by minimizing a leave-one-out cross-validation criterion, in which the B-splines are used to approximate the nonparametric function and the least-squares estimation is applied to estimate the unknown parameters in each candidate model. The resulting model-averaging estimator, CPLMA estimator, is shown to be asymptotically optimal. A simulation study and two real data examples indicate that our method posseses some advantages over other model-selection and model-averaging methods.

Based on the results in this paper, we can further explore the optimal model averaging for the semiparametric partially linear quantile regression models with censored data. In addition, it is worthwhile to apply other optimal model-averaging methods, such as the model-averaging method based on Kullback–Leibler distance, to generalized partially linear models with censored responses.

Author Contributions

Conceptualization, W.C.; methodology, G.H., W.C. and J.Z.; software, G.H. and J.Z.; supervision, W.C. and J.Z.; writing—original draft, G.H.; writing—review and editing, J.Z. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The PBC dataset is available in the R package “survival”, and the MCL dataset is available at http://llmpp.nih.gov/MCL(accessed on 18 January 2023).

Acknowledgments

The authors would like to thank the reviewers and editors for their careful reading and constructive comments.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. The curves of the mean of risk with the number of knots over 500 replications.

Figure 2. Risk comparisons for Design 1 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 3. Risk comparisons for Design 1 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 4. Risk comparisons for Design 1 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 5. Risk comparisons for Design 1 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 6. Risk comparisons for Design 2 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 7. Risk comparisons for Design 2 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 8. Risk comparisons for Design 2 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 9. Risk comparisons for Design 2 when [Forumla omitted. See PDF.] and the censoring rate is about [Forumla omitted. See PDF.].

Figure 10. Boxplots for MSPEs of five methods for the MCL data.

Table 1

The basic notations used in this paper.

Notations	Descriptions
$T_i$	The survival time of the ith subject
$Y_i$	The response variable, a transformation of $T_i$
$X_i$	The covariate vector of the ith subject
${\delta }_i$	The censoring indicator of the ith subject
$C_i$	The last follow up time of the ith subject
$Z_i$	The observed time, equal to $min(Y_i,C_i)$
$G(·)$	The cumulative distribution function of $C_i$
$Z_G$	$n\times 1$ synthetic response vector
$\mu$	$n\times 1$ conditional mean vector of the response
$B_{\left(s\right)}$	$n\times k_n$ B-spline basis matrix for the sth model
${\beta }_{\left(s\right)}$	$p_s\times 1$ linear regression coefficient vector for the sth model
${\alpha }_{\left(s\right)}$	$k_n\times 1$ spline coefficient vector for the sth model
${\widehat{G}}_n(·)$	The Kaplan–Meier estimator of $G(·)$
${\widehat{\beta }}_{G,\left(s\right)}$	The estimator of ${\beta }_{\left(s\right)}$ with $G(·)$
${\widehat{\alpha }}_{G,\left(s\right)}$	The estimator of ${\alpha }_{\left(s\right)}$ with $G(·)$
${\widehat{\mu }}_{G,\left(s\right)}$	The estimator of $\mu$ for the sth model with $G(·)$
${\widehat{\mu }}_G\left(\omega \right)$	The model-averaging estimator of $\mu$ with $G(·)$
${\tilde{\mu }}_{G,\left(s\right)}$ or ${\tilde{\mu }}_{{\widehat{G}}_n,\left(s\right)}$	the sth jackknife estimator of $\mu$ with $G(·)$ or ${\widehat{G}}_n(·)$
${\tilde{\mu }}_G\left(\omega \right)$ or ${\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)$	the jackknife model-averaging estimator of $\mu$ with $G(·)$ or ${\widehat{G}}_n(·)$

Table 2

The mean and median of RMSPE across 1000 repetitions.

${\mathit{n}}_0$	Method	BIC	SAIC	SBIC	CPLMA
140	mean	1.005	0.987	0.983	0.979
	median	1.014	0.992	0.995	0.987
160	mean	1.005	0.987	0.983	0.982
	median	1.006	0.986	0.988	0.984
180	mean	1.011	0.991	0.991	0.986
	median	1.011	0.990	0.990	0.981
200	mean	1.014	0.993	0.995	0.989
	median	1.012	0.984	0.990	0.976
220	mean	1.012	0.995	0.997	0.994
	median	1.020	0.994	1.003	0.993
240	mean	1.008	0.995	0.998	0.993
	median	1.017	0.996	0.999	0.988

Table 3

Diebold–Mariano test results for the differences in MSPE.

${\mathit{n}}_0$	Method	${\displaystyle \frac{\mathbf{AIC}}{\mathbf{BIC}}}$	${\displaystyle \frac{\mathbf{AIC}}{\mathbf{SAIC}}}$	${\displaystyle \frac{\mathbf{AIC}}{\mathbf{SBIC}}}$	${\displaystyle \frac{\mathbf{AIC}}{\mathbf{CPLMA}}}$	${\displaystyle \frac{\mathbf{BIC}}{\mathbf{SAIC}}}$	${\displaystyle \frac{\mathbf{BIC}}{\mathbf{SBIC}}}$	${\displaystyle \frac{\mathbf{BIC}}{\mathbf{CPLMA}}}$	${\displaystyle \frac{\mathbf{SAIC}}{\mathbf{SBIC}}}$	${\displaystyle \frac{\mathbf{SAIC}}{\mathbf{CPLMA}}}$	${\displaystyle \frac{\mathbf{SBIC}}{\mathbf{CPLMA}}}$
140	DM	−3.013	15.130	12.335	14.858	16.743	25.341	16.633	4.973	6.361	2.834
	p-value	0.003	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.005
160	DM	−3.014	16.607	12.490	14.862	17.196	30.331	18.474	4.995	5.538	0.942
	p-value	0.003	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.347
180	DM	−8.082	12.355	7.874	11.238	22.554	34.561	21.914	−0.348	4.439	5.679
	p-value	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.728	0.000	0.000
200	DM	−11.473	12.320	4.744	11.393	23.962	32.286	22.550	−3.721	5.288	8.690
	p-value	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
220	DM	−9.509	8.500	2.308	5.587	19.004	22.316	17.152	−4.011	1.085	5.059
	p-value	0.000	0.000	0.021	0.000	0.000	0.000	0.000	0.000	0.278	0.000
240	DM	−5.484	7.332	1.441	5.848	12.901	15.427	12.998	−4.175	2.110	6.561
	p-value	0.000	0.000	0.150	0.000	0.000	0.000	0.000	0.000	0.035	0.000

Table 4

The mean and median of RMSPE across 1000 repetitions.

${\mathit{n}}_0$	Method	BIC	SAIC	SBIC	CPLMA
55	mean	1.011	0.982	0.988	0.923
	median	0.951	0.965	0.937	0.918
65	mean	0.992	0.976	0.987	0.947
	median	0.982	0.970	0.973	0.939

Appendix A

To prove Theorem 1, we first give some notations which are used in the following context. Similar to $L_G\left(\omega \right)$ and $R_G\left(\omega \right)$, we define the loss function of ${\tilde{\mu }}_G\left(\omega \right)$ as ${\tilde{L}}_G\left(\omega \right)={\parallel {\tilde{\mu }}_G\left(\omega \right)-\mu \parallel }^2$, and the risk function as ${\tilde{R}}_G\left(\omega \right)=E\left({\tilde{L}}_G\left(\omega \right)\right|X,U)$. Straightforward calculation yields (A1)$${\tilde{R}}_G\left(\omega \right)=\parallel \tilde{A}\left(\omega \right)\mu {\parallel }^2+\mathrm{tr}\left[\tilde{P}\left(\omega \right){\Omega }_G{\left\{\tilde{P}\left(\omega \right)\right\}}^{\top }\right],$$ where $\tilde{P}\left(\omega \right)={\displaystyle \sum _{s=1}^S}{\omega }_s{\tilde{P}}_{\left(s\right)}$, and $\tilde{A}\left(\omega \right)=I-\tilde{P}\left(\omega \right)$.

Lemma A1.

Under Conditions (C.2)–(C.7), we have the following results: (A2) $$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{R_G\left(\omega \right)}{{\tilde{R}}_G\left(\omega \right)}}-1|=o_p\left(1\right),\end{array}$$ (A3) $$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{{\tilde{L}}_G\left(\omega \right)}{{\tilde{R}}_G\left(\omega \right)}}-1|=o_p\left(1\right),\end{array}$$ (A4) $$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{L_G\left(\omega \right)}{R_G\left(\omega \right)}}-1|=o_p\left(1\right).\end{array}$$

Proof. 

According to the proof of (A.45), (A.48) and (A.44) in [45], we know that Equations (A2)–(A4) are satisfied. □

Lemma A2.

If Conditions (C.4)–(C.7) hold, then (A5) $$\begin{array}{c}\hfill \underset{\omega \in W}{sup}\overline{\lambda }\left\{P\left(\omega \right)\right\}\le 2,\end{array}$$ (A6) $$\begin{array}{c}\hfill \underset{\omega \in W}{sup}\overline{\lambda }\left\{{\Lambda }\left(\omega \right)\right\}=o_p\left(1\right),\end{array}$$ and (A7) $$\begin{array}{c}\hfill \underset{\omega \in W}{sup}\overline{\lambda }\left\{\tilde{P}\left(\omega \right)\right\}=O_p\left(1\right).\end{array}$$

Proof. 

By the inequalities for maximum singular value, we obtain (A8)$$\begin{array}{cc}\hfill \underset{\omega \in W}{sup}\overline{\lambda }\left\{P\left(\omega \right)\right\}& =\underset{\omega \in W}{sup}\overline{\lambda }\left\{\sum _{s=1}^S{\omega }_s{P}_{\left(s\right)}\right\}\le \underset{1\le s\le S}{\max }\overline{\lambda }\left\{P_{\left(s\right)}\right\}\hfill \\ \hfill & \le \underset{1\le s\le S}{\max }\left[\overline{\lambda }\left\{Q_{\left(s\right)}\right\}+\overline{\lambda }\left\{{\tilde{X}}_{\left(s\right)}{\left({\tilde{X}}_{\left(s\right)}^{\top }{\tilde{X}}_{\left(s\right)}\right)}^{-1}{\tilde{X}}_{\left(s\right)}^{\top }\right\}\right]\le 2,\hfill \end{array}$$ where the last inequality holds because the matrices $Q_{\left(s\right)}$ and ${\tilde{X}}_{\left(s\right)}{\left({\tilde{X}}_{\left(s\right)}^{\top }{\tilde{X}}_{\left(s\right)}\right)}^{-1}{\tilde{X}}_{\left(s\right)}^{\top }$ are all symmetric and idempotent.

By the definition of ${{\Lambda }}_{\left(s\right)}$, we have (A9)$$\begin{array}{cc}\hfill \underset{\omega \in W}{sup}\overline{\lambda }\left\{{\Lambda }\left(\omega \right)\right\}& =\underset{\omega \in W}{sup}\overline{\lambda }\left\{\sum _{s=1}^S{\omega }_s{{\Lambda }}_{\left(s\right)}\right\}\le \underset{1\le s\le S}{\max }\overline{\lambda }\left({{\Lambda }}_{\left(s\right)}\right)\hfill \\ \hfill & \le \underset{1\le s\le S}{\max }\overline{\lambda }\left(D_{\left(s\right)}{A}_{\left(s\right)}\right)\le \underset{1\le s\le S}{\max }\overline{\lambda }\left(D_{\left(s\right)}\right)\overline{\lambda }\left(I-P_{\left(s\right)}\right)\hfill \\ \hfill & \le \underset{1\le s\le S}{\max }\overline{\lambda }\left(D_{\left(s\right)}\right)\left\{1+\overline{\lambda }\left(P_{\left(s\right)}\right)\right\}\le \underset{1\le s\le S}{\max }3\overline{\lambda }\left(D_{\left(s\right)}\right)\hfill \\ \hfill & \le \underset{1\le s\le S}{\max }\underset{1\le i\le n}{\max }{\displaystyle \frac{3{\varphi }_{ii}^{\left(s\right)}}{1-{\varphi }_{ii}^{\left(s\right)}}}=O_p\left({\displaystyle \frac{\overline{p}+\overline{q}}{n}}\right),\hfill \end{array}$$ where the last equality is obtained based on Condition (C.6). Then, Equation (A6) is derived by Condition (C.4).

From (A5) and (A6), it can be shown that (A10)$$\begin{array}{cc}\hfill \underset{\omega \in W}{sup}\overline{\lambda }\left\{\tilde{P}\left(\omega \right)\right\}& \le \underset{1\le s\le S}{\max }\overline{\lambda }\left({\tilde{P}}_{\left(s\right)}\right)=\underset{1\le s\le S}{\max }\overline{\lambda }\left(P_{\left(s\right)}-{{\Lambda }}_{\left(s\right)}\right)\hfill \\ \hfill & \le \underset{1\le s\le S}{\max }\overline{\lambda }\left(P_{\left(s\right)}\right)+\underset{1\le s\le S}{\max }\overline{\lambda }\left({{\Lambda }}_{\left(s\right)}\right)\hfill \\ \hfill & \le 2+O_p\left({\displaystyle \frac{\overline{p}+\overline{q}}{n}}\right)=O_p\left(1\right).\hfill \end{array}$$

The proof of Lemma A2 is completed. □

Lemma A3.

Assuming that Conditions (C.1), (C.2), and (C.5) are satisfied, we obtain (A11) $$\begin{array}{c}\hfill \parallel Z_G-Z_{{\widehat{G}}_n}{\parallel }^2=O_p\left(1\right).\end{array}$$

Proof. 

This result is from Lemma $6.2$ in [37] directly; we omit the proof procedure. □

Proof of Theorem 1. 

According to [44,45], Theorem 1 is valid if we can prove (A12)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{R_G\left(\omega \right)}{{\tilde{R}}_G\left(\omega \right)}}-1|=o_p\left(1\right),\end{array}$$ (A13)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{{\tilde{L}}_G\left(\omega \right)}{{\tilde{R}}_G\left(\omega \right)}}-1|=o_p\left(1\right),\end{array}$$ (A14)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{L_{{\widehat{G}}_n}\left(\omega \right)}{R_G\left(\omega \right)}}-1|=o_p\left(1\right),\end{array}$$ and (A15)$$\begin{array}{c}\hfill {\displaystyle \frac{{\tilde{L}}_G\left(\widehat{\omega }\right)}{{inf}_{\omega \in W}{\tilde{L}}_G\left(\omega \right)}}-1=o_p\left(1\right).\end{array}$$

By Lemma A1 and Conditions (C.2)–(C.6), Equations (A12) and (A13) are satisfied. Next, we present the proofs of (A14) and (A15), which completes the proof of Theorem 1.

For (A14), because (A16)$$\begin{array}{cc}\hfill \underset{\omega \in W}{sup}& |{\displaystyle \frac{L_{{\widehat{G}}_n}\left(\omega \right)}{R_G\left(\omega \right)}}-1|=\underset{\omega \in W}{sup}|{\displaystyle \frac{\parallel {\widehat{\mu }}_{{\widehat{G}}_n}{\left(\omega \right)-\mu \parallel }^2}{R_G\left(\omega \right)}}-1|\hfill \\ \hfill & =\underset{\omega \in W}{sup}|{\displaystyle \frac{\parallel \mu -{\widehat{\mu }}_G\left(\omega \right)+{\widehat{\mu }}_G\left(\omega \right)-{\widehat{\mu }}_{{\widehat{G}}_n}{\left(\omega \right)\parallel }^2}{R_G\left(\omega \right)}}-1|\hfill \\ \hfill & \le \underset{\omega \in W}{sup}|{\displaystyle \frac{L_G\left(\omega \right)}{R_G\left(\omega \right)}}-1|+2\underset{\omega \in W}{sup}\left|{\displaystyle \frac{{\left\{L_G\left(\omega \right)\right\}}^{1/2}\parallel {\widehat{\mu }}_G\left(\omega \right)-{\widehat{\mu }}_{{\widehat{G}}_n}\left(\omega \right)\parallel }{R_G\left(\omega \right)}}\right|\hfill \\ \hfill & +\underset{\omega \in W}{sup}\left|{\displaystyle \frac{\parallel {\widehat{\mu }}_G\left(\omega \right)-{\widehat{\mu }}_{{\widehat{G}}_n}{\left(\omega \right)\parallel }^2}{R_G\left(\omega \right)}}\right|,\hfill \end{array}$$ it is sufficient to verify that (A17)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}|{\displaystyle \frac{L_G\left(\omega \right)}{R_G\left(\omega \right)}}-1|=o_p\left(1\right),\end{array}$$ and (A18)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}\left|{\displaystyle \frac{\parallel {\widehat{\mu }}_G\left(\omega \right)-{\widehat{\mu }}_{{\widehat{G}}_n}{\left(\omega \right)\parallel }^2}{R_G\left(\omega \right)}}\right|=o_p\left(1\right).\end{array}$$

Equation (A17) can be directly obtained by Lemma A1. As for (A18), by Cauchy–Schwarz inequality, we have (A19)$$\begin{array}{cc}\hfill \underset{\omega \in W}{sup}& \left|{\displaystyle \frac{\parallel {\widehat{\mu }}_G\left(\omega \right)-{\widehat{\mu }}_{{\widehat{G}}_n}{\left(\omega \right)\parallel }^2}{R_G\left(\omega \right)}}\right|=\underset{\omega \in W}{sup}\left|{\displaystyle \frac{\parallel P\left(\omega \right)Z_G-P\left(\omega \right)Z_{{\widehat{G}}_n}{\parallel }^2}{R_G\left(\omega \right)}}\right|\hfill \\ \hfill & \le {\xi }_G^{-1}{\overline{\lambda }}^2\left\{P\left(\omega \right)\right\}\parallel Z_G-Z_{{\widehat{G}}_n}{\parallel }^2\le 4{\xi }_G^{-1}\parallel Z_G-Z_{{\widehat{G}}_n}{\parallel }^2=o_p\left(1\right),\hfill \end{array}$$ where the last equality follows Lemma A3 and ${\xi }_G\to \infty$, which is implied by Condition (C.3). Then Equation (A14) is obtained.

A simple calculation yields (A20)$$\begin{array}{cc}\hfill C{V}_{{\widehat{G}}_n}\left(\omega \right)& =\parallel Z_{{\widehat{G}}_n}-{\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right){\parallel }^2=\parallel Z_{{\widehat{G}}_n}-\mu +\mu -{\tilde{\mu }}_G\left(\omega \right)+{\tilde{\mu }}_G\left(\omega \right)-{\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right){\parallel }^2\hfill \\ \hfill & =\parallel Z_{{\widehat{G}}_n}-\mu {\parallel }^2+{\tilde{L}}_G\left(\omega \right)+\Phi \left(\omega \right),\hfill \end{array}$$ where the term $\parallel Z_{{\widehat{G}}_n}{-\mu \parallel }^2$ is unrelated to $\omega$, and (A21)$$\begin{array}{cc}\hfill \Phi \left(\omega \right)=& \parallel {\tilde{\mu }}_G\left(\omega \right)-{\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right){\parallel }^2+2{\left(Z_{{\widehat{G}}_n}-Z_G\right)}^{\top }\left\{\mu -{\tilde{\mu }}_G\left(\omega \right)\right\}\hfill \\ \hfill & +2e_G^{\top }\tilde{A}\left(\omega \right)\mu -2e_G^{\top }\tilde{P}\left(\omega \right)e_G+2{\left(Z_{{\widehat{G}}_n}-Z_G\right)}^{\top }\left\{{\tilde{\mu }}_G\left(\omega \right)-{\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)\right\}\hfill \\ \hfill & +2e_G^{\top }\tilde{P}\left(\omega \right)\left(Z_G-Z_{{\widehat{G}}_n}\right)+2{\left\{\mu -{\tilde{\mu }}_G\left(\omega \right)\right\}}^{\top }\left\{{\tilde{\mu }}_G\left(\omega \right)-{\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)\right\}.\hfill \end{array}$$

Considering (A13), (A15) is implied by (A22)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}{\displaystyle \frac{|\Phi (\omega \left)\right|}{{\tilde{R}}_G\left(\omega \right)}}=o_p\left(1\right).\end{array}$$

Using Cauchy–Schwarz inequality, Lemma A1, Lemma A3, and [45] to establish (A22), we only need to show (A23)$$\begin{array}{c}\hfill \underset{\omega \in W}{sup}{\displaystyle \frac{\parallel {\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)-{\tilde{\mu }}_G{\left(\omega \right)\parallel }^2}{{\tilde{R}}_G\left(\omega \right)}}=o_p\left(1\right).\end{array}$$

By Equation (A7) and Lemma A3, we observe (A24)$$\begin{array}{cc}\hfill \underset{\omega \in W}{sup}& {\displaystyle \frac{\parallel {\tilde{\mu }}_{{\widehat{G}}_n}\left(\omega \right)-{\tilde{\mu }}_G{\left(\omega \right)\parallel }^2}{{\tilde{R}}_G\left(\omega \right)}}=\underset{\omega \in W}{sup}{\displaystyle \frac{\parallel \tilde{P}\left(\omega \right)(Z_{{\widehat{G}}_n}-Z_G){\parallel }^2}{{\tilde{R}}_G\left(\omega \right)}}\hfill \\ \hfill & \le \underset{\omega \in W}{sup}{\displaystyle \frac{{\overline{\lambda }}^2\left\{\tilde{P}\left(\omega \right)\right\}{\parallel Z_{{\widehat{G}}_n}-Z_G\parallel }^2}{{\tilde{R}}_G\left(\omega \right)}}=o_p\left(1\right).\hfill \end{array}$$

Thus, we can obtain (A15). This concludes the proof. □

References

1. Engle, R.F.; Granger, C.W.J.; Rice, J.; Weiss, A. Semiparametric estimates of the relation between weather and electricity sales. J. Am. Stat. Assoc.; 1986; 81, pp. 310-320. [DOI: https://dx.doi.org/10.1080/01621459.1986.10478274]

2. Speckman, P. Kernel smoothing in partial linear models. J. R. Stat. Soc. Ser. B Stat. Methodol.; 1988; 50, pp. 413-436. [DOI: https://dx.doi.org/10.1111/j.2517-6161.1988.tb01738.x]

3. Heckman, N.E. Spline smoothing in a partly linear model. J. R. Stat. Soc. Ser. B Stat. Methodol.; 1986; 48, pp. 244-248. [DOI: https://dx.doi.org/10.1111/j.2517-6161.1986.tb01407.x]

4. Shi, J.; Lau, T. Empirical likelihood for partially linear models. J. Multivar. Anal.; 2000; 72, pp. 132-148. [DOI: https://dx.doi.org/10.1006/jmva.1999.1866]

5. Härdle, W.; Liang, H.; Gao, J. Partially Linear Models; Springer Science & Business Media: Berlin, Germany, 2000.

6. Claeskens, G.; Hjort, N.L. Model Selection and Model Averaging; Cambridge University Press: Cambridge, UK, 2008.

7. Hansen, B.E.; Racine, J.S. Jackknife model averaging. J. Econom.; 2012; 167, pp. 38-46. [DOI: https://dx.doi.org/10.1016/j.jeconom.2011.06.019]

8. Racine, J.S.; Li, Q.; Yu, D.; Zheng, L. Optimal model averaging of mixed-data kernel-weighted spline regressions. J. Bus. Econ. Stat.; 2022.in press [DOI: https://dx.doi.org/10.1080/07350015.2022.2118126]

9. Akaike, H. Statistical predictor identification. Ann. Inst. Statist. Math.; 1970; 22, pp. 203-217. [DOI: https://dx.doi.org/10.1007/BF02506337]

10. Schwarz, G. Estimating the dimension of a model. Ann. Statist.; 1978; 6, pp. 461-464. [DOI: https://dx.doi.org/10.1214/aos/1176344136]

11. Claeskens, G.; Hjort, N.L. The focused information criterion. J. Am. Stat. Assoc.; 2003; 98, pp. 900-916. [DOI: https://dx.doi.org/10.1198/016214503000000819]

12. Ni, X.; Zhang, H.; Zhang, D. Automatic model selection for partially linear models. J. Multivar. Anal.; 2009; 100, pp. 2100-2111. [DOI: https://dx.doi.org/10.1016/j.jmva.2009.06.009]

13. Raheem, S.E.; Ahmed, S.E.; Doksum, K.A. Absolute penalty and shrinkage estimation in partially linear models. Comput. Stat. Data Anal.; 2012; 56, pp. 874-891. [DOI: https://dx.doi.org/10.1016/j.csda.2011.09.021]

14. Xie, H.; Huang, J. SCAD-penalized regression in high-dimensional partially linear models. Ann. Statist.; 2009; 37, pp. 673-696. [DOI: https://dx.doi.org/10.1214/07-AOS580]

15. Peng, J.; Yang, Y. On improvability of model selection by model averaging. J. Econom.; 2022; 229, pp. 246-262. [DOI: https://dx.doi.org/10.1016/j.jeconom.2020.12.003]

16. Hoeting, J.A.; Madigan, D.; Raftery, A.E.; Volinsky, C.T. Bayesian model averaging: A tutorial. Statist. Sci.; 1999; 14, pp. 382-417.

17. Hansen, B.E. Least squares model averaging. Econometrica; 2007; 75, pp. 1175-1189. [DOI: https://dx.doi.org/10.1111/j.1468-0262.2007.00785.x]

18. Zhang, X.; Zou, G.; Carroll, R.J. Model averaging based on Kullback-Leibler distance. Stat. Sin.; 2015; 25, pp. 1583-1598. [DOI: https://dx.doi.org/10.5705/ss.2013.326] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/27761098]

19. Liu, Q.; Okui, R.; Yoshimura, A. Generalized least squares model averaging. Economet. Rev.; 2016; 35, pp. 1692-1752. [DOI: https://dx.doi.org/10.1080/07474938.2015.1092817]

20. Gao, Y.; Zhang, X.; Wang, S.; Zou, G. Model averaging based on leave-subject-out cross-validation. J. Econom.; 2016; 192, pp. 139-151. [DOI: https://dx.doi.org/10.1016/j.jeconom.2015.07.006]

21. Zhang, X.; Liu, C. Model averaging prediction by K-fold cross-validation. J. Econom.; 2022.in press

22. Lu, X.; Su, L. Jackknife model averaging for quantile regressions. J. Econom.; 2015; 188, pp. 40-58. [DOI: https://dx.doi.org/10.1016/j.jeconom.2014.11.005]

23. Zhang, X.; Wang, W. Optimal model averaging estimation for partially linear models. Stat. Sin.; 2019; 29, pp. 693-718. [DOI: https://dx.doi.org/10.5705/ss.202015.0392]

24. Zhu, R.; Wan, A.T.K.; Zhang, X.; Zou, G. A Mallows-type model averaging estimator for the varying-coefficient partially linear model. J. Am. Stat. Assoc.; 2019; 114, pp. 882-892. [DOI: https://dx.doi.org/10.1080/01621459.2018.1456936]

25. Xie, J.; Yan, X.; Tang, N. A model-averaging method for high-dimensional regression with missing responses at random. Stat. Sin.; 2021; 31, pp. 1005-1026. [DOI: https://dx.doi.org/10.5705/ss.202018.0297]

26. Wei, Y.; Wang, Q.; Liu, W. Model averaging for linear models with responses missing at random. Ann. Inst. Statist. Math.; 2021; 73, pp. 535-553. [DOI: https://dx.doi.org/10.1007/s10463-020-00759-y]

27. Zhang, X.; Chiou, J.; Ma, Y. Functional prediction through averaging estimated functional linear regression models. Biometrika; 2018; 105, pp. 945-962. [DOI: https://dx.doi.org/10.1093/biomet/asy041]

28. Zhang, X.; Ma, Y.; Carroll, R.J. MALMEM: Model averaging in linear measurement error models. J. R. Stat. Soc. Ser. B Stat. Methodol.; 2019; 81, pp. 763-779. [DOI: https://dx.doi.org/10.1111/rssb.12317]

29. Ando, T.; Li, K.C. A model-averaging approach for high-dimensional regression. J. Am. Stat. Assoc.; 2014; 109, pp. 254-265. [DOI: https://dx.doi.org/10.1080/01621459.2013.838168]

30. Ando, T.; Li, K.C. A weight-relaxed model averaging approach for high-dimensional generalized linear models. Ann. Statist.; 2017; 45, pp. 2654-2679. [DOI: https://dx.doi.org/10.1214/17-AOS1538]

31. Zeng, D.; Lin, D. Efficient estimation for the accelerated failure time model. J. Am. Stat. Assoc.; 2007; 102, pp. 1387-1396. [DOI: https://dx.doi.org/10.1198/016214507000001085]

32. Wang, H.J.; Wang, L. Locally weighted censored quantile regression. J. Am. Stat. Assoc.; 2009; 104, pp. 1117-1128. [DOI: https://dx.doi.org/10.1198/jasa.2009.tm08230]

33. Hjort, N.L.; Claeskens, G. Focused information criteria and model averaging for the Cox hazard regression model. J. Am. Stat. Assoc.; 2006; 101, pp. 1449-1464. [DOI: https://dx.doi.org/10.1198/016214506000000069]

34. Du, J.; Zhang, Z.; Xie, T. Focused information criterion and model averaging in censored quantile regression. Metrika; 2017; 80, pp. 547-570. [DOI: https://dx.doi.org/10.1007/s00184-017-0616-1]

35. Sun, Z.; Sun, L.; Lu, X.; Zhu, J.; Li, Y. Frequentist model averaging estimation for the censored partial linear quantile regression model. J. Statist. Plann. Inference; 2017; 189, pp. 1-15. [DOI: https://dx.doi.org/10.1016/j.jspi.2017.04.001]

36. Yan, X.; Wang, H.; Wang, W.; Xie, J.; Ren, Y.; Wang, X. Optimal model averaging forecasting in high-dimensional survival analysis. Int. J. Forecast.; 2021; 37, pp. 1147-1155. [DOI: https://dx.doi.org/10.1016/j.ijforecast.2020.12.004]

37. Liang, Z.; Chen, X.; Zhou, Y. Mallows model averaging estimation for linear regression model with right censored data. Acta Math. Appl. Sin. E.; 2022; 38, pp. 5-23. [DOI: https://dx.doi.org/10.1007/s10255-022-1054-z]

38. Koul, H.; Susarla, V.; Ryzin, J.V. Regression analysis with randomly right-censored data. Ann. Statist.; 1981; 9, pp. 1276-1288. [DOI: https://dx.doi.org/10.1214/aos/1176345644]

39. Xia, X. Model averaging prediction for nonparametric varying-coefficient models with B-spline smoothing. Stat. Pap.; 2021; 62, pp. 2885-2905. [DOI: https://dx.doi.org/10.1007/s00362-020-01218-9]

40. De Boor, C. A Practical Guide to Splines; Springer: New York, NY, USA, 2001.

41. Kaplan, E.L.; Meier, P. Nonparametric estimation from incomplete observations. J. Am. Stat. Assoc.; 1958; 53, pp. 457-481. [DOI: https://dx.doi.org/10.1080/01621459.1958.10501452]

42. Hu, G.; Cheng, W.; Zeng, J. Model averaging by jackknife criterion for varying-coefficient partially linear models. Comm. Statist. Theory Methods; 2020; 49, pp. 2671-2689. [DOI: https://dx.doi.org/10.1080/03610926.2019.1580736]

43. Turlach, B.A.; Weingessel, A.; Moler, C. Quadprog: Functions to Solve Quadratic Programming Problems. R Package Version 1.5-8. 2019; Available online: https://CRAN.R-project.org/package=quadprog (accessed on 16 December 2022).

44. Wei, Y.; Wang, Q. Cross-validation-based model averaging in linear models with response missing at random. Stat. Probab. Lett.; 2021; 171, 108990. [DOI: https://dx.doi.org/10.1016/j.spl.2020.108990]

45. Zhang, X.; Wan, A.T.K.; Zou, G. Model averaging by jackknife criterion in models with dependent data. J. Econom.; 2013; 174, pp. 82-94. [DOI: https://dx.doi.org/10.1016/j.jeconom.2013.01.004]

46. Wan, A.T.; Zhang, X.; Zou, G. Least squares model averaging by Mallows criterion. J. Econom.; 2010; 156, pp. 277-283. [DOI: https://dx.doi.org/10.1016/j.jeconom.2009.10.030]

47. Fan, J.; Ma, Y.; Dai, W. Nonparametric independence screening in sparse ultra-high-dimensional varying coefficient models. J. Am. Stat. Assoc.; 2014; 109, pp. 1270-1284. [DOI: https://dx.doi.org/10.1080/01621459.2013.879828] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25309009]

48. Bates, D.M.; Venables, W.N. Splines: Regression Spline Functions and Classes. R Package Version 3.6-1. 2019; Available online: https://CRAN.R-project.org/package=splines (accessed on 15 December 2022).

49. Therneau, T.M.; Lumley, T.; Elizabeth, A.; Cynthia, C. Survival: Survival Analysis. R Package Version 3.4-0. 2022; Available online: https://CRAN.R-project.org/package=survival (accessed on 15 December 2022).

50. Tibshirani, R. The lasso method for variable selection in the Cox model. Stat. Med.; 1997; 16, pp. 385-395. [DOI: https://dx.doi.org/10.1002/(SICI)1097-0258(19970228)16:4<385::AID-SIM380>3.0.CO;2-3]

51. Shows, J.H.; Lu, W.; Zhang, H.H. Sparse estimation and inference for censored median regression. J. Statist. Plann. Inference; 2010; 140, pp. 1903-1917. [DOI: https://dx.doi.org/10.1016/j.jspi.2010.01.043] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/20607110]

52. Diebold, F.X.; Mariano, R.S. Comparing predictive accuracy. J. Bus. Econ. Stat.; 1995; 13, pp. 253-263.

53. Rosenwald, A.; Wright, G.; Wiestner, A.; Chan, W.C.; Connors, J.M.; Campo, E.; Gascoyne, R.D.; Grogan, T.M.; Muller-Hermelink, H.K.; Smeland, E.B. et al. The proliferation gene expression signature is a quantitative integrator of oncogenic events that predicts survival in mantle cell lymphoma. Cancer Cell; 2003; 3, pp. 185-197. [DOI: https://dx.doi.org/10.1016/S1535-6108(03)00028-X] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/12620412]

54. Ma, S.; Du, P. Variable selection in partly linear regression model with diverging dimensions for right censored data. Stat. Sin.; 2012; 22, pp. 1003-1020. [DOI: https://dx.doi.org/10.5705/ss.2010.267]