纵向数据广泛应用在生物医药、流行病学和经济学等领域中,它通过对同一个体在不同的时间点重复测量得到, 具有组内相关、组间独立的特性。在实际的统计建模中,如果忽略这种特性,将会降低估计的效率,所以如何处理组内相关性是纵向数据研究不可避免的问题。Liang等[1]提出的广义估计方程(GEE)方法可以有效地处理组内相关问题,该方法利用含少量讨厌参数的工作相关矩阵建立估计方程,即使工作矩阵错误也可以得到回归系数的相合估计。但是GEE方法首先需要得到讨厌参数的相合估计,同时讨厌参数的相合估计可能不存在。为了提高估计的效率,Qu等[2]提出了二次推断函数法(QIF),利用某些基矩阵的线性组合来逼近工作相关矩阵的逆矩阵和避免估计讨厌参数,既使所得估计具有相合性,又提高了参数估计的效率。已有许多学者将GEE和QIF方法应用于不同的模型[3-5],但是这些应用方法大多类似于加权最小二乘估计,不具有稳健性。
在纵向数据中,由于重复测量,一个个体的突变往往会产生一系列的离群点,因此在纵向数据的研究中稳健估计是一个十分重要的问题[6-8]。为了获得更好的稳健性和提高估计的效率,Wang等[9]提出了基于指数得分函数的稳健回归估计方法, 通过引入一个调节参数来提高估计的稳健性和有效性。Lv等[10]将加权指数得分函数和GEE方法相结合,得到了纵向数据下广义线性模型参数的有效且稳健的估计(ERGEE)。但是估计讨厌参数必然会降低估计的效率。本文将指数得分函数和QIF方法相结合, 利用加权的指数得分函数有效控制协变量和因变量中离群点的影响,同时基于QIF方法建立稳健二次推断函数,得到了纵向数据下广义线性模型参数的ERQIF估计,并进一步提高了估计的效率。最后通过模拟计算进一步验证了所得估计的有限样本性质。
1 参数估计方法 1.1 ERQIF方法假定纵向数据中包含n个个体,对每个个体重复观测mi次, 总的观测次数
$ {\mu _{ij}} = E\left( {{y_{ij}}\left| {{\mathit{\boldsymbol{x}}_{ij}}} \right.} \right) = g\left( {{\mathit{\boldsymbol{\beta }}^{\rm{T}}}{\mathit{\boldsymbol{x}}_{ij}}} \right) $ |
$ Var\left( {{y_{ij}}\left| {{\mathit{\boldsymbol{x}}_{ij}}} \right.} \right) = \phi v\left( {{\mu _{ij}}} \right)\;\;\;\;\;\left( {i = 1, \cdots ,n;j = 1, \cdots ,{m_i}} \right) $ |
其中β为未知的p×1回归系数,g(·)和v(·)分别是已知的逆连接函数和方差函数,ϕ为尺度参数。不失一般性,假设不同个体间的观测是相互独立的,且同一个个体内部的观测具有一定的相关关系。
Lv等[10]基于指数得分函数提出ERGEE方法,通过求解如式(1)所示的估计方程得到了参数β的基于GEE有效且稳健的估计
$ \mathit{\boldsymbol{U}}_{\rm{n}}^{{\rm{ERGEE}}}\left( {\mathit{\boldsymbol{\beta }},\alpha } \right) = \sum\limits_{i = 1}^n {\mathit{\boldsymbol{D}}_i^{\rm{T}}\mathit{\boldsymbol{V}}_i^{ - 1}h_i^\gamma \left( {{\mu _i}\left( \mathit{\boldsymbol{\beta }} \right)} \right)} = 0 $ | (1) |
其中Di=∂μi/∂β是ni×p矩阵;Vi=Ri(α)Ai1/2,Ri(α)是依赖于参数α的mi×mi工作相关矩阵,且Ai=φdiag(v(μi1), …, v(μimi)), μi=(μi1, …, μimi)T; hiγ(μi(β))=Wi[ψγ(μi(β))-Ci(μi(β))], ψγ(μi)= ψγ(Ai-1/2(Yi-μi)), Ci(μi)=E[ψγ(μi)],ψγ(t)=
权重矩阵Wi=diag(wi1, …, wimi)用于限制协变量中离群点的影响,可取Mallow s-type函数wij为
$ {w_{ij}} = w\left( {{\mathit{\boldsymbol{x}}_{ij}}} \right) = \min \left\{ {1,{{\left( {\frac{{{b_0}}}{{{{\left( {{\mathit{\boldsymbol{x}}_{ij}} - {\mathit{\boldsymbol{m}}_x}} \right)}^{\rm{T}}}\mathit{\boldsymbol{S}}_x^{ - 1}\left( {{\mathit{\boldsymbol{x}}_{ij}} - {\mathit{\boldsymbol{m}}_x}} \right)}}} \right)}^{\rho /2}}} \right\} $ |
其中ρ≥1,b0是自由度为p的卡方分布的0.95分位点,mx和Sx分别表示关于Xij的位置参数和尺度参数的稳健估计,可以利用minimum covariance determinant(MCD)进行估计。
为了简单起见,假定重复观测次数是相同的,即mi=m<∞, 并令R(α)为公共的工作相关矩阵。当重复观测次数不同时,可采用Xue等[11]的方法加以处理。在实际计算中,为了得到
$ {\mathit{\boldsymbol{R}}^{ - 1}} \approx {a_1}{\mathit{\boldsymbol{M}}_1} + {a_2}{\mathit{\boldsymbol{M}}_2} + \cdots + {a_k}{\mathit{\boldsymbol{M}}_k} $ | (2) |
其中M1是单位矩阵,M2,…, Mk是某些已知的对称的基矩阵[2],a1, a2,…, ak是未知常数。
将式(2)代入式(1)构造扩展得分向量为
$ {{\mathit{\boldsymbol{\bar U}}}_n}\left( \mathit{\boldsymbol{\beta }} \right) = \frac{1}{n}\sum\limits_{i = 1}^n {{\mathit{\boldsymbol{U}}_i}\left( \mathit{\boldsymbol{\beta }} \right)} = \frac{1}{n}\sum\limits_{i = 1}^n \\ {\left( {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{D}}_i^{\rm{T}}\mathit{\boldsymbol{A}}_i^{ - 1/2}{\mathit{\boldsymbol{M}}_1}h_i^\gamma \left( {{\mathit{\boldsymbol{\mu }}_i}\left( \mathit{\boldsymbol{\beta }} \right)} \right)}\\ \vdots \\ {\mathit{\boldsymbol{D}}_i^{\rm{T}}\mathit{\boldsymbol{A}}_i^{ - 1/2}{\mathit{\boldsymbol{M}}_k}h_i^\gamma \left( {{\mathit{\boldsymbol{\mu }}_i}\left( \mathit{\boldsymbol{\beta }} \right)} \right)} \end{array}} \right)} $ |
由于Un(β)的维数为kp,显然高于β的维数p,因此求解方程Un(β)=0会出现过度识别问题。基于GMM估计的思想,定义有效且稳健的ERQIF估计方程如式(3)
$ {Q_n}\left( \mathit{\boldsymbol{\beta }} \right) = n\mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( \mathit{\boldsymbol{\beta }} \right)C_n^{ - 1}\left( \mathit{\boldsymbol{\beta }} \right){{\mathit{\boldsymbol{\bar U}}}_n}\left( \mathit{\boldsymbol{\beta }} \right) $ | (3) |
其中,
通过极小化目标函数Qn(·)得到β的估计
$ {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{ERQIF}}}} = \arg {\min _\mathit{\boldsymbol{\beta }}}{Q_n}\left( \mathit{\boldsymbol{\beta }} \right) $ |
利用如式(4)的Newto n-Raphson迭代算法得到估计
$ {{\mathit{\boldsymbol{\hat \beta }}}^{\left( {k + 1} \right)}} = {{\mathit{\boldsymbol{\hat \beta }}}^{\left( k \right)}} - {\left[ {{\nabla ^2}{Q_n}\left( {{{\mathit{\boldsymbol{\hat \beta }}}^{\left( k \right)}}} \right)} \right]^{ - 1}}\nabla {Q_n}\left( {{{\mathit{\boldsymbol{\hat \beta }}}^{\left( k \right)}}} \right) $ | (4) |
其中,∇Qn(β)和∇2Qn(β)分别为Qn(β)对β的一阶和二阶导数,且满足
$ \left\{ \begin{array}{l} \nabla {Q_n}\left( \mathit{\boldsymbol{\beta }} \right) = 2n\nabla \mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( \mathit{\boldsymbol{\beta }} \right)C_n^{ - 1}\left( \mathit{\boldsymbol{\beta }} \right){{\mathit{\boldsymbol{\bar U}}}_n}\left( \mathit{\boldsymbol{\beta }} \right) - \\ \;\;\;\;\;\;n\mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( \mathit{\boldsymbol{\beta }} \right)\nabla C_n^{ - 1}\left( \mathit{\boldsymbol{\beta }} \right){{\mathit{\boldsymbol{\bar U}}}_n}\left( \mathit{\boldsymbol{\beta }} \right)\\ {\nabla ^2}{Q_n}\left( \mathit{\boldsymbol{\beta }} \right) = 2n\nabla \mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( \mathit{\boldsymbol{\beta }} \right)C_n^{ - 1}\left( \mathit{\boldsymbol{\beta }} \right)\nabla {{\mathit{\boldsymbol{\bar U}}}_n}\left( \mathit{\boldsymbol{\beta }} \right) + {O_p}\left( 1 \right) \end{array} \right. $ |
设β0为参数β的真值,β所处的参数空间记为Θ。首先给出一些假设条件[4](本文中Eβ0(·)均表示随机变量关于真值β0的数学期望):
1) 参数空间Θ是可识别的,若β≠β0, 则Eβ0[U1(β)]≠0;
2) 对任意β∈Θ,Eβ0[U1(β)]存在且有限,同时关于β连续;
3) Θ是紧空间,且β0为Θ的一个内点;
4) 存在β0的某个邻域O0,当β∈O0时,有
$ {C_n}\left( \mathit{\boldsymbol{\beta }} \right)\xrightarrow{{a.s.}}{E_{{\mathit{\boldsymbol{\beta }}_0}}}\left[ {{\mathit{\boldsymbol{U}}_1}\left( \mathit{\boldsymbol{\beta }} \right)\mathit{\boldsymbol{U}}_1^{\rm{T}}\left( \mathit{\boldsymbol{\beta }} \right)} \right] \triangleq {\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_{{\mathit{\boldsymbol{\beta }}_0}}}\left( \mathit{\boldsymbol{\beta }} \right) $ |
其中,Σβ0(β)为正定矩阵,且关于β∈O0连续;
5) Un(β)关于β的一阶导数存在且连续,且当
$ \frac{{\partial {{\mathit{\boldsymbol{\bar U}}}_n}\left( {{\mathit{\boldsymbol{\beta }}^ * }} \right)}}{{\partial \mathit{\boldsymbol{\beta }}}}\xrightarrow{P}\mathit{\boldsymbol{B}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right) = {E_{{\mathit{\boldsymbol{\beta }}_0}}}\left[ {\frac{{\partial {{\mathit{\boldsymbol{\bar U}}}_n}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)}}{{\partial \mathit{\boldsymbol{\beta }}}}} \right] $ |
基于条件1)~4),根据文献[12]和文献[2]中相关定理的证明,可得定理1。
定理1 在条件1)~4)下,通过式(3)得到的
定理2 在条件1)~5)下,若Un(β)和Cn(β)的二阶导数有有限的均值和方差,则当n→∞时,有
$ \sqrt n \left( {{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{ERQIF}}}} - {\mathit{\boldsymbol{\beta }}_0}} \right)\xrightarrow{D}{N_p}\left( {0,{\mathit{\boldsymbol{J}}^{ - 1}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)} \right) $ |
即证明了所得估计的渐进正态性。其中,J(β0)=BT(β0)Σβ0-1(β0)B(β0)。
证明 因
$ {\bf{0}} = \nabla {Q_n}\left( {{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{ERQIF}}}}} \right) = \nabla {Q_n}\left( {{\mathit{\boldsymbol{\beta }}_{\rm{0}}}} \right) + {\nabla ^2}{Q_n}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)\left( {{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{ERQIF}}}} - {\mathit{\boldsymbol{\beta }}_0}} \right) $ |
其中
$ \sqrt n \left( {{{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{ERQIF}}}} - {\mathit{\boldsymbol{\beta }}_0}} \right) = - {\left( {\frac{1}{{2n}}{\nabla ^2}{Q_n}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)} \right)^{ - 1}}\frac{1}{{2\sqrt n }}\nabla {Q_n}\left( {{\mathit{\boldsymbol{\beta }}_{\rm{0}}}} \right) $ |
其中
$ \nabla {Q_n}\left( {{\mathit{\boldsymbol{\beta }}_{\rm{0}}}} \right) = 2n\nabla \mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)C_n^{ - 1}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)\nabla {{\mathit{\boldsymbol{\bar U}}}_n}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right) + {O_p}\left( 1 \right) $ |
$ {\nabla ^2}{Q_n}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right) = 2n\nabla \mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)C_n^{ - 1}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)\nabla {{\mathit{\boldsymbol{\bar U}}}_n}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right) + {O_p}\left( 1 \right) $ |
由条件4)和5)以及
$ \begin{gathered} \frac{1}{{2\sqrt n }}{Q_n}\left( {{\mathit{\boldsymbol{\beta }}_{\rm{0}}}} \right) = \nabla \mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)C_n^{ - 1}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)\sqrt n {{\mathit{\boldsymbol{\bar U}}}_n}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right) + \hfill \\ {Q_p}\left( {{n^{ - 1/2}}} \right) = {\mathit{\boldsymbol{B}}^{\rm{T}}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_{{\mathit{\boldsymbol{\beta }}_0}}^{ - 1}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)\sqrt n {{\mathit{\boldsymbol{\bar U}}}_n}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right) + {Q_p}\left( {{n^{ - 1/2}}} \right)\xrightarrow{D} \hfill \\ {N_p}\left( {{\bf{0}},\mathit{\boldsymbol{J}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right)} \right)\frac{1}{{2n}}{\nabla ^2}{Q_n}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right) = \nabla \mathit{\boldsymbol{\bar U}}_n^{\rm{T}}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)C_n^{ - 1}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)\nabla {{\mathit{\boldsymbol{\bar U}}}_n} \hfill \\ \left( {\mathit{\boldsymbol{\tilde \beta }}} \right) + {O_p}\left( 1 \right)\xrightarrow{P}\mathit{\boldsymbol{J}}\left( {{\mathit{\boldsymbol{\beta }}_0}} \right) \hfill \\ \end{gathered} $ |
因此可得当n→∞时,
推论1 在定理2的条件下,可得
$ {{\mathit{\boldsymbol{\hat B}}}_n} = \sum\limits_{i = 1}^n {{{\left( {{{\mathit{\boldsymbol{\hat B}}}_{i1}}, \cdots ,{{\mathit{\boldsymbol{\hat B}}}_{ik}}} \right)}^{\rm{T}}}} $ |
$ {{\mathit{\boldsymbol{ \boldsymbol{\hat \varSigma} }}}_n} = \sum\limits_{i = 1}^n {{{\left( {{{\mathit{\boldsymbol{ \boldsymbol{\hat \varSigma} }}}_{i1}}, \cdots ,{{\mathit{\boldsymbol{ \boldsymbol{\hat \varSigma} }}}_{ik}}} \right)}^{\rm{T}}}} $ |
$ \mathit{\boldsymbol{ \boldsymbol{\hat \varSigma} }}_{is}^{\rm{T}} = \mathit{\boldsymbol{D}}_i^{\rm{T}}\mathit{\boldsymbol{A}}_i^{ - 1/2}{\mathit{\boldsymbol{M}}_s}h_i^\gamma \left( {{\mu _i}\left( {\mathit{\boldsymbol{\hat \beta }}} \right)} \right)h_i^\gamma {\left( {{\mu _i}\left( {\mathit{\boldsymbol{\hat \beta }}} \right)} \right)^{\rm{T}}}{\mathit{\boldsymbol{M}}_s}\mathit{\boldsymbol{A}}_i^{ - 1/2}{\mathit{\boldsymbol{D}}_i} $ |
$ \mathit{\boldsymbol{\hat B}}_{is}^{\rm{T}} = \mathit{\boldsymbol{D}}_i^{\rm{T}}\mathit{\boldsymbol{A}}_i^{ - 1/2}{\mathit{\boldsymbol{M}}_s}\dot h_i^\gamma \left( {{\mu _i}\left( {\mathit{\boldsymbol{\hat \beta }}} \right)} \right){\mathit{\boldsymbol{D}}_i}\left( {s = 1, \cdots ,k} \right) $ |
$ \dot h_i^\gamma \left( {{\mathit{\boldsymbol{\mu }}_i}\left( {\mathit{\boldsymbol{\hat \beta }}} \right)} \right) = \partial h_i^\gamma \left( {{\mathit{\boldsymbol{\mu }}_i}\left( {\mathit{\boldsymbol{\hat \beta }}} \right)} \right)/\partial {\mathit{\boldsymbol{\mu }}_i}\left| {_{{\mathit{\boldsymbol{\mu }}_i} = {\mathit{\boldsymbol{\mu }}_i}\left( {\mathit{\boldsymbol{\hat \beta }}} \right)}} \right. $ |
$ \mathit{\boldsymbol{\hat \beta }} = {{\mathit{\boldsymbol{\hat \beta }}}_{{\rm{ERQIF}}}} $ |
要通过式(3)得到ERQIF估计
$ \hat \phi = {\left\{ {1.483{\rm{medi}}\left\{ {\left| {{{\hat \xi }_{ij}} - {\rm{medi}}\left( {{{\hat \xi }_{ij}}} \right)} \right|} \right\}} \right\}^2} $ |
其中
$ {\gamma _{{\rm{opt}}}} = \arg \mathop {\min }\limits_\gamma \left( {\det \left( {{\rm{\hat Cov}}\left( {\mathit{\boldsymbol{\tilde \beta }}} \right)} \right)} \right) $ |
其中(
通过一个迭代算法得到ERQIF估计,具体步骤如下:
1) 给定β的初始估计
2) 利用当前估计
3) 利用公式(4)迭代求解
4) 若
在数值计算中,利用独立工作相关矩阵的GEE方法[1]得到初始估计
为了得到估计的有限样本性质,本文通过随机模拟考察数据出现离群点时估计的稳健性和工作相关矩阵被正确或错误指定对估计效率的影响。用本文提出的ERQIF方法与GEE、ERGEE方法分别计算3种参数估计的样本偏差(Bias)和标准差(SD)并进行比较。设样本容量n=100,重复观测次数m=5, 模拟次数为500。在模拟计算中,分别指定组内工作相关矩阵Ri(α)为Exch和AR(1)两种结构。
3.1 线性模型当响应变量是连续型变量时,假设模型满足:
$ {y_{ij}} = \sum\limits_{k = 1}^3 {x_{ij}^{\left( k \right)}{\mathit{\boldsymbol{\beta }}_{0k}}} + {\varepsilon _{ij}}\left( {i = 1, \cdots ,n;j = 1, \cdots ,m} \right) $ | (5) |
其中,β01=0.7,β02=0.7,β03=-0.4;xij(k)~N(0, 1);Cov(xij(k), xij(l))=0.5|k-l|,k, l=1, 2, 3。随机误差εi=(εi1, …, εi5)T的真实工作相关矩阵是Exch结构,相关系数α=0.7。
分别考虑两种情况:
① εij~N(0, 1);
② εij~0.9N(0, 1)+0.1N(0, 100)。
用3种方法计算得到两种情况下两种结构的估计值,结果分别如表 1、表 2所示。
从表 1可以看出,不管相关结构是否正确指定,3种方法的Bias差异不大,说明3种方法均能得到参数的相合估计。由GEE和ERGEE方法相比较可知,ERQIF方法的SD均较小,特别是当工作相关结构被错误指定时SD有较大降低,说明ERQIF方法可避免对讨厌参数的估计,同时提高估计的效率。
从表 2可以看出,GEE方法的SD值受离群点的影响很大,说明该方法不能得到参数的稳健估计。通过引入加权指数得分函数得到的ERGEE和ERQIF方法可以有效控制离群点的影响,从而保证了估计的稳健性。同时从SD的结果可以看出,不管相关结构是否正确指定,ERQIF方法均优于ERGEE方法,其估计效率有了显著提高。
3.2 逻辑回归模型当响应变量是0-1型变量时,假设模型满足:
$ \ln \left( {{\mu _{ij}}/\left( {1 - {\mu _{ij}}} \right)} \right) = {x_{ij}}{\beta _0}\left( {i = 1, \cdots ,n;j = 1, \cdots ,m} \right) $ | (6) |
其中xij独立地服从[0, 1]上的均匀分布,β0=0.4,Yi的真实相关矩阵Ri(α)为可交换结构,且α=0.2。在模拟计算中,分别指定Ri(α)为Exch和AR(1)两种结构;参考文献[14]生成不同结构下的yij;为了研究稳健性,对部分数据进行扰动,即将Xij(i=1, …, 10)的每个值都减去5,计算得到3种方法下参数的偏差和标准差如表 3所示。
从表 3可以看出,3.2节和3.1节的模拟结果相似。在不同情况下,ERQIF方法在Bias与SD上均比其他方法占优势,特别当数据有扰动时,ERGEE与ERQIF方法比非稳健GEE方法的Bias有显著降低;且ERQIF较ERGEE的SD更小,且在工作相关结构被错误指定的情况下,ERQIF也能得到更稳健的估计,估计效率更高。
总之,从3.1节和3.2节的模拟结果可以得出,当数据出现离群点或工作相关矩阵被错误指定时,ERQIF方法均有更好的稳健性和有效性。
4 结束语本文将加权指数得分函数和QIF方法相结合,得到了纵向数据下广义线性模型的ERQIF估计,证明了估计的渐近性质,并通过数值计算得到了稳健的模拟结果。由于本文仅研究了模型参数估计问题,未来可望进一步研究各种半参数模型的稳健估计和变量选择问题,并对估计的稳健性进行理论分析。
[1] |
Liang K Y, Zeger S L. Longitudinal data analysis using generalized linear models[J]. Biometrika, 1986, 73(1): 13-22. DOI:10.1093/biomet/73.1.13 |
[2] |
Qu A, Lindsay B G, Li B. Improving generalised estimating equations using quadratic inference functions[J]. Biometrika, 2000, 87(4): 823-836. DOI:10.1093/biomet/87.4.823 |
[3] |
Lian H, Liang H, Wang L. Generalized additive partial linear models for clustered data with diverging number of covariates using GEE[J]. Statistica Sinica, 2014, 24: 173-196. |
[4] |
Qu A, Li R Z. Quadratic inference functions for varying-coefficient models with longitudinal data[J]. Biometrics, 2006, 62(2): 379-391. DOI:10.1111/j.1541-0420.2005.00490.x |
[5] |
Tian R Q, Xue L G, Liu C L. Penalized quadratic inference functions for semiparametric varying coefficient partially linear models with longitudinal data[J]. Journal of Multivariate Analysis, 2014, 132: 94-110. DOI:10.1016/j.jmva.2014.07.015 |
[6] |
Fan Y L, Qin G Y, Zhu Z Y. Variable selection in robust regression models for longitudinal data[J]. Journal of Multivariate Analysis, 2012, 109: 156-167. DOI:10.1016/j.jmva.2012.03.007 |
[7] |
Zheng X Y, Fung W K, Zhu Z Y. Robust estimation in joint mean-covariance regression model for longitudinal data[J]. Annals of the Institute of Statistical Mathematics, 2013, 65(4): 617-638. DOI:10.1007/s10463-012-0383-8 |
[8] |
Qin G Y, Zhu Z Y, Fung W K. Robust estimation of generalized partially linear model for longitudinal data with dropouts[J]. Annals of the Institute of Statistical Mathematics, 2016, 68(5): 977-1000. DOI:10.1007/s10463-015-0519-8 |
[9] |
Wang X Q, Jiang Y L, Huang M, et al. Robust variable selection with exponential squared loss[J]. Journal of the American Statistical Association, 2013, 108(502): 632-643. DOI:10.1080/01621459.2013.766613 |
[10] |
Lv J, Yang H, Guo C H. An efficient and robust variable selection method for longitudinal generalized linear models[J]. Computational Statistics & Data Analysis, 2015, 82: 74-88. |
[11] |
Xue L, Qu A N, Zhou J H. Consistent model selection for marginal generalized additive model for correlated data[J]. Journal of the American Statistical Association, 2010, 105(492): 1518-1530. DOI:10.1198/jasa.2010.tm10128 |
[12] |
Hansen L P. Large sample properties of generalized method of moments estimators[J]. Econometrica, 1982, 50(4): 1029-1054. DOI:10.2307/1912775 |
[13] |
Yao W X, Lindsay B G, Li R Z. Local modal regression[J]. Journal of Nonparametric Statistics, 2012, 24(3): 647-663. DOI:10.1080/10485252.2012.678848 |
[14] |
Oman S D. Easily simulated multivariate binary distributions with given positive and negative correlations[J]. Computational Statistics & Data Analysis, 2009, 53(4): 999-1005. |