・166・ 中国卫生统计2017年2月第34卷第1期 混合效应模型中的方差成分检验 曾 平 赵杨 陈峰。
在很多科学问题中需要在混合效应模型框架下对 随机效应方差成分(暂记为 )进行假设检验 j,也 即检验 :r =0。除直接科学兴趣外,许多间接医学 问题也能转化为对方差成分的检验。例如,为判断在 惩罚样条回归中是参数模型还是非参数模型更合适, Claeskens 首先建立混合效应模型,将模型选择问题 转化为对随机效应方差成分是否等于零的假设检验问 题,最后通过限制性似然比检验 :丁 =0。其他研究 者也用同样的方法处理过类似问题Ⅲ8 。然而,方差 成分为非负参数,对方差成分的假设检验是非标准的: 在 下 位于参数空间边缘。由于这种限制,常用 的渐近 无效分布对似然比统计量不再成立¨ ’ 。 混合效应模型中的方差成分检验吸引了广泛研究兴 趣l卜8. 。本文主要介绍方差成分检验的似然比方 法,方差成分得分检验或Wald检验可参考其他相关 文献 13,20,21]。 线性混合效应模型的方差成分检验 设y为连续变量, =(X --, )和Z=(z ,…, Zq)分别为,z×P和n×q的矩阵, 为样本量。在广 泛意义上】,、 和Z既适于非独立数据结构 I2 ,如 重复测量数据或聚群数据 ,也适于独立数据情 形 ¨,如惩罚样条回归。更加详细的例子可参 考上述文献。y、 和z之间的关系由线性混合效应 模型 刻画 p q Y=∑Xjflj+∑Zkb + , (1) 用矩阵的形式模型(1)可表达为 Y= +z易+s,b~Ⅳ(0,r K),占~Ⅳ(0, In) (2) 式中 为固定效应,b为随机效应, 为未知方差成 分, 为已知q×q矩阵描述随机效应间相互关系,L :国家自然科学基金项目(81402765,81473070,81373102);国家统计 局全国统计科学研究项目(2014LY1 12);江苏省教育厅高校哲学社会 科学研究基金项目(2013SJD790032,2013SJB790059);南京医科大学公 共卫生学院优秀博士论文培育项目 1.徐州医科大学公共卫生学院流行病与卫生统计学教研室(221004) 2.南京医科大学公共卫生学院生物统计学系 为n维单位阵, 为残差方差,与b不相关。模型中 y的期望和方差分别为 E(y)= ,Var(1,)=.1-2ZKZ + I =or ^ (3) 其中 ^=AZKZ +J ,A=,/2/o-。。设参数为0=( ,A, or ),模型(2)的对数似然函数为 ’ r 1 (0):一 1{,z log(o" )+ (y一 ) (y—
L )+logl 1) (4) 似然比统计量定义为
=2sup[ (0)一 (0)] (5) 1.渐近混合卡方分布 显然检验/40:b=0等价于检验/40: =0。如前 文所述虽然定义似然比统计量并不难,但是 : =0 是非标准形式检验,渐近 无效分布不再成立。这 里,与非标准形式检验对应是标准形式检验,后者指参 数位于参数空间内部 。Self和Liang…、Liang和 Selfl3 以及Stram和Lee_2 假设:如y能够划分为许多 独立同分布的亚向量、亚向量维度趋于无穷且研究为 均匀设计,那么职 渐近服从0.5Q +0.5ox , 为 退化到零的点分布, 是自由度为1的卡方分布。然 而,Crainiceanu和Ruppert 认为,Stram和Lee的假设 条件在实际中很难满足,并显示0.5ox +0.5ox 会导 致保守结果。为验证这一点,图1给出了模拟例子。 对特殊纵向数据,Pinheiro和Bates 发现非等比例卡 方混合分布更加合适,如0.65x +0.35x 。但对更一 般的似然比方差成分检验,Crainiceanu和Ruppert_8 证 明 和 的混合比例取决于具体数据,没有满足所 有情形的混合比例参数H 。 数据根据模型(1)模拟 ,设丁 =0,然后对比 L尺 和0.5ox +0.5ox 的分位数。0.5ox +0.5ox 的分位数总是大于对应 的分位数(散点和虚线 重叠表示两分布一致),意味着在0.50Xo+0.5ox 提 供的临界值过大。因此如果采用0.5 +0.5 作为 似然比统计量的理论无效分布将产生保守结果,降低 统计效能和导致假阴性结果。 2.基于谱分解的模拟分布 为了获得似然比统计量的无效分布,Crainiceanu Chinese Jouma1 ofHealth Statistics.Feb.2017.Vo1.34.No.1 模拟计算的似然比统计量 模拟计算的似然比统计量 模拟计算的似然比统计量模拟计算的似然比统计量 图1模拟计算的似然比统计量和0.5 +0.50x 的分位数。 和Ruppert在谱分解基础上提出了一种基于模拟抽样 残差,由 条件下的模型计算获得;h。 为矩阵H= 算法 。假设模型(2)中A已知,可获得卢和 的极x(x X) X 对角线上元素。在 : =0条件下,模型 大似 估计值 (1)退化为简单线性模型Y:∑7 ,q-Eo假设台和 卢=( lX) y, =(y P 1p ̄y)/ 分别是JB和 z在 模型下的 值。基于重抽样
, 获得似然比统计量无效分布可执行为:①根据原始数 带人对数似然函数,获得剖面似然函数 据j[)获得似然比统计量的观察值,记做职 如;②获得
(A)=÷{,z 1og( l,)+logf f} (7) 重抽样数据 ,采用相同的模型和估计算法获得关 于D 似然比统计量,记做职 ;③重复第二步 次, 一f RT r 一 得到B个 ;则似然比统计量统计量的无效分布由
~ 职 组成,其对应的P值由LRn 7 等于或大于职 , 、 …,:LL, ,7 /’J H Lti L J n 、 J ’~/ J …n
{,zlog[¨ g ¨ )) 比例 p: 。具体而言在 是矩阵 z ZK 的特征根,其中, 第二步中如果通过分布N( , )获得重抽样】, 、同 (A) 2,Dn(A) 丁_ 2+ 时保持 和Z不变,则为参数b。。tstrap检验;如果Y
n叩 =X +台 、且台 通过有放回方式从台。随机获得,那 T, ,0、 厶+1“ 么为非参数bootstrap检验;如果Y = o+ 、且
是 z PoZK 特征根,P。=, 一X(X X) X 。则 LRT 以分布收敛于.厂LK 。 基于上述谱分解可建立模拟算法获得职 的无 效分布。需要注意的是,该分布是有限样本下的精确 分布而非大样本近似分布。谱分解模拟算法只与随机 效应维度有关、与样本量无关,因此在实际中具有较高 计算速度 。。。 。该算法获得的无效分布能够保证似 然比检验有效控制I型错误率,与其他方差成分检验 方法相比具有更高统计效能 , 。’ 。当Stram和 Lee假设条件满足时,基于模拟得到的精确分布弱收 敛于0.50 ̄ +0.50 ̄ 。 3.重抽样方法 统计分析中当面对复杂假设检验时诉诸于重抽样 方法是一种常用策略 。重抽样避免了数学推导, 概念上简单并且容易操作,属于计算密集型统计技术。 Fitzmaurice等 以及Lee和Braun_4 提出了似然比 方差成分的permutation检验,Faraway 建立了似然 比方差成分的参数bootstrap检验。其他研究者对似 然比方差成分提出了类似的重抽样检验方法¨ J。
设 。 =  ̄/1一h。 称为调整残差,称台。为原始
通过无放回方式从台。随机获得,那么为permutation检 验。由此可见,参数bootstrap检验、非参数bootstrap 检验和permutation检验的主要差别在于如何产生新 的数据。这三种重抽样方法往往获得相似的结果,但 相对参数和非参数bootstrap检验,permutation检验假 设条件更少 , 。另外进行重抽样的对象是调整 残差 而非原始残差 ,这是因为 。是同方差,而 是异方差 。一般而言,实际数据分析 取值1000 左右就足够 ’ .3 51 。
广义线性混合效应模型的方差成分检验 1.1ogistic混合效应模型和PQL算法 目前几乎所有似然比方差成分检验都是在线性混 合效应模型框架下完成,即都是针对连续应变量的;而 对离散应变量似然比方差成分检验的研究还很少。以 logistic回归为例,接下来讨论广义混合效应模型的似 然比方差成分检验。设y为二分类应变量, 和z同 前。y、 和z之间的关系通过logistic混合效应模 型 卜 描述
(】,),log【 )= + +Zb, ・168・ b—N(O,T2K) (10) 为y期望。感兴趣的检验是14o:西=0,同样等价于 Hn:丁 =0。模型(1O)的对数似然函数可表示为:L n (/3,丁 )=∑J.Ef(Yi l/3,b)f(b l丁 )] 。似然比检验
i=l 统计量
2[L 。(J8, )一L(/3, =0)] (11) 然而,由于对数似然函数涉及积分运算,实际中除极有 限和简单情况外,大多数时候都无法精确计算 L 。(卢, )和获得精确的对数似然值。因此不能直 接通过式(11)进行似然比检验。 logistic混合模型的参数估计算法包括数值积 分 ,惩罚拟似然(penalized quasi—likelihood,PQL)和 边际拟似然(marginal quasi—likelihood,MQL) ,伪似 然估计 川,Monte Carlo估计 和Gibbs抽样 。 数值积分只适用于低维度的积分运算,一般很难处理 超过五维的积分 ;Monte Carlo估计和Gibbs抽样的 不足在于计算繁重。PQL和MQL能够通过近似方法 有效处理高维积分问题,计算相对简单,且能够通过重 复调用已有的线性混合模型程序执行。对logistic混 合模型进行方差成分的似然比检验至少包含以下的三 个困难 ]:与线性混合模型不同,通常不能获得关于 logistic混合模型对数函数的闭型表达式,因此精确计 算logistic混合模型的对数似然值和似然比统计量变 得很困难;由于涉及高维计算运算,很难获得logistic 混合模型参数的精确估计值;即使能够有效估计logis— cic混合模型和计算其对数似然比统计量,如何获得该 统计量的无效分布也不清楚。基于logistic混合模型 已有的理论发展,本文仅仅专注其方差成分的似然比 检验。 2.基于重抽样方法的伪似然比方差成分检验 为了处理上面的问题,可以通过基于重抽样的伪 似然比方法进行方差成分检验 。具体步骤如下:① 采用PQL算法估计logistic混合模型参数;PQL算法 可通过R软件中的函数glmmPQL进行,该函数包含 在MASS软件包中 ;②当PQL算法收敛时,有工作 应变量 Y = +Z西+8,b~N(0,T2K), —N(0,R) (12) 和伪似然函数 1 . . ( )=一 一[1og l D I+(y 一X卢) D一 (y 一X )] 厶 +C (13) c为一个常数,D=r ZKZ +R, =( D ) D y 。③将 (.r )看做y 的对数似然函数,计算伪 似然比统计量LRT =2sup o[ vH!(7- )一 ,H0 ( )];④通过重抽样的方法获得伪似然比统计量 LRT 的无效分布。 中国卫生统计2017年2月第34卷第1期 其他问题和可能的研究方向 目前,由于以下几方面的原因,似然比方差成分检 验主要集中在线性混合效应模型并且模型中只包含一 个方差成分 ’他’" :①能够较为容易估计线性混合 效应模型的参数和计算其对数似然值,进而计算似然 比统计量;②多个方差成分共存时,不存在简单的对数 似然函数谱分解形式,因此很难获得似然比统计量的 无效分布;③由于可能的高维积分运算,一般难以精确 估计广义混合效应模型的参数和计算其对数似然值, 以及获得似然比统计量的无效分布;④重抽样方法属 于计算密集型统计方法,难以大规模运用。 在标准的参数假设检验中似然比检验被证明是在 小样本时最优的 。事实上,由于似然比检验充分利 用了H。条件下的模型和H 条件下的模型,相对得分 检验和Wald检验,似然比检验被证明在非标准参数 假设检验也具有更高统计效能 .3 U'如J。因此,推广和 发展似然比检验以适合更加广泛的情形具有理论和实 际意义。具体而言,以下几方面问题值得进一步探索: ①在线性混合效应模型中当存在多个方差成分而仅对 其中某个成分感兴趣时,Greven 等建议可以首先从 应变量中消除冗余随机效应;此时,尚包括应变量的残 差和待检验方差成分,然后按照前文只包含单个方差 成分的似然比检验进行统计分析。然而,这样做的缺 点在于,没有考虑到随机效应之间可能存在的相关,而 且在模型估计多个方差成分时还可能存在数值问题, 如模型不收敛或计算不稳定。②在线性混合效应模型 中假设残差是独立的;当存在相关协方差结构,即 ~ N(0, ( )), 为一般形式协方差矩阵、包含未知参 数向量 ;在这种情况下可先估计残差协方差,在此基 础上建立伪数据,然后将伪数据当做原始数据按照线 性混合效应模型框架下相似的方法进行似然比方差成 分检验 。但这种方法在协方差矩阵 的选择以及 其适用性等方面需要更多的研究。③前文在广义线性 混合效应模型框架下建立的伪似然比方差成分检验其 有效性取决于PQL算法;PQL算法本身属于有偏估 计;虽然提出了校正PQL偏倚的算法 卜 J,有偏和无