Stability and $L_{2}$-penalty in Model Averaging

We assume that $S=\{z_{i}=(y_{i},x_{i}^{{}^{\prime}})^{{}^{\prime}}\in\mathcal{Z},i=1,...,n\}$ is a simple random sample from distribution $\mathcal{D}$, where $y_{i}$ is the $i$-th observation of response variable, and $x_{i}$ is the $i$-th observation of covariates. Let $z^{*}=(y^{*},x^{*^{\prime}})^{{}^{\prime}}$ be an observation from distribution $\mathcal{D}$ and independent of $S$.

In model averaging, $M$ approximating models are selected first in order to describe the relationship between response variable and covariates. We assume that the hypothesis spaces of $M$ approximating models are

where $x^{*}_{m}$ consists of some elements of $x^{*}$, and $\mathcal{F}_{m}$ is a given function set. For example, in MMA, to estimate $E(y^{*}|x^{*})$, we take

where $dim(\cdot)$ represents the dimension of the vector. For the $m$-th approximating model, a proper estimation method $A_{m}$ is selected, and $\hat{h}_{m}$, the estimator of $h_{m}$, is obtained based on $S$ and $A_{m}$. Then, the hypothesis space of model averaging is defined as follows:

where $W$ is a given weight space, and $H(\cdot)$ is a given function of weight vector and estimators of $M$ approximating models. In MMA, we take

as the model average estimator of $E(y^{*}|x^{*})$. An important problem with model averaging is the choice of model weights. Here, the estimator $\hat{w}$ of weight vector is obtained based on $S$ and a proper weight selection criterion $A(w)$ that makes $\hat{w}$ be optimal in a certain sense.

The selection of $A_{m},m=1,...,M$ and $A(w)$ is closely related to the definition of the loss function. Let $L[\hat{h}(x^{*},w),y^{*}]$ be a real value loss function which is defined in $\mathcal{H}\times\mathcal{Y}$, where $\mathcal{Y}$ is the value space of $y^{*}$. Then, the risk function is defined as follows:

which is an MSPE under the sample set S and weight vector $w$.

In this paper, we mainly discuss whether $F(\hat{w},S)$ can approximate smallest possible risk $inf_{w\in W}F(w,S)$. If $A(w)$ has such a property, we say that $A(w)$ is consistent. For fixed $m$, Shalev-Shwartz et al. (2010) defined the stability and consistency of $A_{m}$, and discussed their relationship. Obviously, for model averaging, we need to pay more attention to the stability and consistency of weight selection. We note that relevant conclusions of Shalev-Shwartz et al. (2010) cannot be directly applied to model averaging because $\mathcal{H}$ depends on $S$. Therefore, we extend the relevant definitions and conclusions to model averaging. The following is the definition of consistency:

In statistical learning theory, stability concerns how much the algorithm’s output varies when $S$ changes a little.“Leave-one-out(Loo)” and “Replace-one(Ro)” are two common tools used to evaluate stability. Loo considers the change in the algorithm’s output after removing an observation from $S$, and Ro considers such a change after replacing an observation in $S$ with an observation that is independent of S. Accordingly, the stability is called Loo stability and Ro stability respectively. Here we will give the formal definitions of Loo stability and Ro stability. To this end, we first give the definition of algorithm symmetry:

Now let $S^{-i}$ be the sample set after removing the $i$-th observation from $S$, $\hat{h}^{-i}_{m}$ be the estimator of $h_{m}$ based on $S^{-i}$ and $A_{m}$, $\hat{w}^{-i}$ be the estimator of weight vector based on $S^{-i}$ and $A(w)$, and $F(w,S^{-i})=E_{z^{*}}\big{\{}L[\hat{h}^{-i}(x^{*},w),y^{*}]\big{\}}$, where $\hat{h}^{-i}(x^{*},w)=H[w,\hat{h}^{-i}_{1}(x^{*}_{1}),...,\hat{h}^{-i}_{M}(x^{*}_{M})]$. We define Loo stability as follows:

Let $S^{i}$ be the sample set S with the $i$-th observation replaced by $z_{i}^{*}=(y_{i}^{*},x_{i}^{*^{\prime}})^{{}^{\prime}}$, $\hat{h}^{i}_{m}$ be the estimator of $h_{m}$ based on $S^{i}$ and $A_{m}$, and $\hat{w}^{i}$ be the estimator of weight vector based on $S^{i}$ and $A(w)$, where $z_{i}^{*}$ is from distribution $\mathcal{D}$ and independent of $S$. Let $F(w,S^{i})=E_{z^{*}}\big{\{}L[\hat{h}^{i}(x^{*},w),y^{*}]\big{\}}$, then

where $\hat{h}^{i}(x^{*},w)=H[w,\hat{h}^{i}_{1}(x^{*}_{1}),...,\hat{h}^{i}_{M}(x^{*}_{M})]$. Therefore, we define Ro stability as follows:

Before discussing the relationship between stability and consistency, we give the definitions of AERM and generalization. The empirical risk function is defined as follows:

Vapnik (1998) proved some theoretical properties of the empirical risk minimization principle. However, when the sample size is small, the empirical risk minimizer tends to produce over-fitting phenomenon. Therefore, the structural risk minimization principle is proposed in Vapnik (1998), and the method to satisfy this principle is usually an AERM. Shalev-Shwartz et al. (2010) also discussed the deficiency of the empirical risk minimization principle and the importance of AERM.

In statistical learning theory, generalization refers to the performance of the concept learned by models on unknown samples. It can be seen from Definition 7 that the generalization of $A(w)$ describes the difference between using $\hat{w}$ to fit the training set $S$ and predict the unknown sample.

Note that for any $i\in\{1,...,n\}$,

so we give the following theorem to illustrate the relationship between Loo stability and Ro stability:

Shalev-Shwartz et al. (2010) emphasized that Ro stability and Loo stability are in general incomparable notions, but Theorem 2.1 shows that they are closely related.

By definitions of generalization and Ro stability, we have

and then we give the following theorem to illustrate the equivalence of Ro stability and generalization:

Theorem 2.2 shows that stability is an important property of weight selection criterion, which can ensure that the corresponding estimator has good generalization performance.

Let $\hat{w}^{*}\in W$ satisfy $F(\hat{w}^{*},S)=inf_{w\in W}F(w,S)$. Note that

so we give the following theorem to illustrate the relationship between stability and consistency:

Since $\hat{w}^{*}$ and $\mathcal{H}$ depend on $S$, unlike Lemma 15 in Shalev-Shwartz et al. (2010), Theorem 2.3 requires $\hat{w}^{*}$ to generalize with rate $\epsilon_{n}$. In the next section, we will propose a $L_{2}$-penalty model averaging method and prove that it has stability and consistency under certain reasonable conditions.

We assume that the response variable $y_{i}$ and covariates $x_{i}=(x_{1i},x_{2i},...)$ satisfy the following data generating process:

and M approximating models are given by

where $b_{mi}=\mu_{i}-\sum_{k=1}^{k_{m}}x_{m(k)i}\theta_{m,(k)}$ is the approximating error of the $m$-th approximating model. We assume that the $M$-th approximating model contains all the considered covariates. To simplify notation, we let $x_{i}=(x_{(1)i},...,x_{(k_{M})i})$ throughout the rest of this article, where $x_{(k)i}=x_{M(k)i}$.

Let $y=(y_{1},...,y_{n})^{{}^{\prime}}$, $\mu=(\mu_{1},...,\mu_{n})^{{}^{\prime}}$, $e=(e_{1},...,e_{n})^{{}^{\prime}}$, $b_{m}=(b_{m1},...,b_{mn})^{{}^{\prime}}$. Then, the corresponding matrix form of the true model is $y=X_{m}\theta_{m}+b_{m}+e$, where $\theta_{m}=(\theta_{m,(1)},...,\theta_{m,(k_{m})})^{{}^{\prime}}$, and $X_{m}$ is the design matrix of the $m$-th approximating model. When there is an approximating model such that $b_{m}=\textbf{0}$, it indicates that the true model is included in the $m$-th approximating model, i.e. the model is correctly specified. Unlike Hansen (2007) , Wan et al. (2010) and Hansen and Racine (2012), we do not require that the infimum of $R_{n}(w)$ (this is defined in section 3.2) tends to infinity, and therefore allow that the model is correctly specified.

Let $\pi_{m}\in R^{K\times k_{m}}$ be the variable selection matrix satisfying $X_{M}\pi_{m}=X_{m}$ and $\pi_{m}^{{}^{\prime}}\pi_{m}=I_{k_{m}}$, $m=1,...,M$. Then, the hypothesis spaces of M approximating models are

The least squares estimator of $\theta_{m}$ is $\hat{\theta}_{m}=(X_{m}^{{}^{\prime}}X_{m})^{-1}X_{m}^{{}^{\prime}}y$, $m=1,...,M$.

Let $P_{m}=X_{m}(X_{m}^{{}^{\prime}}X_{m})^{-1}X_{m}^{{}^{\prime}}$, $P(w)=\sum_{m=1}^{M}w_{m}P_{m}$, $L_{n}(w)=\|\mu-P(w)y\|_{2}^{2}$, and $R_{n}(w)=E_{e}[L_{n}(w)]$. When $\sigma_{i}^{2}\equiv\sigma^{2}$, Hansen (2007) and Wan et al. (2010) used Mallows criterion $C_{n}(w)=\|y-\hat{\Omega}w\|_{2}^{2}+2\sigma^{2}w^{{}^{\prime}}\kappa$ to select a model weight vector from hypothesis space $\mathcal{H}^{0}=\big{\{}x^{*^{\prime}}\hat{\theta}(w),w\in W^{0}\big{\}}$ and proved that the estimator of weight vector asymptotically minimizes $L_{n}(w)$, where $\hat{\Omega}=(P_{1}y,...,P_{M}y)$, $\kappa=(k_{1},...,k_{M})^{{}^{\prime}}$, and $\hat{\theta}(w)=\sum_{m=1}^{M}w_{m}\pi_{m}\hat{\theta}_{m}$. Hansen and Racine (2012) used Jackknife criterion $J_{n}(w)=\|y-\bar{\Omega}w\|_{2}^{2}$ to select a model weight vector from hypothesis space $\mathcal{H}^{0}$ and proved that the estimator of weight vector asymptotically minimizes $L_{n}(w)$ and $R_{n}(w)$, where $\bar{\Omega}=\big{[}y-D_{1}(I-P_{1})y,...,y-D_{M}(I-P_{M})y\big{]}$ with $D_{m}=diag[(1-h_{ii}^{m})^{-1}]$ and $h_{ii}^{m}=x_{i}^{{}^{\prime}}\pi_{m}(X_{m}^{{}^{\prime}}X_{m})^{-1}\pi_{m}^{{}^{\prime}}x_{i},i=1,...,n$.

Different from Hansen (2007), Wan et al. (2010) and Hansen and Racine (2012), we focus on whether the model averaging estimator can asymptotically minimize MSPE when the model weights are not restricted. Let $\hat{\gamma}=(x^{*^{\prime}}\pi_{1}\hat{\theta}_{1},...,x^{*^{\prime}}\pi_{M}\hat{\theta}_{M})$. Then, the risk function and the empirical risk function are defined as:

and

respectively. Since Hansen (2007), Wan et al. (2010) and Hansen and Racine (2012) restrict $w\in W^{0}$, the corresponding estimators of weight vector do not necessarily asymptotically minimize $F(w,S)$ on $R^{M}$. An intuitive way that enables the estimator of weight vector to asymptotically minimize $F(w,S)$ on $R^{M}$ is to remove the restriction $w\in W^{0}$ directly.

Let $\hat{P}$ and $\bar{P}$ be the orthogonal matrices satisfying $\hat{P}^{{}^{\prime}}\hat{\Omega}^{{}^{\prime}}\hat{\Omega}\hat{P}=diag(\hat{\zeta}_{1},...,\hat{\zeta}_{M})$, $\bar{P}^{{}^{\prime}}\bar{\Omega}^{{}^{\prime}}\bar{\Omega}\bar{P}=diag(\bar{\zeta}_{1},...,\bar{\zeta}_{M})$, where $\hat{\zeta}_{1}\leq...\leq\hat{\zeta}_{M}$ and $\bar{\zeta}_{1}\leq...\leq\bar{\zeta}_{M}$ are the eigenvalues of $\hat{\Omega}^{{}^{\prime}}\hat{\Omega}$ and $\bar{\Omega}^{{}^{\prime}}\bar{\Omega}$ respectively. We assume that $E_{z^{*}}(\hat{\gamma}^{{}^{\prime}}\hat{\gamma})$, $\hat{\Omega}^{{}^{\prime}}\hat{\Omega}$ and $\bar{\Omega}^{{}^{\prime}}\bar{\Omega}$ are invertible( this is reasonable under Assumption 3 ), then

and

From this, we can see that, in order to satisfy the consistency, $\hat{w}^{0}$ and $\bar{w}^{0}$ should be good estimators of $\hat{w}^{*}$. However, when approximating models are highly correlated, the minimum eigenvalues of $\hat{\Omega}^{{}^{\prime}}\hat{\Omega}$ and $\bar{\Omega}^{{}^{\prime}}\bar{\Omega}$ may be small so that $\|\hat{w}^{0}\|_{2}^{2}=\sum_{m=1}^{M}\frac{a_{m}^{2}}{\hat{\zeta}_{m}^{2}}\geq\frac{a_{1}^{2}}{\hat{\zeta}_{1}^{2}}$ and $\|\bar{w}^{0}\|_{2}^{2}=\sum_{m=1}^{M}\frac{b_{m}^{2}}{\bar{\zeta}_{m}^{2}}\geq\frac{b_{1}^{2}}{\bar{\zeta}_{1}^{2}}$ are too large, which usually result in extreme weights, where $(a_{1},a_{2},...,a_{M})^{{}^{\prime}}=\hat{P}^{{}^{\prime}}\hat{\Omega}y$, $(b_{1},b_{2},...,b_{M})^{{}^{\prime}}=\bar{P}^{{}^{\prime}}\bar{\Omega}y$. Therefore, similar to ridge regression in Hoerl and Kennard (1970), we make the following correction to $C_{n}(w)$ and $J_{n}(w)$:

where $\lambda_{n}>0$ is a tuning parameter. The above corrections are actually $L_{2}$-penalty for weight vector. Let $\hat{Z}=(\hat{\Omega}^{{}^{\prime}}\hat{\Omega}+\lambda_{n}I)^{-1}\hat{\Omega}^{{}^{\prime}}\hat{\Omega}$ , $\bar{Z}=(\bar{\Omega}^{{}^{\prime}}\bar{\Omega}+\lambda_{n}I)^{-1}\bar{\Omega}^{{}^{\prime}}\bar{\Omega}$. Then

and

In the next subsection, we discuss the theoretical properties of $C(w,S)$ and $J(w,S)$.

Let $\lambda_{max}(\cdot)$ and $\lambda_{min}(\cdot)$ be the maximum and minimum eigenvalues of a square matrix respectively, and $\chi$ be the value space of K covariates, where $K=k_{M}$. In order to discuss the stability and consistency of the proposed method, we need the following assumptions:

Assumption 1 is mild and similar conditions can be found in Zou and Zhang (2009) and Zhao et al. (2020). From $X_{m}^{{}^{\prime}}X_{m}=\pi_{m}^{{}^{\prime}}X_{M}^{{}^{\prime}}X_{M}\pi_{m}$, we see that, under Assumption 1, for any m $\in\{1,...,M\}$, we have

Shalev-Shwartz et al. (2010) assumed that the loss function is bounded, which is usually not satisfied in traditional regression analysis. We replace this assumption with Assumption 2. Tong and Wu (2017) assumed that $\chi\times\mathcal{Y}$ is a compact subset of $R^{K+1}$, under which Assumption 2 is obviously true. Assumption 3 requires the minimum eigenvalue of $\hat{\Omega}^{{}^{\prime}}\hat{\Omega}$ to have a lower bound away from 0 and the maximum eigenvalue to have an order $O(nM)$ a.s.. A similar assumption is used in Liao and Zou (2020). Lemma 3 guarantees the rationality of the assumptions about the eigenvalues of $E_{z^{*}}(\hat{\gamma}^{{}^{\prime}}\hat{\gamma})$ and $\bar{\Omega}^{{}^{\prime}}\bar{\Omega}$. Assumption 4 is a mild assumption for tuning parameter in order to avoid excessive penalty. In Section 3.4, we look for tuning parameters only from $[0,M\log n]$.

Let $\hat{V}(\lambda_{n})=\|\hat{Z}\hat{w}^{0}-\hat{Z}\hat{w}^{*}\|_{2}^{2}$, $\hat{B}(\lambda_{n})=\|\hat{Z}\hat{w}^{*}-\hat{w}^{*}\|_{2}^{2}$, $\bar{V}(\lambda_{n})=\|\bar{Z}\bar{w}^{0}-\bar{Z}\hat{w}^{*}\|_{2}^{2}$, and $\bar{B}(\lambda_{n})=\|\bar{Z}\hat{w}^{*}-\hat{w}^{*}\|_{2}^{2}$. We define