超高维数据下半参数工具变量模型的双惩罚正则估计

Double Penalized Regularization Estimation for Semiparametric Instrumental Variable Models with Ultrahigh Dimensional Data

摘要: 在超高维数据下，考虑一类含内生性协变量的半参数工具变量模型的估计。结合惩罚估计技术以及SIS特征筛选方法，提出了一种识别有效工具变量的双惩罚正则估计方法，并得到了模型参数分量和非参数分量的工具变量调整估计。在一些正则条件下，证明了参数分量和非参数分量的估计均是相合的。最后通过一些数值仿真模拟研究了所提方法的有限样本性质。

Abstract: The estimation of a semiparametric instrumental variable model with endogenous covariates is considered under ultra-high dimensional data. Combined with penalty estimation technology and SIS feature screening method, a double penalized regularization estimation method for identifying optimal instrumental variables is proposed, and an optimal instrumental variable adjusted esti-mator of model parameter and nonparametric components are obtained. Under some regular con-ditions, it is proved that the estimators of parametric and nonparametric components are con-sistent. Finally, the finite sample properties of the proposed method are studied by some numerical simulations.

1. 引言

设 $\left\{X_i,Z_i,U_i,Y_i\right\}$， $i=1,\cdots ,n$ 表示一组容量为n的样本，其中 $X_i$， $Z_i$ 和 $U_i$ 是协变量， $Y_i$ 是响应变量。本文考虑如下一类带变系数的半参数模型：

$Y_i=X_i^{\text{T}}\beta +Z_i^{\text{T}}\theta \left(U_i\right)+{\epsilon }_i$, $i=1,\cdots ,n$ (1)

其中 $\beta ={\left({\beta }_1,\cdots ,{\beta }_p\right)}^{\text{T}}$ 是一个p维的未知参数向量， $\theta \left(U_i\right)={\left({\theta }_1\left(U_i\right),\cdots ,{\theta }_k\left(U_i\right)\right)}^{\text{T}}$ 是一个k维的未知非参数函数系数向量， ${\epsilon }_i$ 表示模型误差。模型(1)同时含有参数分量和非参数函数分量，因此在对实际问题统计建模过程中具有较好的灵活性和广泛的适应性。该模型近来引起了许多学者的广泛关注，并且对模型中的参数和非参数分量提出了大量的估计方法(参见文 [1] [2] [3])。但是以上文献是在假定所有协变量均为外生性协变量的前提下对模型的估计问题进行讨论，而对高维数据进行统建模中假定所有协变量均为外生性协变量是不切实际的。如果某些协变量具有内生性，上述文献提出的估计方法给出的估计则不是相合的，而是会产生一定的内生性偏差。

对于具有内生性变量的模型，工具变量调整技术则可以提供一个消除内生性偏差的方法，进而得到模型的相合估计。由于半参数模型具有较好的灵活性和广泛的适应性，近年来对含内生性变量的半参数模型的统计推断也越来越受到广泛的关注。例如，Cai和Xiong [4] 对一类含内生协变量的变系数部分线性模型的估计问题，提出了一种基于工具变量调整的三阶段估计方法。Yuan等 [5] 对一类含内生协变量的变系数部分线性模型的变量选择问题，提出了一种基于工具变量调整的变量选择方法。Yang等 [6] 则对一类含内生性变量的部分线性单指标模型的经验似然统计推断问题，提出了一种工具变量调整的经验似然估计方法。

注意到上述文献所提出工具变量调整方法均是在工具变量的维数为固定的低维数据下进行的，对超高维数据下的半参数工具变量模型的估计问题，上述给出的基于工具变量调整的估计方法将不能直接使用。为此，本文讨论当某些协变量为内生性协变量并且工具变量为超高维数据的情况下，研究模(1)的估计问题。具体地，假定模型(1)中的 $Z_i$ 和 $U_i$ 为外生性协变量，而 $X_i$ 为内生性协变量，并满足：

$X_i=\Gamma {\xi }_i+e_i$, $i=1,2,\cdots ,n$ (2)

其中 ${\xi }_i$ 是一个由工具变量构成的q维向量， $\Gamma$ 是一个 $p\times q$ 维的未知参数矩阵。 $e_i$ 是模型误差。注意到 $Z_i$ 和 $U_i$ 为外生性协变量，而 $X_i$ 为内生性协变量，那么本文可假定 $Z_i$ 和 $U_i$ 直接影响响应变量 $Y_i$，这表明 $Z_i$ 和 $U_i$ 均与 $X_i$ 相互独立。另外，在接下来的估计过程中假定协变量 $X_i$ 的维数p是固定的，而工具变量 ${\xi }_i$ 为超高维工具变量，即 ${\xi }_i$ 的维数q远远大于样本容量n。进一步，本文假设未知参数矩阵 $\Gamma$ 是稀疏的，即 $\Gamma$ 中只有一小部分参数是非零的，这也意味着只有小部分工具变量是有效工具变量。

本文的目标是识别并估计 $\Gamma$ 中的非零参数，进而识别有效工具变量。并利用识别出的有效工具变量以及工具变量调整技术对模型(1)中的参数和非参数分量给出一个估计过程。具体地，首先结合内生协变量 $X_i$ 和工具变量 ${\xi }_i$ 的相关关系构建一个辅助的回归模型，基于构建的辅助回归模型，并利用SIS (Sure Independence Screening)特征筛选方法提出了一种有效工具变量的识别方法。值得一提的是，注意到 $X_i$ 是一个多维向量以及 $\Gamma$ 是一个超高维矩阵，因此经典的SIS方法(参见文 [7])将不能再直接应用。本文提出了一种基于辅助回归模型的SIS特征筛选方法，并给出了一种基于双惩罚的迭代算法，改进并推广了已有的SIS变量筛选过程。其次，利用识别出的有效工具变量并采用工具变量调整技术，对模型(1)中的参数和非参数分量给出了一种基于工具变量调整的估计方法，并在一些正则条件下证明了所得估计是渐近相合估计量。

2. 估计过程及迭代算法

注意到 $X_i$ 为内生性协变量， ${\xi }_i$ 为对应的工具变量，所以有 $E\left({\epsilon }_i|X_i\right)\ne 0$ 和 $E\left({\epsilon }_i|{\xi }_i\right)=0$。由第1节中的讨论可知外生性协变量 $Z_i$ 和 $U_i$ 均与内生性协变量 $X_i$ 相互独立。那么结合 $X_i=\Gamma {\xi }_i+e_i$，可以证明工具变量 ${\xi }_i$ 也独立于外生性协变量 $Z_i$ 和 $U_i$。记 ${\xi }_{is}$ 为工具变量 ${\xi }_i={\left({\xi }_{i1},\cdots ,{\xi }_{iq}\right)}^{\text{T}}$ 的第s个分量，那经计算可得：

$Cov\left(Y_i,{\xi }_{is}\right)={\displaystyle \underset{j=1}{\overset{p}{\sum }}Cov\left(X_{ij},{\xi }_{is}\right){\beta }_j}+{\displaystyle \underset{j=1}{\overset{k}{\sum }}Cov\left(Z_{ij},{\xi }_{is}\right){\theta }_j}\left(U_i\right)+Cov\left({\epsilon }_i,{\xi }_{is}\right)={\displaystyle \underset{j=1}{\overset{p}{\sum }}Cov\left(X_{ij},{\xi }_{is}\right){\beta }_j}$ (3)

另外，Bound等 [8] 以及Didelez等 [9] 指出，当 ${\xi }_{is}$ 是有效工具变量时， ${\xi }_{is}$ 应与内生性变量 $X_i$ 高度相关。即如果对所有 $j\in \left\{1,\cdots ,p\right\}$ 有 $Cov\left(X_{ij},{\xi }_{is}\right)=0$，那么 ${\xi }_{is}$ 是一个无效的工具变量。另外注意到 $\beta \ne 0$，因此(3)式表明如果 $Cov\left(Y_i,{\xi }_{is}\right)\ne 0$ 显著成立，那么 ${\xi }_{is}$ 则是一个有效工具变量。基于该结论并利用SIS特征筛选方法，本文可以根据响应变量 $Y_i$ 和工具变量 ${\xi }_{is}$， $s\in \left\{1,\cdots ,q\right\}$ 之间的边际相关性选择有效工具变量。具体地，首先定义响应变量 $Y_i$ 和工具变量 ${\xi }_{is}$， $s\in \left\{1,\cdots ,q\right\}$ 的一个辅助回归模型如下：

$Y_i={\omega }_1{\xi }_{i1}+{\omega }_2{\xi }_{i2}+\cdots +{\omega }_q{\xi }_{iq}+v_i$, $i=1,\cdots ,n$ (4)

其中 $v_i$ 是模型误差满足 $E\left(v_i|{\xi }_{i1},\cdots ,{\xi }_{iq}\right)=0$。对于任给定的 $s\in \left\{1,\cdots ,q\right\}$，结合模型(4)并利用类似张等 [10] 的估计方法，定义关于分量的边际回归模型为 $Y_i={\omega }_s{\xi }_{is}+v_i$，并且 ${\omega }_s$ 的边际最小二乘估计可通过最小化以下目标函数给出。

${\displaystyle \underset{i=1}{\overset{n}{\sum }}Q\left(Y_i,{\omega }_s{\xi }_{is}\right)}\equiv {\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(Y_i-{\omega }_s{\xi }_{is}\right)}^2}$ (5)

通过最小化(5)式,可以得到边际最小二乘估计量为

${\hat{\omega }}_s^M={\left({\displaystyle \underset{i=1}{\overset{n}{\sum }}{\xi }_{is}^2}\right)}^{-1}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\xi }_{is}{Y}_i}$ (6)

结合(4)和(6)式，则可识别出有效工具变量的子集如下：

${\hat{M}}_{v_n}=\left\{1\le s\le q:\left|{\hat{\omega }}_s^M\right|\ge v_n\right\}$ (7)

其中 $v_n$ 是预先给定的一个门限值。另外，利用类似张等 [10] 的讨论可知，基于边际回归所识别出的有效工具变量可能会存在一定的选择性偏差。比如，某一有效工具变量与内生性变量 $X_i$ 的某些分量存在高度相关，但该工具变量可能与响应变量 $Y_i$ 相关性不显著，进而导致无法识别出该有效工具变量。为了克服这个问题，本文提出了一个新的基于双惩罚的迭代算法。在如下迭代算法中，我们还利用识别出的有效工具变量提出了内生性变量的调整过程，并给出了基于有效工具变量调整的模型参数分量 $\beta$ 和非参数分量 $\theta \left(u\right)$ 的估计过程。具体算法如下：在此算法中，本文首先预先确定一个稀疏参数大小 $d=O\left(n/\ln n\right)$。然后，算法如下。

第一阶段：有效工具变量识别算法

第1步：对所有的 $s\in \left\{1,\cdots ,q\right\}$，由(7)式识别出一个包含 $S_1$ 个有效工具变量的子集 ${\mathcal{A}}_1$，其中 $S_1=\lfloor 2d/3\rfloor$ 以及 $d=n/\ln n$。

第2步：记 ${\xi }_{{\mathcal{A}}_1}$ 为集合 ${\mathcal{A}}_1$ 对应的有效工具变量，并且模型(2)中的对应系数矩阵则为 ${\Gamma }_{{\mathcal{A}}_1}$。然后基于模型(2)，通过最小化如下惩罚目标函数来选择 ${\mathcal{A}}_1$ 的子集 ${\mathcal{M}}_1$ ：

$\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\Vert X_i-{\Gamma }_{{\mathcal{A}}_1}{Z}_{i,{\mathcal{A}}_1}\Vert }_2^2}+p_{\lambda }\left({\Gamma }_{{\mathcal{A}}_1}\right)$

其中 $p_{\lambda }(\cdot )$ 是一个惩罚函数。在实际应用中，Lasso惩罚，SCAD惩罚以及MCP惩罚等惩罚函数均可以使用。令 ${\hat{\Gamma }}_{{\mathcal{A}}_1}$ 为 ${\Gamma }_{{\mathcal{A}}_1}$ 的惩罚估计矩阵，有 ${\mathcal{M}}_1$ 是 ${\hat{\Gamma }}_{{\mathcal{A}}_1}$ 的非零列对应的指标集。

第3步：对任给定的 $s\in {\mathcal{M}}_1^c=\left\{1,\cdots ,q\right\}\backslash {\mathcal{M}}_1$，基于模型(4)，最小化如下惩罚目标函数：

${\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left[Y_i-\left({\xi }_{i,{\mathcal{M}}_1}^{\text{T}}{\omega }_{{\mathcal{M}}_1}+{\xi }_{is}{\omega }_s\right)\right]}^2}+p_{\lambda }(\; \omega \; s\; )$

对任意 $s\in {\mathcal{M}}_1^c$，令 ${\omega }_s$ 的估计量为 ${\hat{\omega }}_s^{\left(2\right)}$，那么通过排序 $\left\{\left|{\hat{\omega }}_s^{\left(2\right)}\right|:s\in {\mathcal{M}}_1^c\right\}$ 选择一个大小为 $k_2=\min \left(5,d-\left|{\mathcal{M}}_1\right|\right)$ 的指数子集 ${\mathcal{A}}_2$。进一步，用第2步的算法从 ${\mathcal{M}}_1\cup {\mathcal{A}}_2$ 中选择一个子集 ${\mathcal{M}}_2$。

第4步：对第3步进行迭代直到有 ${\mathcal{M}}_l\ge d$ 或 ${\mathcal{M}}_l={\mathcal{M}}_{l-1}$，所选的最终子集表示为 ${\mathcal{M}}_{final}$。注意 ${\mathcal{M}}_{final}$ 满足 $\left|{\mathcal{M}}_{final}\right|\le d<n$，因此，提出的迭代算法可实现降低维数的效果。

第二阶段：基于有效工具变量的内生性协变量的工具变量调整算法

首先用 ${\xi }^{*}$ 表示 ${\mathcal{M}}_{final}$ 集合中对相应的工具变量。那么 ${\xi }^{*}$ 则是识别出的有效工具变量组成的一个向量。另外，从模型(2)中可知 $E\left\{X{\xi }^{\text{T}}\right\}=\Gamma E\left\{\xi {\xi }^{\text{T}}\right\}$，那么 $\Gamma$ 中有效工具变量对应的参数矩阵的矩估计量可表示为：

${\hat{\Gamma }}^{\ast }=\left({\displaystyle \underset{i=1}{\overset{n}{\sum }}{X}_i{\xi }_i^{*\text{T}}}\right){\left({\displaystyle \underset{i=1}{\overset{n}{\sum }}{\xi }_i^{*}{\xi }_i^{*\text{T}}}\right)}^{-1}$

因此，基于有效工具变量调整的协变量可定义为 $X^{\ast }={\hat{\Gamma }}^{\ast }{\xi }^{\ast }$。

第三阶段：模型参数和非参数分量的估计

记 $B\left(u\right)={\left(B_1\left(u\right),\cdots ,B_L\left(u\right)\right)}^{\text{T}}$ 为表示阶数为m的B样条基数，其中 $L=\kappa +m$， $\kappa$ 为内部结点的个数。那么 ${\theta }_s\left(u\right)$ 可以渐近表示为 ${\theta }_s\left(u\right)\approx B{\left(u\right)}^{\text{T}}{\gamma }_s$，其中 ${\gamma }_s={\left({\gamma }_{s1},\cdots ,{\gamma }_{sL}\right)}^{\text{T}}$ 是基函数系数向量。记 $W_i=I_k\otimes B\left(U_i\right)\cdot Z_i$ 以及 $\gamma ={\left({\gamma }_1^{\text{T}},\cdots ,{\gamma }_k^{\text{T}}\right)}^{\text{T}}$，那么通过最小化如下目标函数来给出 $\beta$ 和 $\gamma$ 的估计。

$M_n\left(\beta ,\gamma \right)=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(Y_i-X_i^{*\text{T}}\beta -W_i^{\text{T}}\gamma \right)}^2}$ (8)

记 $\hat{\beta }$ 和 $\hat{\gamma }$ 为最小化(8)式的解，那么 $\hat{\beta }$ 则是基于工具变量调整的参数分量 $\beta$ 估计量，并且非参数分量 ${\theta }_s\left(u\right)$ 的估计可以定义为 ${\hat{\theta }}_s\left(u\right)\approx B{\left(u\right)}^{\text{T}}{\hat{\gamma }}_s$。

另外在上述算法过程中对非参数函数系进行样条展开时，内部结点取为等距结点，并且内部结点个数 $\kappa$ 通过最小化如下交叉验证得分函数进行得到。

$CV\left(\kappa \right)=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left\{Y_i-X_i^{\ast \text{T}}{\hat{\beta }}_{\left[i\right]}-Z_i^{\text{T}}{\hat{\theta }}_{\left[i\right]}\left(U_i\right)\right\}}^2}$

其中 ${\hat{\beta }}_{\left[i\right]}$ 和 ${\hat{\theta }}_{\left[i\right]}\left(u\right)$ 表示删除第i个样本个体后，基于(8)式得到的 $\beta$ 和 $\theta \left(u\right)$ 估计量。

3. 渐近性质

令 $S\left(Y_i,{\omega }_s{\xi }_{is}\right)=\partial Q\left(Y_i,{\omega }_s{\xi }_{is}\right)/\partial {\omega }_s$，其中 $Q\left(.,.\right)$ 由式(5)定义，那么由(6)式定义的估计量 ${\hat{\omega }}_s^M$ 也是估计方程 ${\displaystyle {\sum }_{i=1}^{n}S\left(Y_i,{\omega }_s{\xi }_{is}\right)}=0$ 的解。另外对每一个 $s\in \left\{1,\cdots ,q\right\}$，定义 ${\omega }_s^M$ 为 $E\left\{S\left(Y_i,{\omega }_s{\xi }_{is}\right)\right\}=0$ 的解。那么结合(5)式并简单计算可得

${\omega }_s^M=E\left({\xi }_{s}Y\right)/E\left({\xi }_s^2\right)$ (9)

进一步，记 ${\mathcal{M}}_{\ast }=\left\{1\le s\le q:{\omega }_s\ne 0\right\}$ 表示所有有效工具变量构成的指标集，假设 $\left|{\mathcal{M}}_{\ast }\right|=m<q$，即有效工具变量只是所有工具变量一个子集。此外令 $\mathcal{B}=\left\{{\omega }_s:\left|{\omega }_s\right|\le B\right\}$ 其中B为一个给定的常数。为了给出有效工具变量识别以及所得估计量的渐近性质，首先给出如下的一些正则性条件。

C1：记 $I\left({\omega }_s^M\right)=E\left\{S\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}$，那么假定 $I{\left({\omega }_s^M\right)}_{\mathcal{B}}={\sup }_{{\omega }_s^M\in \mathcal{B}}I\left({\omega }_s^M\right)$ 存在上界 $c_0$，其中 $c_0$ 表示某一正常数。

C2：对 $\mathcal{B}$ 中任意的 ${\omega }_s$ 和 ${{\omega }^{\prime }}_s$，函数 $Q\left(\cdot ,\cdot \right)$ 满足以下利普席茨条件：

$\left|Q\left(Y_i,{\omega }_s{\xi }_{is}\right)-Q\left(Y_i,{{\omega }^{\prime }}_s{\xi }_{is}\right)\right|\le c_1\left|{\xi }_{is}\right|\left|{\omega }_s-{{\omega }^{\prime }}_k\right|$

其中 $c_1$ 表示某一正常数。

C3：对所有的 ${\omega }_s\in \mathcal{B}$，函数 $Q\left(Y_i,{\omega }_s{\xi }_{is}\right)$ 关于 ${\omega }_s$ 是凸函数且满足：

$\left|E\left\{Q\left(Y_i,{\omega }_s{\xi }_{is}\right)-Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\right|\ge c_2{\left|{\omega }_s-{\omega }_s^M\right|}^2$

其中 $c_2$ 表示某一正常数。

C4：对任意的 $s\in \left\{1,\cdots ,q\right\}$，概率趋于1， $\left|{\xi }_{is}\right|\le c_3$ 和 $\left|Y_i\right|\le c_3$ 以概率趋于1成立，其中 $c_3$ 表示某一正常数。

C5：存在一个正常数 $c_4$，使得 ${\sup }_{s\in {\mathcal{M}}_{\ast }}\left|{\omega }_s^M\right|\le c_4$ 成立。此外，对任给定的 $s\in {\mathcal{M}}_{\ast }$，存在一个正常数 $c_5$ 有：

$\left|E\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)-Q\left(Y_i,0\right)\right\}\right|\ge c_5{n}^{-k}$

其中 $k\in \left(0,1/2\right)$。

在这些正则性条件下，如下定理1给出了边际估计量 ${\hat{\omega }}_s^M\left(s=1,\cdots ,q\right)$ 的一致收敛性。

定理1. 假定正则性条件C1~C5成立。那对任给定的 $c_6>0$，存在一个正常数 $c_7$ 有：

$P\left\{\underset{1\le \text{s}\le q}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge c_6{n}^{-k}\right\}\le q\exp \left(-c_7{n}^{1-2k}\right)$

其中k在条件C5中定义。

定理1表明本文提出的有效工具变量识别方法允许工具变量的维数为 $q=o\left(\exp \left(n^{1-2k}\right)\right)$，即在一些正则性条件下，提出的有效工具变量识别方法可以处理超高维的工具变量。另外，如下定理2表明本文提出的有效工具变量识别方法可以渐近以概率1识别出所有的有效工具变量。

定理2. 令 $v_n=c_8{n}^{-2k}$，其中 $c_8$ 是一个正常数并且满足 $0<c_8\le b_1/2$， $b_1=c_5/\left(c_1{c}_3\right)$，那么在正则性条件C1~C5下，可得

$P\left\{{\mathcal{M}}_{\ast }\subset {\hat{\mathcal{M}}}_{v_n}\right\}\ge 1-m\exp \left(-c_7{n}^{1-2k}\right)$，

其中 $m=\left|{\mathcal{M}}_{\ast }\right|$ 表示 ${\mathcal{M}}_{\ast }$ 中的指标个数。

显然，如果 ${\max }_{s\notin {\mathcal{M}}_{\ast }}\left|{\omega }_s^M\right|=o\left(n^{-k}\right)$ 成立，那么从定理1有，对于任给定的 $b_2>0$， ${\max }_{s\notin {\mathcal{M}}_{\ast }}\left|{\hat{\omega }}_s^M\right|\le b_2{n}^{-k}$ 以概率趋于1成立。因此通过定理2中 $v_n$ 的选择，可以保证有效工具变量的识别具有一致相合性，即 $P\left\{{\mathcal{M}}_{\ast }={\hat{\mathcal{M}}}_{v_n}\right\}=1-o\left(1\right)$。

接下来研究模型参数分量估计 $\hat{\beta }$ 和非参数函数系数估计 $\hat{\theta }\left(u\right)$ 的渐进性质。记 ${\Gamma }^{*}$ 表示所识别出的有效工具变量 ${\xi }^{*}$ 对应的系数矩阵。再结合定理2可知，内生性协变量 $X_i$ 以概率趋于1可写为 $X_i={\Gamma }^{\ast }{Z}_i^{\ast }+e_i$， $i=1,\cdots ,n$。为了给出 $\hat{\beta }$ 和 $\hat{\theta }\left(u\right)$ 的渐进性质，首先给出一些关于模型(1)和(2)的正则性条件如下：

C6：非参数函数系数 $\theta \left(u\right)$ 在区间 $\left[0,1\right]$ 内r阶连续可微，其中 $r\ge 2$。

C7：记 $c_1,\cdots ,c_k$ 是 $\left[0,1\right]$ 的内部结点。那么存在常数c，使得 $\max \left\{\left|{\Delta }_{i+1}-{\Delta }_i\right|\right\}=o\left({\kappa }^{-1}\right)$ 以及 $\max \left\{{\Delta }_i\right\}/\min \left\{{\Delta }_i\right\}\le c$，其中 ${\Delta }_i=c_i-c_{i-1}$， $c_0=0$， $c_{\kappa +1}=1$ 和 $\kappa$ 表示内部结点的个数。

C8：有效工具变量对应的系数矩阵 ${\Gamma }^{\ast }$ 是一个行满秩的矩阵。另外记 ${\tau }_1$ 和 ${\tau }_2$ 分别成为矩阵 $E\left({\xi }^{\ast }{\xi }^{\ast \text{T}}\right)$ 的最小和最大的特征值，并且 ${\rho }_1$ 和 ${\rho }_2$ 分别是矩阵 ${\Gamma }^{\ast }{\Gamma }^{\ast \text{T}}$ 的最小和最大的特征值，那么存在常数 $0<{\rho }_1^{\ast }<{\rho }_2^{\ast }<\infty$ 和 $0<{\tau }_1^{\ast }<{\tau }_2^{\ast }<\infty$，使得 ${\rho }_1^{\ast }<{\rho }_1<{\rho }_2<{\rho }_2^{\ast }$ 和 ${\tau }_1^{\ast }<{\tau }_1<{\tau }_2<{\tau }_2^{\ast }$ 成立。

在正则性条件C1~C8下，如下定理3表明估计量 $\hat{\beta }$ 是一致相合的。

定理3. 假设正则性条件C1~C8成立，并且内部结点的数量 $\kappa$ 满足 $\kappa =O\left(n^{1/\left(2r+1\right)}\right)$。那么有

$\Vert \hat{\beta }-{\beta }_0\Vert =O_p\left(\sqrt{\kappa /n}\right)$，

其中 ${\beta }_0$ 是 $\beta$ 的参数真值，r由条件C6定义。

另外如下定理4表明非参数函数系数的估计量 $\hat{\theta }\left(u\right)$ 也是一致相合的，并达到了最优的非参数收敛速度。

定理4. 假设正则性条件C1~C8成立，并且内部结点的数量 $\kappa$ 满足 $\kappa =O\left(n^{1/\left(2r+1\right)}\right)$。那么有

$\Vert {\hat{\theta }}_s\left(u\right)-{\theta }_s\left(u\right)\Vert =O_p\left(n^{-r/\left(2r+1\right)}\right)$，

其中r由条件C6定义。

4. 数值模拟研究

在本节中，我们通过一些蒙特卡洛数值模拟实验来研究本文估计方法的限样本性质。为实施模拟，数据从如下模型产生：

$\{\begin{array}{l}{Y}_i=X_i\beta +Z_i\theta \left(U_i\right)+{\epsilon }_i\\ X_i=\Gamma {\xi }_i+\alpha {\epsilon }_i\end{array}$

其中 $\beta =2$ 和 $\theta \left(u\right)=0.5\sin \left(u\right)+0.5\cos \left(u\right)$。我们取 $\Gamma =\left(2,1.5,1,0.5,0,\cdots ,0\right)$ 是一个 $1\times q$ 维矩阵，并且工具变量由 ${\xi }_{is}~N\left(1,1.5\right)$， $s=1,\cdots ,q$ 产生。显然，从工具变量的生成机制可以看出 ${\xi }_{is}$， $s=1,\cdots ,4$ 是有效工具变量，且其他工具变量是无效工具变量。内生性变量 $X_i$ 根据模型产生，其中模型误差取为 $\epsilon ~N\left(0,0.5\right)$，并且 $\alpha$ 可分别取为 $\alpha =0.2$ 和 $\alpha =0.8$ 两种情况用来表示协变量内生性的不同水平。另外，外生性变量 $U_i$ 来自 $U_i~U\left(0,1\right)$， $Z_i$ 来自 $Z_i~N\left(0,1\right)$，响应变量 $Y_i$ 根据模型生成。

在我们的迭次算法中有许多惩罚性函数可供使用，例如SCAD惩罚(参见Tibshirani [11])、Lasso惩罚(参见Fan和Li [12])和MCP惩罚(参见Zhang [13])等。注意到Zhao和Xue [14] 指出变量选择的模拟结果对惩罚函数的选择不是太敏感，因此在接下来的模拟中，我们仅给出基于SCAD惩罚的数值模拟结果。另外样本大小选择分别为 $n=200$ 、400以及600三种情况，工具变量的维数也分别为 $q=500$ 、1000和1500三种情况，并且对于每一种情况，数值模拟实验进行重复1000次模拟运行。

接下来首先模拟有效工具变量识别的有限样本性质。为此我们选取如下指标进行模拟：包含所有有效工具变量的最小模型大小(MMS)、把真实的有效工具变量正确识别为有效工具变量的个数(TP)，把无效工具变量错误识别为有效工具变量的个数(FP)以及错误选择率(FSR)，其中FSR定义为 $\text{FSR}=\text{FP}/\left(\text{TP}+\text{FP}\right)$。显然，FSR代表了错误识别成有效工具变量的无效工具变量个数在所有识别成有效工具变量个数中所占的比例。在我们的模拟中，所有这些评价指标都是基于1000次重复模拟实验的平均值。基于1000次模拟运行，模拟结果如表1所示。从表1中，我们可以得出以下观察结果：

1) 对模拟过程中所设定的每一种情况，所提出的有效工具变量识别方法均可以给出相对较小的MMS。这表明所提出的有效工具变量识别方法可以有效地控制有效工具变量的边际信号，对于工具变量维数的变化是比较稳定的。

2) 对于任意给定的样本量n和工具变量维数q，在内生性协变量的两种内生性水平下得到的TP、FP和FSR等指标的表现是非常相似的。这表明当内生性水平较低时，所提出的有效工具变量识别方法仍可以识别出有效的工具变量。

3) 对于任意给定的协变量内生性水平，模拟的FP和FSR指标均随着样本量n的增加而减小，而TP指标随着样本量的增加趋于真实有效工具变量的个数。这表明所提出的有效工具变量识别方法可以渐近相合地识别出有效工具变量。

Table 1. The simulation results of the optimal instrumental variable identification

表1. 有效工具变量识别的模拟结果

接下来，我们模拟参数分量 $\beta$ 和非参数函数系数 $\theta \left(u\right)$ 的估计效果。为了进行比较，我们把本文所提出的基于工具变量调整的估计(IVE)方法与忽略协变量的内生性的朴素估计(NE)方法进行比较。这里的朴素估计方法是指忽略协变量 $X_i$ 的内生性，并通过直接最小化如下目标函数来得到 $\beta$ 和 $\theta \left(u\right)$ 的估计量：

${\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(Y_i-X_i^{\text{T}}\beta -W_i^{\text{T}}\gamma \right)}^2}$，

其中 $W_i=B\left(U_i\right)\cdot Z_i$， $\theta \left(u\right)$ 由 $\theta \left(u\right)\approx B{\left(u\right)}^{\text{T}}\gamma$ 近似逼近。

对于参数分量 $\beta$，本文仅给出当工具变量维数为 $q=1000$ 时的模拟结果。在工具变量维数为其他情况下的模拟结果是类相似的，因此不再一一展示。基于1000次重复模拟实验，关于参数分量 $\beta$ 估计量的1000个绝对偏差值 $\left|\hat{\beta }-\beta \right|$ 的箱线图见图1和图2，其中图1是在内生性水平为 $\alpha =0.02$ 下的模拟结果，图2是在内生性水平为 $\alpha =0.08$ 下的模拟结果。由图1和图2可以得到如下结论：

1) 随着样本量的增加，基于IVE估计方法给出的绝对偏差值则明显减小，而基于NE方法给出的绝对偏差值即使样本量增加的情况下仍然较大。这就表明所提出的IVE方法优于NE方法，而且忽略变量内生性的NE方法将会产生一定的内生性偏差。

2) 对于任意给定的样本量n，基于IVE估计方法所给出的模拟结果在不同的内生性水平下是非常相似的。这也表明本文提出的工具变量调整机制可以减弱协变量内生性对估计精度带来的影响。

Figure 1. Under the endogenous level $\alpha =0.2$, the box-plots of absolute bias values for the estimator $\hat{\beta }$

图1. 在内生性水平 $\alpha =0.2$ 下，估计量 $\hat{\beta }$ 的绝对偏差箱线图

Figure 2. Under the endogenous level $\alpha =0.8$, the box-plots of absolute bias values for the estimator $\hat{\beta }$

图2. 在内生性水平 $\alpha =0.8$ 下，估计量 $\hat{\beta }$ 的绝对偏差箱线图

另外，图3给出了在样本量 $n=400$ 的情况下，非参数函数系数 $\theta \left(u\right)$ 的模拟结果，其中短虚线表示基于IVE方法得到的估计曲线，点虚线表示基于NE方法得到的估计曲线，实线表示 $\theta \left(u\right)$ 的真实曲线。从图3可以看出，NE方法给出的估计曲线有相对较大的估计偏差，而本文提出的IVE方法则可以减弱估计偏差，给出的估计曲线非常接近真实曲线。另外还可以看出，在不同的内生性水平下，基于IVE方法给出的估计曲线是非常相似的。综上，上述模拟结果表明，协变量的内生性将会影响参数分量和非参数分量估计的一致相合性，而本文所提出的IVE方法可有效地削弱协变量的内生性对模型参数分量及非参数分量估计的影响。

Figure 3. The estimator of $\theta \left(u\right)$ under the the endogenous level $\alpha =0.2$ (a) and $\alpha =0.8$ (b), respectively

图3. $\theta \left(u\right)$ 的估计曲线，其中(a)是在内生性水平为 $\alpha =0.2$ 时的模拟结果，(b)是在 $\alpha =0.8$ 时的模拟结果

5. 定理的证明

在本节中，我们将给出定理1~4的详细证明过程。令 $E_n\left\{m\left(Y,{\xi }_s\right)\right\}=n^{-1}{\displaystyle {\sum }_{i=1}^{n}m\left(Y_i,{\xi }_{is}\right)}$ 表示数学期望对应的经验测度。另外，定义 ${\mathcal{B}}_s\left(N\right)=\left\{{\omega }_s\in \mathcal{B}:\left|{\omega }_s-{\omega }_s^M\right|\le N\right\}$ 其中 $N>0$ 为某一正常数，并进一步定义 $F_s\left(N\right)={\sup }_{{\omega }_s\in {\mathcal{B}}_s\left(N\right)}\left|\left(E_n-E\right)\left\{Q\left(Y_i,{\omega }_s{\xi }_{is}\right)-Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\right|$。首先给出一些定理证明过程中需要的引理。

引理1. 假设正则条件C1~C5成立，那么有

$\underset{s\in {\mathcal{M}}_{\ast }}{\min }\left|{\omega }_s^M\right|\ge b_1{n}^{-k}$，

其中 $k\in \left(0,1/2\right)$， $b_1=c_5/\left(c_1{c}_3\right)$ 和 ${\omega }_s^M$ 由(9)式定义。

证明：根据条件C2可知

$\left|Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)-Q\left(Y_i,0\right)\right|\le c_1\left|{\omega }_s^M\right|\left|{\xi }_{is}\right|$ (10)

对(10)式两边同时取期望可得

$\left|E\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)-Q\left(Y_i,0\right)\right\}\right|\le E\left\{\left|Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)-Q\left(Y_i,0\right)\right|\right\}\le c_1\left|{\omega }_s^M\right|E\left\{\left|{\xi }_{is}\right|\right\}\le c_1{c}_3\left|{\omega }_s^M\right|$ (11)

另外由条件C5可知

$\left|E\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)-Q\left(Y_i,0\right)\right\}\right|\ge c_5{n}^{-k}$ (12)

从(11)和(12)式可知，对任意 $s\in {\mathcal{M}}_{\ast }$，则有 $\left|{\omega }_s^M\right|\ge b_1{n}^{-k}$。即引理1的证明完毕。

引理2. 假设正则条件C1~C5成立，那么对任给定的 $\zeta >0$ 和 $s\in \left\{1,\cdots ,q\right\}$，则有：

$P\left\{F_s\left(N\right)\ge 4N{b}_3{n}^{-1/2}\left(1+\zeta \right)\right\}\le \exp \left(-2{\zeta }^2\right)$

其中 $b_3=c_1{c}_3$。

证明：结合条件C2，并利用对称化定理(参见Vaart和Wellner [15])可得

$\begin{array}{c}E\left\{F_s\left(N\right)\right\}\le 2E\left\{\underset{{\omega }_s\in {\mathcal{B}}_s\left(N\right)}{\sup }\left|E_n{\delta }_i\left[Q\left(Y_i,{\omega }_s{\xi }_{is}\right)-Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right]\right|\right\}\\ \le 4c_{1}E\left\{\underset{{\omega }_s\in {\mathcal{B}}_s\left(N\right)}{\sup }\left|E_n\left[{\delta }_i{\xi }_{is}\left({\omega }_s-{\omega }_s^M\right)\right]\right|\right\}\end{array}$ (13)

其中 ${\delta }_1,\cdots ,{\delta }_n$ 表示分别以概率为1/2取值为±1的独立且同分布随机变量序列。另外再通过利用Jensen不等式可以得到

$\begin{array}{c}E\left\{\underset{{\omega }_s\in {\mathcal{B}}_s\left(N\right)}{\sup }\left|E_n\left[\delta {\xi }_{is}\left({\omega }_s-{\omega }_s^M\right)\right]\right|\right\}\le \underset{{\omega }_s\in {\mathcal{B}}_s\left(N\right)}{\sup }\left|{\omega }_s-{\omega }_s^M\right|\times E\left[\left|E_n\delta {\xi }_{is}\right|\right]\\ \le N{\left(E{\left|E_n\delta {\xi }_{is}\right|}^2\right)}^{1/2}\\ =\frac{N}{\sqrt{n}}{\left[E\left\{{\delta }^2{\xi }_{is}^2\right\}\right]}^{1/2}\\ \le N{c}_3{n}^{-1/2}\end{array}$ (14)

结合(13)和(14)式则可以得到 $E\left\{F_s\left(N\right)\right\}\le 4N{b}_1{n}^{-1/2}$。另外再由条件C2可知在集合 ${\mathcal{B}}_s\left(N\right)$ 上有

$\left|Q\left(Y_i,{\omega }_s{\xi }_{is}\right)-Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right|\le b_{1}N$

因此，令 $L=2N{b}_1$ 并利用集中定理(参见Massart [16])可得

$\begin{array}{c}\mathrm{Pr}\left\{F_s\left(N\right)\ge 4N{b}_1{n}^{-1/2}\left(1+\zeta \right)\right\}\le \mathrm{Pr}\left[F_s\left(N\right)\ge E\left\{F_s\left(N\right)\right\}+4N{b}_1{n}^{-1/2}\zeta \right]\\ \le \exp \left[-n{\left(4N{b}_1{n}^{-1/2}\zeta \right)}^2/\left(2L^2\right)\right]\\ =\exp \left(-2{\zeta }^2\right)\end{array}$

即引理2证明完毕。

定理1证明. 令 ${\omega }_{\lambda s}=\lambda {\hat{\omega }}_s^M+\left(1-\lambda \right){\omega }_s^M$，其中 $\lambda ={\left(1+\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|/N\right)}^{-1}\le 1$，且N满足条件(C1)和(C2)。则一个简单的计算有 $\left|{\omega }_{\lambda s}-{\omega }_s^M\right|=\lambda \left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\le N$，其表明 ${\omega }_{\lambda s}\in {\mathcal{B}}_s\left(N\right)$。因为 ${\hat{\omega }}_s^M$ 是 $E_n\left\{Q\left(Y_i,{\omega }_s{\xi }_{is}\right)\right\}$ 的最小取值，然后调用凸性，我们得到

$\begin{array}{c}{E}_n\left\{Q\left(Y_i,{\omega }_{\lambda s}{\xi }_{is}\right)\right\}\le \lambda E_n\left\{Q\left(Y_i,{\hat{\omega }}_s^M{\xi }_{is}\right)\right\}+\left(1-\lambda \right)E_n\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\\ \le \lambda E_n\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}+\left(1-\lambda \right)E_n\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\\ =E_n\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\end{array}$ (15)

另外，注意到 ${\omega }_s^M$ 是 $E\left\{S\left(Y_i,{\omega }_s{\xi }_{is}\right)\right\}=0$ 的解，因此有 ${\omega }_s^M$ 是 $E\left\{Q\left(Y_i,{\omega }_s{\xi }_{is}\right)\right\}$ 的最小值点。因此有

$E\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\le E\left\{Q\left(Y_i,{\omega }_{\lambda s}{\xi }_{is}\right)\right\}$ (16)

进而结合(15)和(16)式，通过简单计算可以得到

$E\left\{Q\left(Y_i,{\omega }_{\lambda s}{\xi }_{is}\right)\right\}-E\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\le \left(E-E_n\right)\left\{Q\left(Y_i,{\omega }_{\lambda s}{\xi }_{is}\right)-Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\le F_s\left(N\right)$ (17)

其中 $F_s\left(N\right)={\sup }_{{\omega }_s\in {\mathcal{B}}_s\left(N\right)}\left|\left(E-E_n\right)\left\{Q\left(Y_i,{\omega }_{\lambda s}{\xi }_{is}\right)-Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\right|$。另外再由条件C5可得

$\left|E\left\{Q\left(Y_i,{\omega }_{\lambda s}{\xi }_{is}\right)\right\}-E\left\{Q\left(Y_i,{\omega }_s^M{\xi }_{is}\right)\right\}\right|\ge c_2{\left|{\omega }_{\lambda s}-{\omega }_s^M\right|}^2$ (18)

因此，结合(17)和(18)式可以得到

$\left|{\omega }_{\lambda s}-{\omega }_s^M\right|\le {\left(F_s\left(N\right)/c_2\right)}^{1/2}$ (19)

对任给定的 $\delta$ 并结合(19)式可得 $P\left(\left|{\omega }_{\lambda s}-{\omega }_s^M\right|\ge \delta \right)\le P\left(F_s\left(N\right)\ge c_2{\delta }^2\right)$。取 $\delta =N/2$，那么有

$P\left(\left|{\omega }_{\lambda s}-{\omega }_s^M\right|\ge N/2\right)\le P\left(F_k\left(N\right)\ge c_2{N}^2/4\right)$ (20)

结合 ${\omega }_{\lambda s}$ 的定义，简单计算可得 ${\omega }_{\lambda s}-{\omega }_s^M=\lambda \left({\hat{\omega }}_s^M-{\omega }_s^M\right)$。因此可以得到

$\begin{array}{c}P\left(\left|{\omega }_{\lambda s}-{\omega }_s^M\right|\ge N/2\right)=P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge {\lambda }^{-1}N/2\right)\\ =P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge \left[1+\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|/N\right]N/2\right)\\ =P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge N\right)\end{array}$ (21)

取 $N=16\left(c_1{c}_3/c_2\right)n^{-1/2}\left(1+\zeta \right)$，则由(20)和(21)可得

$P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge 16\left(c_1{c}_3/c_2\right)n^{-1/2}\left(1+\zeta \right)\right)\le P\left(F_s\left(N\right)\ge 4N{c}_1{c}_3{n}^{-1/2}\left(1+\zeta \right)\right)$ (22)

因此结合(22)式并利用引理2可得

$P\left\{\sqrt{n}\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge 16\left(c_1{c}_3/c_2\right)\left(1+\zeta \right)\right\}\le \exp \left(-2{\zeta }^2\right)$

取 $1+\zeta =\left(c_2{c}_6\right)/\left(c_1{c}_3\right)n^{1/2-k}/16$，那么则有

$P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge c_6{n}^{-k}\right)\le \exp \left(-2{\left[\left(c_2{c}_6\right)/\left(c_1{c}_3\right)n^{1/2-k}/16-1\right]}^2\right)$ (23)

注意到 ${\left[\left(c_2{c}_6\right)/\left(c_1{c}_3\right)n^{1/2-k}/16-1\right]}^2=O\left(n^{1-2k}\right)$，因此由(23)式可得

$P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge c_6{n}^{-k}\right)\le \exp \left(-c_7{n}^{1-2k}\right)$ (24)

其中 $c_7$ 为一个正常数。另外简单计算可得

$P\left(\underset{1\le s\le q}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge c_6{n}^{-k}\right)\le {\displaystyle \underset{s=1}{\overset{q}{\sum }}P\left(\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge c_6{n}^{-k}\right)}$ (25)

进而结合(24)和(25)式可得

$P\left(\underset{1\le s\le q}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge c_6{n}^{-k}\right)\le q\exp \left(-c_7{n}^{1-2k}\right)$

这就完成了定理1的证明。

定理2证明. 我们首先定义 ${\mathcal{A}}_n$ 如下

${\mathcal{A}}_n=\left\{\underset{s\in {\mathcal{M}}_{\ast }}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\le b_1{n}^{-k}/2\right\}$

其中 $b_1$ 在引理1中定义。由引理1可知对任意 $s\in {\mathcal{M}}_{\ast }$ 有 ${\left|{\omega }_s^M\right|}^2\ge b_1{n}^{-k}$，因此在 ${\mathcal{A}}_n$ 上利用三角不等式可以得到

$\left|{\hat{\omega }}_s^M\right|\ge \left|{\omega }_s^M\right|-\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge \left|{\omega }_s^M\right|-\underset{s\in {\mathcal{M}}_{\ast }}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge b_1{n}^{-k}/2$

取 $v_n=c_8{n}^{-k}$，其中 $0<c_8<b_1/2$，那么可得

$P\left({\mathcal{M}}_{\ast }\subset {\hat{\mathcal{M}}}_{v_n}\right)\ge 1-P\left({\mathcal{A}}_n^c\right)=1-P\left(\underset{s\in {\mathcal{M}}_{\ast }}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge b_1{n}^{-k}/2\right)$ (26)

另外由定理1可得

$P\left(\underset{s\in {\mathcal{M}}_{\ast }}{\max }\left|{\hat{\omega }}_s^M-{\omega }_s^M\right|\ge b_1{n}^{-k}/2\right)\le m\exp \left(-c_7{n}^{1-2k}\right)$ (27)

其中m是 ${\mathcal{M}}_{\ast }$ 中的元素数。那么结合(26)和(27)式则完成了定理2的证明。

定理3证明. 由第3节中的讨论可知 $X_i={\Gamma }^{\ast }{\xi }_i^{\ast }+e_i$ 以概率趋于1成立。另外注意 ${\hat{\Gamma }}^{\ast }$ 是 ${\Gamma }^{\ast }$ 的矩估计量，因此可得 ${\hat{\Gamma }}^{\ast }={\Gamma }^{\ast }+O_p\left(\sqrt{1/n}\right)$。再结合 $E\left(e_i\right)=0$，并通过简单计算可得

$\begin{array}{c}\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(X_i-X_i^{*}\right)}^{\text{T}}\beta }=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left({\Gamma }^{\ast }{\xi }_i^{\ast }+e_i-{\hat{\Gamma }}^{\ast }{\xi }_i^{\ast }\right)}^{\text{T}}\beta }\\ =\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\xi }_i^{\ast \text{T}}{\left({\Gamma }^{\ast }-{\hat{\Gamma }}^{\ast }\right)}^{\text{T}}\beta }+\frac{1}{\sqrt{n}}\left(\frac{1}{\sqrt{n}}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{e}_i^{\text{T}}\beta }\right)\\ =O_p\left(\Vert {\Gamma }^{\ast }-{\hat{\Gamma }}^{\ast }\Vert \right)+O_p\left(n^{-1/2}\right)\\ =O_p\left(\sqrt{1/n}\right)\end{array}$ (28)

记 $\beta ={\beta }_0+{\delta }_n{\eta }_1$， $\gamma ={\gamma }_0+{\delta }_n{\eta }_2$ 和 $R\left(U_i\right)=Z_i^{\text{T}}\theta \left(U_i\right)-W_i^{\text{T}}{\gamma }_0$，其中 ${\delta }_n=\sqrt{\kappa /n}$， ${\beta }_0$ 和 ${\gamma }_0$ 是参数真值，且 ${\eta }_1$ 和 ${\eta }_2$ 表示常数向量。那么由(28)式可得到

$\begin{array}{c}{M}_n\left(\beta ,\gamma \right)=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(Y_i-X_i^{\ast \text{T}}\beta -W_i^{\text{T}}\gamma \right)}^2}\\ =\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left[e_i^{\text{T}}{\beta }_0+{\epsilon }_i+R\left(U_i\right)-{\delta }_n{\xi }_i^{\ast \text{T}}{\Gamma }^{\ast \text{T}}{\eta }_1-{\delta }_n{W}_i^{\text{T}}{\eta }_2\right]}^2+O_p\left({\Vert {\Gamma }^{\ast }-{\hat{\Gamma }}^{\ast }\Vert }^2\right)}\\ =\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left[e_i^{\text{T}}{\beta }_0+{\epsilon }_i+R\left(U_i\right)-{\delta }_n{\psi }_i^{\text{T}}\eta \right]}^2+O_p\left(1/n\right)}\end{array}$ (29)

其中 ${\psi }_i={\left({\xi }_i^{\ast \text{T}}{\Gamma }^{\ast \text{T}},W_i^{\text{T}}\right)}^{\text{T}}$ 和 $\eta ={\left({\eta }_1^{\text{T}},{\eta }_2^{\text{T}}\right)}^{\text{T}}$。类似地，我们可得到

$M_n\left({\beta }_0,{\gamma }_0\right)=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(Y_i-X_i^{\ast \text{T}}{\beta }_0-W_i^{\text{T}}{\gamma }_0\right)}^2}=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left[e_i^{\text{T}}{\beta }_0+{\epsilon }_i+R\left(U_i\right)\right]}^2+O_p\left(1/n\right)}$ (30)

进一步记 ${\Delta }_n\left(\beta ,\gamma \right)=M_n\left(\beta ,\gamma \right)-M_n\left({\beta }_0,{\gamma }_0\right)$，那么由(29)和(30)可得

$\begin{array}{c}{\Delta }_n\left(\beta ,\gamma \right)=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\delta }_n^2{\psi }_i^{\text{T}}\eta {\eta }^{\text{T}}{\psi }_i}-\frac{2}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\delta }_n\left(e_i^{\text{T}}{\beta }_0+{\epsilon }_i\right){\psi }_i^{\text{T}}\eta }\\ -\frac{2}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\delta }_{n}R\left(U_i\right){\psi }_i^{\text{T}}\eta }+o_p\left({\delta }_n^2\right)\\ \equiv J_1+J_2+J_3+O_p\left({\delta }_n^2\right)\end{array}$ (31)

结合 $E\left(e_i^{\text{T}}{\beta }_0+{\epsilon }_i|{\xi }_i^{\ast },W_i\right)=0$，简单计算可得

$J_2=-\frac{2{\delta }_n}{\sqrt{n}}\left[\frac{1}{\sqrt{n}}{\displaystyle \underset{i=1}{\overset{n}{\sum }}\left(e_i^{\text{T}}{\beta }_0+{\epsilon }_i\right){\psi }_i^{\text{T}}\eta }\right]=O_p\left(\frac{{\delta }_n}{\sqrt{n}}\right)\Vert \eta \Vert =o_p\left({\delta }_n^2\right)\Vert \eta \Vert$ (32)

另外由Schumaker [17] 可知 $\Vert R\left(U_i\right)\Vert =O_p\left(k^{-r}\right)=O_p\left({\delta }_n\right)$。因此简单计算可得

$J_3=O_p\left({\delta }_n\Vert R\left(U_i\right)\Vert \right)\Vert \eta \Vert =O_p\left({\delta }_n^2\right)\Vert \eta \Vert$ (33)

注意到 $J_1=O_p\left({\delta }_n^2\right){\Vert \eta \Vert }^2$，那么对于一个充分大的常数C，在 $\Vert \eta \Vert =C$ 条件下 $J_2$ 和 $J_3$ 完全可由 $J_1$ 控制。另外注意 $J_1$ 是一个正数，因此结合(31)~(33)式可得，对任给的 $\upsilon >0$，存在一个大常数C使得

$P\left\{\underset{\Vert \eta \Vert =C}{\inf }{M}_n\left(\beta ,\gamma \right)>M_n\left({\beta }_0,{\gamma }_0\right)\right\}\ge 1-\upsilon$ (34)

另外由 $M_n\left(.,.\right)$ 的凸性可知，(34)式表明存在一个局部最小值点 $\hat{\beta }$ 和 $\hat{\gamma }$，使得 $\hat{\beta }-{\beta }_0=O_p\left({\delta }_n\right)$ 和 $\hat{\gamma }-{\gamma }_0=O_p\left({\delta }_n\right)$ 以不小于 $1-\upsilon$ 的概率成立。再由 $\upsilon$ 的任意性，则完成定理3的证明。

定理4证明. 利用类似Zhao和Xue [14] 的讨论可得

${\Vert {\hat{\theta }}_s\left(u\right)-{\theta }_s\left(u\right)\Vert }^2\le 2{\left({\hat{\gamma }}_s-{\gamma }_{s0}\right)}^{\text{T}}H\left({\hat{\gamma }}_s-{\gamma }_{s0}\right)+{\displaystyle {\int }_0^1{R}_s{\left(u\right)}^2\text{d}u}$ (35)

其中 $R_s\left(u\right)={\theta }_s\left(u\right)-B^{\text{T}}\left(u\right){\gamma }_{s0}$ 和 $H={\displaystyle {\int }_0^{1}B\left(u\right)B^{\text{T}}\left(u\right)\text{d}u}$。另外，由定理3的证明过程可知

$\Vert {\hat{\gamma }}_s-{\gamma }_{s0}\Vert =O_p\left(\sqrt{\kappa /n}\right)$。结合条件C6和 $\kappa =O\left(n^{1/\left(2r+1\right)}\right)$，则可以得到 $O_p\left(\sqrt{\kappa /n}\right)=O_p\left(n^{-r/\left(2r+1\right)}\right)$。因此再结合 $\Vert H\Vert =O\left(1\right)$ 可得

${\left({\hat{\gamma }}_s-{\gamma }_{s0}\right)}^{\text{T}}H\left({\hat{\gamma }}_s-{\gamma }_{s0}\right)=O_p\left(n^{\frac{-2r}{2r+1}}\right)$ (36)

另外注意到 $\Vert R_s\left(u\right)\Vert =O\left({\kappa }^{-r}\right)=O\left(n^{\frac{-r}{2r+1}}\right)$，则可以得到

${\displaystyle {\int }_0^1{R}_s{\left(u\right)}^2\text{d}u}=O_p\left(n^{\frac{-2r}{2r+1}}\right)$ (37)

结合(35)~(37)式则完成了定理4的证明。

6. 结论

本文在工具变量为超高维的情况下，对一类含内生性协变量的半参数模型的估计问题进行了研究。结合构造辅助回归模型和惩罚估计技术，提出了一种基于双惩罚的正则估计方法。所提出的估计方法不仅可以给出模型参数分量和非参数分量的相合估计，还可以识别出有效工具变量。模拟研究表明，即使工具变量为超高维数据的情况下，本文提出的方法仍可以相合的识别出所有有效工具变量。另外，模拟结果还表明本文提出的基于工具变量调整的估计方法对模型参数和非参数函数系数的估计也是行之有效的。

注意到本文考虑的半参数模型具有较好的灵活性和广泛的适应性，经典的线性模型、部分线性模型以及部分线性可加模型均是本文所研究模型的特殊情形。因此，本文提出的双惩罚正则估计方法可以方便地推广到超高维数据下部分线性模型以及部分线性可加模型等半参数模型的统计推断研究上。另外，在本文的研究中，工具变量 ${\xi }_i$ 假定为超高维数据，而内生性协变量 $X_i$ 则假定为维数固定的低维数据。因此，接下来一个值得进一步研究的课题是在工具变量 ${\xi }_i$ 和内生性协变量 $X_i$ 均为超高维数据的情况下，研究超高维半参数工具变量模型的统计推断问题。

基金项目

本文研究成果得到国家社会科学基金项目《高维内生协变量的半参数建模及其在环境治理绩效测度中的应用》(编号：18BTJ035)资助。

NOTES

^*通讯作者。