高维广义线性模型的惩罚拟似然SCAD估计

高维广义线性模型的惩罚拟似然SCAD估计 Penalized Quasi-Likelihood SCAD Estimator in High-Dimensional Generalized Linear Models first-author 闫莉，女，讲师，主要从事复杂数据统计分析方面的研究. E-mail: lyan@snnu.edu.cn 闫莉，女，讲师，主要从事复杂数据统计分析方面的研究. E-mail: lyan@snnu.edu.cn lyan@snnu.edu.cn corresp 陈夏, E-mail: xchen80@snnu.edu.cn 陈夏, E-mail: xchen80@snnu.edu.cn xchen80@snnu.edu.cn 陕西师范大学 数学与信息科学学院，陕西 西安 710119 陕西师范大学 数学与信息科学学院，陕西 西安 710119 School of Mathematics and Information Science, Shaanxi Normal University, Xi'an 710119, Shaanxi, China School of Mathematics and Information Science, Shaanxi Normal University, Xi'an 710119, Shaanxi, China 基于光滑切片绝对偏差(smoothly clipped absolute deviation，SCAD)惩罚的拟似然方法，研究高维广义线性模型的变量选择和参数估计问题.所提方法只需要对响应变量期望函数和方差函数的正确设定.在适当的正则条件下，证明了拟似然SCAD估计具有相合性和Oracle性质.最后通过数据模拟和实例分析，验证了所提方法的有限样本性质. This paper considers the problem of variable selection and estimation in high-dimensional generalized linear models based on the quasi-likelihood method with SCAD(smoothly clipped absolute deviation) penalty. This method requires only the correct specification of the expectation and variance function of the response variable. Under reasonable conditions, the consistency and the Oracle property of the quasi-likelihood SCAD estimator is proved. Some simulations and a real data analysis are given to illustrate the and good performance of the proposed method. 广义线性模型 高维数据 变量选择 拟似然方法 光滑切片绝对偏差 generalized linear models high-dimensional data variable selection quasi-likelihood method SCAD (smoothly clipped absolute deviation) 引言 1 广义线性模型(generalized linear models, GLMs)的理论是对经典线性模型理论的重要推广，自Nelder与Wedderburn [ 1 1 ] 正式引入GLMs以来，已经成功应用到了很多领域.起初，研究者们只局限于响应变量服从指数分布族的情形.1974年，Wedderburn [ 2 2 ] 提出了GLMs拟似然函数的概念，该方法在建模时可以只要求响应变量期望函数与方差函数有正确设定，不要求响应变量的分布为指数型分布族.后续研究表明，只对期望函数有正确设定但方差函数未知的情况下，此方法依然有效 [ 3 3 , 4 4 ] .这种研究方法称为拟似然方法，得到的参数估计称为拟似然估计. 在经济和生物医学等领域，常常会遇到多变量大样本的数据，因此变量选择在高维数据分析中发挥着重要作用.在高维数据变量选择问题的研究中，应用较为广泛的是基于惩罚函数的估计方法.它可以同时实现变量选择和参数估计.常见的惩罚函数有SCAD(smoothly clipped absolute deviation) [ 5 5 ] 、Bridge [ 6 6 ] 、Lasso [ 7 7 ] 以及Lasso的衍生版本，如Adaptive Lasso [ 8 8 ] 、Elastic Net [ 9 9 ] 、Adaptive Elastic Net [ 10 10 ] 及针对群组数据的Group Lasso [ 11 11 ] 等.对于经典线性模型，基于惩罚函数的变量选择方法已经取得了众多研究成果.例如，Huang和Xie [ 12 12 ] 研究了高维线性模型下SCAD估计的渐近性质.Knight和Fu [ 13 13 ] 研究了参数个数有限时，Bridge估计的渐近性质.Huang等 [ 14 14 ] 又将此结果拓展到参数无穷维的情形.在此基础上，Wang等 [ 15 15 ] 研究了广义线性模型中参数发散时，Bridge估计的渐近性质.Xie和Huang [ 16 16 ] 研究了SCAD惩罚下，高维部分线性模型参数估计及其渐近性质.Fan和Li [ 5 5 ] 以及Fan和Peng [ 17 17 ] 证明了参数发散时，SCAD估计仍具有Oracle性质.Zou [ 8 8 ] 证明了Lasso方法不具有Oracle性质，并提出了Adaptive Lasso方法.更多关于高维变量选择的研究参见文献[ 18 18 ]. 在文献中尚未见到利用惩罚的拟似然SCAD方法，研究高维广义线性模型的参数估计和变量选择问题.拟似然方法不需要知道响应变量的分布，只需均值函数和方差函数的正确设定，这大大扩大了其应用范围.基于此，本文考虑在高维广义线性模型中，利用惩罚拟似然SCAD方法进行变量选择和参数估计.在适当的条件下，证明了拟似然SCAD估计的相合性和Oracle性质.最后通过数据模拟和实例分析，说明了所提方法的有限样本性质. GLMs的惩罚拟似然SCAD估计 1 考虑如下广义线性模型 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257361&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257361&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257361&type=middle 其中 y i 是响应变量， X i =( X i 1 , …, X ip n ) T 是 p n ×1维协变量， β n 是 p n ×1维的未知参数向量， g (·)的逆是联系函数. ε 1 , ε 2 , …, ε n 是相互独立的随机误差序列项，满足 E ( ε i )=0，Var( ε i )= σ ( μ i )，这里 μ i = g ( X i T β n ). 1974年，Wedderburn提出如下拟似然函数 [ 2 2 ] http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257366&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257366&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257366&type=middle 其中 μ i ( β n )= g ( X i T β n ).参数 β n 的极大拟似然估计定义为方程 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257371&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257371&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257371&type=middle 的根.这里 Q n * ( β n )和指数型情况下所定义的对数似然函数的作用是相同的，而方程(1)的作用也与指数情形下的对数似然方程一致. 惩罚拟似然目标函数定义为 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257375&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257375&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257375&type=middle 其中 p λ n (·)是惩罚函数， λ n 是调整参数，SCAD惩罚函数 p λ n (·)的一阶导函数为 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257380&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257380&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257380&type=middle 其中 λ n ＞0， a ＞2.这里取Fan和Li所建议的固定调整参数 a =3.7 [ 5 5 ] .最大化惩罚拟似然目标函数 Q n ( β n )，得到参数 β n 的拟似然SCAD估计 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257522&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257522&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257522&type=middle . 记参数真值为 β n 0 ，简记为 β 0 .令 β 0 =( β 10 T , β 20 T ) T , 其中 β 10 是 k n 维非零向量， β 20 是 m n 维零向量，这里 k n 为非零系数的个数， m n 为零系数的个数.同样地，令 X =( w T , z T ) T , 相应的样本为 X i =( w i T , z i T ) T , 其中 w i =( X i 1 , …， X ik n ) T , z i =( X i ( k n +1) , …， X ip n ) T . 主要结果 1 首先给出本文中的一些记号. http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257385&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257385&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257385&type=middle 令 λ 1 n 和 λ 2 n 分别是 Σ n 的最小和最大特征值， τ 1 n 和 τ 2 n 分别是 Σ 1 n 的最小和最大特征值.下面列出本文所需要的条件： 条件1    g 定义于[ A , B ]=[ μ － ε , μ + ε ]，对于某个 ε ＞0，其中 μ =inf{ X i T β : i ≥1, β ∈ Θ n }， μ =sup{ X i T β : i ≥1, β ∈ Θ n }， g 的三阶导数存在且在[ A , B ]上有界， g 的导数[ A , B ]上处处不为0， σ 定义于[ g ( A ), g ( B )]，在其上有有界的二阶导数，且 σ 在[ g ( A ), g ( B )]上有非0下界. 条件2   { X i , i ≥1}有界. 条件3    p n →∞, p n 3 / n →0, ( nk n )1/2 λ n →0. 条件4   存在常数 b 0 , b 1 ，且0 < b 0 < b 1 ，使得 b 0 ≤min{| β 10 j |, 1≤ j ≤ k n }≤max{| β 10 j |, 1≤ j ≤ k n }≤ b 1 ，其中 β 10 j 是 β 10 的第 j 个元素. 条件5   存在常数 λ ＞0和00 < τ 1 0 < τ 2 0 < ∞使得 λ ≤ λ 1 n ， τ 1 ≤ τ 1 n ≤ τ 2 n ≤ τ 2 . 下面考虑惩罚拟似然估计的渐近性质.在本文中 C 表示一个正常数，它可以在不同的位置表示不同的数. 定理1 (相合性)  在条件1~5下，如果 nλ n 2 = O (1)，则(2)式存在一个局部最大值点 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257526&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257526&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257526&type=middle ，使得 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257390&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257390&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257390&type=middle 证   令 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257530&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257530&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257530&type=middle ，下面要证明对任给的 ε ＞0，有一个充分大的常数 C ，对足够大的 n ，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257393&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257393&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257393&type=middle 即以概率趋向于1在球 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257534&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257534&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257534&type=middle 内存在一个局部最大值点 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257537&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257537&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257537&type=middle ，使得 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257540&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257540&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257540&type=middle ，其中 u 是一个 p n ×1的向量. 易知 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257396&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257396&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257396&type=middle 首先考虑第一项 I 1 .利用泰勒展开，得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257399&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257399&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257399&type=middle 其中 β n * 是介于 β 0 + a n u 和 β 0 之间的向量.对于 I 11 ，在条件1和条件2下，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257543&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257543&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257543&type=middle .根据Cauchy-Schwarz不等式，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257402&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257402&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257402&type=middle 对于 I 12 ，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257404&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257404&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257404&type=middle 在条件1下，可以得到 E [ d 1 ( β 0 , X 1 )] 2 ≤ C ，再由切比雪夫不等式，得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257407&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257407&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257407&type=middle 因此，由Cauchy-Schwarz不等式，知 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257409&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257409&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257409&type=middle 对于 I 13 ，由条件1~3，可以得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257411&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257411&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257411&type=middle 现在考虑 I 2 ，如果 nλ n 2 = O (1)，由条件3和条件4，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257413&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257413&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257413&type=middle 因此对于一个充分大的常数 C ，(3)式的符号完全由 I 122 决定，而这一项是负的.定理证毕. 定理2 (Oracle性质)  令在定理1中的局部最大值点为 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257547&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257547&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257547&type=middle ，其中 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257549&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257549&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257549&type=middle 和 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257551&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257551&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257551&type=middle 分别是 β 1 n 和 β 2 n 的估计.假设 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257554&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257554&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257554&type=middle 且条件1~5都成立，我们有 1) 稀疏性： http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257557&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257557&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257557&type=middle . 2) 渐近正态性：令 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257561&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257561&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257561&type=middle ，其中 α n 为任意的 k n ×1维向量且满足‖ α n ‖ 2 ≤1.则 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257415&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257415&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257415&type=middle 证   首先证明定理的第一个部分.令 ε n = C http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257563&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257563&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257563&type=middle , 只需证明当 n →∞时以趋于1的概率对于任意给的 β 1 n 满足 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257566&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257566&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257566&type=middle , 有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257418&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257418&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257418&type=middle 其中 j = k n +1, …， p n .由于 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257568&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257568&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257568&type=middle , 由泰勒展开，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257422&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257422&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257422&type=middle 其中 β 0 * 是介于 β n 和 β 0 之间的向量，对于第一项 J n 1 ，由条件1, 我们得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257423&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257423&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257423&type=middle 对于 J n 2 ，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257426&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257426&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257426&type=middle 根据Cauchy-Schwarz不等式和(4)式，得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257433&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257433&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257433&type=middle 对于 J n 22 , 有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257436&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257436&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257436&type=middle 因此，得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257439&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257439&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257439&type=middle 对于 J n 3 ，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257442&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257442&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257442&type=middle 由条件1，得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257570&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257570&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257570&type=middle .利用Cauchy-Schwarz不等式，知 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257448&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257448&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257448&type=middle 且 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257453&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257453&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257453&type=middle 在条件1到条件3下，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257457&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257457&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257457&type=middle 由 J n 31 和 J n 32 ，得 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257461&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257461&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257461&type=middle 由(5), (6)和(7)式，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257465&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257465&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257465&type=middle 由 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257572&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257572&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257572&type=middle , 易得 β nj 的符号决定了 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257575&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257575&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257575&type=middle http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257579&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257579&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257579&type=middle . 现在来证明定理的第二部分.由定理1知，存在一个局部最大值点 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257582&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257582&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257582&type=middle 使得 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257468&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257468&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257468&type=middle 由泰勒展开，有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257472&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257472&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257472&type=middle 其中 β 10 * 是介于 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257583&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257583&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257583&type=middle 和 β 1 n 之间的向量.所以 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257475&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257475&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257475&type=middle 根据Cauchy-Schwarz不等式和(4)式，易知 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257479&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257479&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257479&type=middle 再利用上述证明 J n 3 的方法，可以得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257483&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257483&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257483&type=middle 因此 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257487&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257487&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257487&type=middle 对于任意的 k n ×1的向量 α n , 由于 y i = μ i ( β 10 , 0)+ ε i , (8)式可改写为 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257491&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257491&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257491&type=middle 令 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257589&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257589&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257589&type=middle , 则有 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257495&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257495&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257495&type=middle 令 v i = n －1/2 S n －1 a i , w i = v i ε i ，现在验证Lindeberg-Feller中心极限定理的条件，也就是对任给的 ε ＞0, http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257499&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257499&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257499&type=middle 只需要证明 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257597&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257597&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257597&type=middle , 或者等价于 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257599&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257599&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257599&type=middle | v i |→0.由条件1到条件3，可以得到 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257503&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257503&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257503&type=middle 因此 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257507&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257507&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257507&type=middle 定理证毕. 数值模拟和实例分析 1 计算 2 本节通过局部二次算法(LQA)来计算惩罚拟似然估计 [ 5 5 ] .调整参数 λ n 的选择方式, 采用Wang等 [ 19 19 ] 所建议的贝叶斯信息准则(BIC). 假设初值 β n (0) , 若 β nj (0) 非常接近于0, 则令 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257602&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257602&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257602&type=middle =0.惩罚函数 p λ n (| β nj |)可由下式进行逼近 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257511&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257511&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257511&type=middle 令 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257606&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257606&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257606&type=middle 是最优化算法第 k 步的迭代值，由牛顿迭代法得到下述迭代公式 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257514&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257514&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257514&type=middle 其中 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257608&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257608&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257608&type=middle 是一个对角矩阵，其第( j , j )个元素是 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257613&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257613&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257613&type=middle .当 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257617&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257617&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257617&type=middle 时迭代终止，本文中我们设定 τ =1.0×10 －4 . 数值模拟 2 本节通过数值模拟验证所提估计的有限样本性质.为了进行比较，在模拟中分别给出SCAD、Lasso、Bridge和Oracle 4种估计，这里Oracle估计是指事先已知道模型参数中哪些是零系数.用( C , I )来评估变量选择的结果，其中 C 表示参数零系数被正确选择为零的平均个数， I 表示参数非零系数错误估计为零的平均个数.同时，利用模型误差 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257621&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257621&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257621&type=middle http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257624&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257624&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257624&type=middle 来比较不同方法的表现. 考虑下面的广义线性模型： http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257518&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257518&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257518&type=middle 其中 g ( x )=2cos x , V ( x )=exp( x ), σ =0.8, β 是一个 p n 维的向量，且 β j = j ，1≤ j ≤4, β j =0, 5≤ j ≤ p n . X i 来自于 p n 维的多元正态分布，均值为零，协方差矩阵的第( i , j )个元素为 r i － j , 其中 r =0.2或0.8, http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257627&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257627&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257627&type=middle .考虑3种样本情况 n =200, p n =17； n =400, p n =21和 n =600, p n =24. 模拟的结果列在 表 1 表 1 和 表 2 表 2 中.由 表 1 表 1 和 表 2 表 2 可以得到：1) 4种方法下非零系数的估计都接近真值，SCAD和Bridge方法得到的估计比较接近，均优于Lasso，是具有Oracle性质的; 2)当样本量增大时，SCAD估计的MRME比Bridge小，但是均比Lasso小很多; 3) 4种方法对于非零系数的选择都比较好，几乎没有出现错误选择的情况. ( n , p n ) r 方法 http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257630&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257630&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257630&type=middle http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257633&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257633&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257633&type=middle http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257636&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257636&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257636&type=middle http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257639&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257639&type=small http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=4257639&type=middle 8 (200, 17) 4 0.2 Lasso 0.987 8 1.983 2 2.984 6 3.971 7 Bridge 0.995 5 1.997 1 3.001 7 4.001 1 SCAD 0.995 4 2.001 7 2.999 7 3.999 2 Oracle 1.000 2 1.999 9 3.002 9 3.999 5 4 0.8 Lasso 0.998 7 1.985 2 2.988 7 3.993 2 Bridge 0.994 6 2.007 6 2.995 4 4.001 0 SCAD 1.000 1 2.001 6 2.997 7 4.000 3 Oracle 1.004 0 2.000 5 2.994 1 3.999 7 8 (400, 21) 4 0.2 Lasso 0.983 8 1.986 5 2.982 5 3.986 0 Bridge 0.995 0 1.998 4 2.995 2 3.999 6 SCAD 1.000 8 1.999 0 2.998 7 4.001 2 Oracle 1.000 4 1.999 1 2.999 8 4.000 6 4 0.8 Lasso 1.001 1 1.997 2 3.000 1 3.990 4 Bridge 0.994 6 1.998 1 3.000 8 4.002 3 SCAD 1.000 6 1.997 7 2.999 8 4.001 7 Oracle 1.001 0 1.999 0 2.998 4 3.999 2 8 (600, 24) 4 0.2 Lasso 0.981 5 1.990 1 2.988 6 3.985 0 Bridge 0.994 9 1.997 5 2.997 5 3.998 1 SCAD 0.999 5 2.000 1 2.999 9 4.001 5 Oracle 1.000 5 2.000 9 3.001 8 3.999 3 4 0.8 Lasso 0.997 0 2.001 2 2.998 7 3.994 7 Bridge 0.997 3 2.002 6 2.998 8 3.999 2 SCAD 0.997 8 2.002 1 3.000 1 3.998 2 Oracle 0.999 9 1.999 8 2.999 6 3.998 9 MRME这一列代表着ME与MEL比值的中位数，其中ME是选择模型的模型误差，MEL是在全模型下无惩罚函数普通极大似然估计的模型误差 2 ( n , p n ) 2 方法 3 r =0.2 2 3 r =0.8 MRME C I MRME C I 4 (200, 17) Lasso 0.649 1 10.53 4 0 0.531 7 10.382 0 Bridge 0.335 7 12.46 0 0.299 7 12.228 0 Oracle 0.221 4 13 0 0.178 9 13 0 SCAD 0.322 8 12.614 0 0.210 6 12.35 0 4 (400, 21) Lasso 0.533 6 14.388 0 0.673 9 13.638 0 Bridge 0.245 2 16.732 0 0.284 6 16.692 0 Oracle 0.175 8 17 0 0.186 8 17 0 SCAD 0.208 7 16.896 0 0.212 8 16.776 0 4 (600, 24) Lasso 0.394 6 17.134 0 0.443 3 17.312 0 Bridge 0.170 3 19.944 0 0.193 2 19.874 0 Oracle 0.125 9 20 0 0.165 4 20 0 SCAD 0.152 4 19.976 0 0.176 7 19.654 0 实例分析 2 这里讨论的是南非心脏病的数据 [ 20 20 ] .这组数据来自南非西开普的三个农村地区.研究的目的在于建立高发病率地区缺血性心脏病风险因素的强度.这组数据包含了462个样本，其中有9个不同的变量和1个响应变量.响应变量是用“presence”来表示的，它表示是否患有心脏病.把响应变量0~1化，利用Logistic模型拟合这组数据.在拟合中不仅采用了前述的3种方法，而且也通过普通极大似然估计(MLE)和最优子集方法对这组数据进行变量选择和参数估计.这里Bridge的上指数取 γ =0.05，调整参数 λ n 由BIC准则进行选择.SCAD，Bridge，Lasso的 λ n 分别为0.482 1，0.024 1，0.007 8.估计结果及标准差在 表 3 表 3 中给出. 2 变量 2 MLE 2 2 Lasso 2 2 Bridge 2 2 SCAD 2 2 最优子集 估计值 标准差 估计值 标准差 估计值 标准差 估计值 标准差 估计值 标准差 sbp 0.004 5 0.006 2 0 0 0 0 0 0 0 0 tob 0.076 4 0.028 4 0.079 0 0.010 0 0.076 9 0.010 6 0.079 5 0.001 0 0.089 9 0.028 6 ldl 0.205 4 0.068 3 0.170 5 0.024 3 0.161 1 0.026 9 0.163 5 0.024 7 0.183 4 0.058 7 adi 0.004 8 0.031 9 0.016 7 0.013 1 0 0 0 0 0 0 fam 0.939 5 0.249 9 0.928 0 0.096 0 0.873 3 0.092 0 0.912 6 0.089 6 0.931 4 0.246 7 typ 0.035 6 0.014 0 0.039 0 0.004 9 0.035 1 0.005 1 0.037 2 0.004 7 0.048 3 0.013 3 obe -0.052 6 0.048 5 -0.058 0 0.02 0 0 0 0 0 0 0 alc 0.005 8 0.005 0 0 0 0 0 0 0 0 0 age 0.050 7 0.013 6 0.045 3 0.005 0 0.051 1 0.004 0 0.049 9 0.004 0 0.051 4 0.011 0 从 表 3 表 3 中可以看到：1) SCAD和Bridge成功实现了最简模型且接近于最优子集的选择结果，而SCAD的标准差略小于Bridge估计, Lasso选择的模型明显较复杂; 2)用于分析数据的5种变量选择和估计方法得到的非零系数的估计也相差不大，且与文献[ 20 20 ]中的估计结果也是一致的. 结论 1 本文主要研究在高维广义线性模型中，利用拟似然SCAD方法进行变量选择和参数估计的问题.在不要求响应变量分布的情况下，给出了参数个数发散时，拟似然SCAD估计的渐近统计性质.通过数值模拟和实例分析，得到了拟似然SCAD方法在参数估计和变量选择方面都具有良好的有限样本表现. 参考文献 NELDER J A, WEDDERBURN R W M. Generalized linear models [J]. Journal of the Royal Statistical Society: Series A, 1972, 135 (3): 370-384. WEDDERBURN R W M. Quasi-likelihood functions, generalized linear models and Gauss-Newton method [J]. Biometrika, 1974, 61 (3): 439-447. DOI: 10.2307/2334725. CHIOU J M, MULLER H G. Nonparametric quasi-likelihood [J]. The Annals of Statistics, 1999, 27 (1): 36-64. CHEN X, CHEN X R. Adaptive quasi-likelihood estimate in generalized linear models [J]. Science in China Series A: Mathematics, 2005, 48 (6): 829-846. FAN J Q, LI R Z. Variable selection via nonconcave penalized likelihood and its Oracle properties [J]. Journal of the American Statistical Association, 2001, 96 (456): 1348-1360. FRANK I E, FRIEDMAN J H. A Statistical view of some chemometrics regression tools [J]. Technometrics, 1993, 35 (2): 109-148. TIBSHIRANI R J. Regression shrinkage and selection via the Lasso [J]. Journal of the Royal Statistical Society Series B-Statistical Methodology, 1996, 58 (1): 267-288. ZOU H. The adaptive Lasso and its Oracle properties [J]. Journal of the American Statistical Association, 2006, 101 (476): 1418-1429. ZOU H, HASTIE T. Regularization and variable selection via the elastic net [J]. Journal of the Royal Statistical Society Series B-Statistical Methodology, 2005, 67 (2): 301-320. ZOU H, ZHANG H H. On the adaptive elastic-net with a diverging number of parameters [J]. The Annals of Statistics, 2009, 37 (4): 1733- 1751. DOI: 10.1214/08-AOS625. YUAN M, LIN Y. Model selection and estimation in regression with grouped variables [J]. Journal of the Royal Statistical Society Series B-Statistical Methodology, 2006, 68 (1): 49-67. HUANG J, XIE H L. Asymptotic Oracle properties of SCAD-penalized least squares estimators [J]. Lecture Notes—Monograph Series, 2007, 55 : 149-166. KNIGHT K, FU W J. Asymptotics for Lasso-type estimators [J]. The Annals of Statistics, 2000, 28 (5): 1356-1378. HUANG J, HOROWITZ J, MA S. Asymptotic properties of Bridge estimators in sparse high-dimensional regression models [J]. The Annals of Statistics, 2008, 36 (2): 587-613. DOI: 10.1214/0090536070000 00875. WANG M Q, SONG L X, WANG X G. Bridge estimation for generalized linear models with a diverging number of parameters [J]. Statistics & Probability Letters, 2010, 80 (21): 1584-1596. XIE H L, HUANG J. SCAD-penalized regression in high-dimensional partially linear models [J]. The Annals of Statistics, 2009, 37 (2): 673-696. DOI: 10.1214/07-AOS580. FAN J Q, PENG H. Nonconcave penalized likelihood with a diverging number of parameters [J]. The Annals of Statistics, 2004, 32 (3): 928-961. DOI:10.1214/009053604000000256. FAN J Q, LÜ J C. A Selective overview of variable selection in high-dimensional feature space [J]. Statistica Sinica, 2010, 20 (1): 101-148. WANG H S, LI R Z, TSAI C. Tuning parameter selectors for the smoothly clipped absolute deviation method [J]. Biometrika, 2007, 94 (3): 553-568. DOI: 10.1093/biomet/asm053. HASTIE T, TIBSHIRANI R, FRIEDMAN J. The Elements of Statistical Learning: Data Mining, Inference and Prediction [M]. New York: Springer-Verlag, 2001. 10.14188/j.1671-8836.2018.06.008 O212.1 whdxxblxb-64-6-533 国家自然科学基金(11801346)；教育部人文社会科学研究青年基金(18YJC910014)；陕西省自然科学基础研究计划(2018JM1024)；中央高校基本科研业务费专项资金(GK201503012, GK201503015) 2018-04-05 2018-12-24 2020-11-30 版权所有© 2018 《武汉大学学报(理学版)》编辑部 Copyright © 2018 Journal of Wuhan University(Natural Science Edition). All rights reserved. 2018 武汉大学学报(理学版) 6 64