稳健Poisson模型:两水平模型与GEE模型在相对危险度或患病率比估计中的应用比较*

周舒冬 郜艳晖△ 李丽霞 张 敏 杨 翌 陈 跃

稳健Poisson模型:两水平模型与GEE模型在相对危险度或患病率比估计中的应用比较*

周舒冬1郜艳晖1△李丽霞1张 敏1杨 翌 陈 跃2

目的 在处理具有层次结构特征的非罕见结局事件资料时，比较基于稳健Poisson模型的两水平模型和GEE模型在估计RR/PR时的应用。方法 将两水平稳健Poisson模型及稳健Poisson-GEE应用到2010年欧洲社会调查资料，估计影响居民生活满意度的各因素相关的PR及95%CI，以说明两模型在理论和应用上的区别和联系。结果 稳健Poisson-GEE模型的PR估计值与稳健Poisson回归模型相同，但置信区间较宽;两水平稳健Poisson模型的PR值较GEE模型为低，显示了随机效应对解释变量的混杂作用。结论 两种方法均可处理具有层次结构特征的非罕见结局事件的RR/PR估计，但两水平模型比GEE可提供更多随机效应的信息，且易于扩展至更高水平或随机系数模型。

层次结构 非罕见结局 稳健Poisson回归 广义估计方程 相对危险度 患病率比

1．广东药学院公共卫生学院流行病与卫生统计学系，广东省分子流行病学重点实验室(510310)

2．Department of Epidemiology and Community Medicine，University of Ottawa

△通信作者:郜艳晖，E-mail:gao_yanhui@163.com

由于logistic回归模型已被广大研究者所熟识，因此很多文献习惯性地将其应用于横断面研究乃至队列研究，计算优势比(odds ratio，OR)，作为相对危险度(relative risk，RR)或患病率比(prevalence ratio，PR)的估计并给予同样解释。但当研究结局频率较高时，OR值作为RR/PR的估计严重地高估暴露因素对结局的影响〔1〕，为此统计学者建议宜用Log-binomial模型或稳健 Poisson回归模型等方法直接计算RR/PR〔2〕。Log-binomial模型属于广义线性模型，采用log链接函数，误差分布为二项分布。但该法的缺陷在于最大似然估计(maximum likelihood estimate，MLE)当参数落在所限制范围的边界，特别是模型中含有连续型协变量时，常得不到似然函数导数为零的极大值，导致模型无法收敛〔3－4〕，后有学者提出采用 COPY 算法〔5〕解决模型不收敛的问题。而稳健 Poisson回归模型〔1，6〕指定误差分布为Poisson分布，应用“三明治”法获得合理的方差估计(sandwich variance)，直接估计RR/PR时不存在收敛困难。在自变量均为分类变量时，估计结果与 Mantel-Haenszel分层分析法非常近似〔2，6〕。

和其他回归模型一样，稳健Poisson回归要求观测单位间独立。但流行病或社会学调查资料中，某些观测单位常根据某些特征聚为一类，甚至具有多水平的特征，如采用多阶段抽样的横断面研究中个体来自同一街道，街道又来自同一社区;或纵向研究中同一个体的多个部位同一指标多次测量，由此获得的数据具有明显的层次结构特征，表现为类间独立、类内相关的特性，从而违背传统回归模型的应用条件。近年来广泛应用的广义估计方程(generalized estimating equation，GEE)和多水平模型(multi-level models)是解决层次结构数据的有力工具，但基于稳健Poisson回归的相关理论和应用研究尚显不足。本研究着力解决结局变量频率较高且存在类内相关时的RR/PR估计，比较稳健Poisson-GEE模型和多水平稳健Poisson模型在该类特征资料中的应用。

模型原理与方法

1.稳健Poisson回归模型

当结局事件频率较高时，为直接估计RR/PR，对独立数据，Zou〔3〕建议使用稳健Poisson模型。设yi和Xi=(xi1，xi2，…，xiP)T分别是第i(i=1，2，…，n)个观测的二分类结局变量和P×1维解释变量向量，其关系可通过Poisson回归模型表示。

式(1)中pi=Pr(yi=1/Xi)，并假设误差分布为Poisson分布。回归系数βp表示当控制其它自变量后，第p个自变量xP每变化一个单位时log(p)的相应变化。因此，与xp相对应的相对危险性为:RR(PR)=exp(βp)。

由于Poisson分布方差等于均数，当应用到二项分布资料时，易出现过度离散(overdispersion)问题，导致参数标准误的高估，产生过宽的置信区间。为此，Cameron〔7〕建议使用稳健方差估计法，如 Huber的“三明治”方差，即:

模型(1)中参数 β =(β0，β1，…，βP)T及“三明治”方差可用准似然(quasi-likelihood)估计，在SAS中可用proc genmod实现，通过在repeated语句中用“subject=”指定个体编号变量(程序见附录)。

2.稳健Poisson-GEE模型

对非独立二分类数据，Zou〔8〕提出仍可用稳健Poisson回归估计RR/PR，这时采用广义估计方程(generalized estimated equation，GEE)的原理进行参数估计，将类内水平的相关作多余参数处理。

设yki和Xki=(xki1，xki2，…，xkiP)T分别为第k(k=1，2，…，K)类内第i(i=1，2，…，nk)个个体的二分类结局变量和P×1维解释变量向量，模型形式同(1):

式(3)中=Pr(yki=1/Xki)，回归系数βp的涵义类似(1)。根据一致性估计方程理论，对参数β=(β0，β1，…，βP)T的“得分(score)”方程为:

方程(4)的解即为参数β的一致估计，方差矩阵形式同(2)，Var()=A－1BA－1

但此时，

式(5)中“三明治”方差的中间项B与(2)不同，先综合各类内个体的“得分”，再根据类别汇总。因此(5)可看作(2)在非独立数据下的自然扩展，同时校正了误差分布的错误指定与类内响应非独立对参数方差估计的影响。因此该模型也可利用SAS中的proc genmod，通过在repeated语句里使用“subject=”指定类别变量来完成(程序见附录)。

3.多水平稳健Poisson模型

与稳健Poisson-GEE模型相比，多水平稳健Poisson模型通过纳入高水平层次上的随机效应来处理层次数据中的类内相关问题。以最简单的两水平稳健Poisson方差分量模型为例，

式(6)中下标和的涵义同(3)，βp为固定效应参数，反映固定效应xp对结局概率对数的影响，且RR(PR)=exp(βp);uk为第k类的随机效应，假定来自正态分布总体，即

如数据结构中含有更高水平层次，或解释变量在高水平上存在随机效应，模型(6)还可扩展为包含更高水平上的随机效应，或随机系数模型。

模型(6)的参数估计可使用SAS中的proc glimmix来完成，调用empirical选项用于指定“三明治”方差的结构，默认缺省时的“classical”即为经典三明治方差;random语句可指定随机截距或其它高水平单位上的随机效应。(程序见附录)

实例分析

实例来源于2010年欧洲社会调查(European Social Survey)的开放数据(www．europeansocialsurvey．org)，该项目是两年一度的多国调查，覆盖欧洲各国。本研究选取2010年欧洲26国49024名居民的“生活满意度”变量作为结局变量，研究自我健康评价及家庭收支情况对结局变量的影响，构建模型时校正了年龄和性别。

49 024名居民中，对生活持满意态度的有21 979人，约占44.83%。表1描述了各因素不同状态下居民生活满意度的分布情况。

表1 各因素不同状态下居民生活满意度的分布情况

考虑到各国居民的生活满意度可能存在国家聚集性，因此建立多水平模型和利用广义估计方程的方法来拟合该数据，国家作为类别指示变量。表2显示了稳健Poisson模型、稳健Poisson-GEE模型及两水平稳健Poisson模型估计的各因素对居民生活满意度影响的PR及95%CI。三个模型结果均显示调整了性别和年龄后，自我健康评价和家庭收支对生活满意度的影响均有统计学意义。但是和稳健Poisson模型相比，GEE模型考虑了各国居民在生活满意度上的国内相关，不仅得到稳健的PR估计值，且估计的PR置信区间较稳健Poisson模型更宽，降低了犯I类错误的风险;而两水平稳健Poisson模型在模型构建时添加随机效应，估计的PR和稳健Poisson-GEE模型的结果不同，调整了性别和年龄后，自我健康评价和家庭收支对生活满意度影响的PR值均低于GEE模型结果，反映了随机效应对解释变量可能存在的混杂效应;从随机效应的方差估计值与其标准误的比值(0.1025/0.0296=3.4628)近似推断随机截距项可能有统计学意义，数据的层次结构不可忽略。此外，本文也拟合了两水平logistic方差分量模型，调整年龄和性别后，得到两解释变量的OR值均高于PR值(表2)。

表2 不同模型估计各因素对居民生活满意度影响的PR和OR及95%CI*

讨 论

当研究结局出现频率较高时，将OR习惯性地解释为RR/PR将严重高估暴露因素对结局的影响，这一问题早已引起统计学者的注意，因此提出各种直接估计RR/PR的模型和方法，如log-binomial模型和稳健Poisson回归模型，其点估计和区间估计均比logistic回归模型的OR解释起来更为合理〔9〕。本文实例使用两水平logistic模型估计自我健康评价和家庭收支对生活满意度影响的OR值均高于稳健Poisson-GEE和两水平稳健Poisson模型估计的PR值。

很多大型的流行病学调查都具有层次结构的特征，广义估计方程或者多水平模型是分析该类资料的两种相对成熟的方法〔10－11〕。本文实例欧洲社会调查项目中，各国居民由于拥有一些共同的属性或国内居民间的相互影响，个体间的研究结局并不独立。本实例拟合两水平稳健Poisson模型时随机效应参数估计结果显示数据的层次结构不可忽略，采用稳健Poisson-GEE模型估计的PR值虽与稳健Poisson模型相等，但置信区间有更宽的变化，除反映GEE模型在处理非独立数据时能有效降低I类错误的能力，某种程度上也体现出该数据具有较为明显的类内相关特征。

多水平稳健Poisson模型和稳健Poisson-GEE模型均可用于非独立数据估计解释变量对常见结局影响的RR/PR。从回归系数的性质看，稳健Poisson-GEE模型中将类内相关作多余参数，主要考虑固定效应，因此回归系数是群体效应(平均效应)的回归系数;而多水平稳健Poisson模型在随机效应条件下，估计的回归系数是个体效应的回归系数，所以两模型参数估计结果有时可能不同。如本研究中两水平模型得到的各因素PR估计值均较GEE模型为低，特别是家庭收支变量，体现了随机效应可能对解释变量的混杂作用，提示引入随机效应项后回归系数变化较大的解释变量也可能存在类内相关，并且和结局变量的随机效应存在某种程度上的关联。本例中根据目前结果提示每个国家的居民家庭收支状况可能存在国家聚集性并对研究的关联产生影响。在进一步分析中可建立稳健Poisson随机系数模型进行考察和比较(我们将另文研究)。此外，从模型的扩展角度来说，稳健Poisson-GEE模型只能处理两水平的资料;而多水平模型理论上可以处理更高水平结构的资料，且可将随机效应分解到解释变量上，进一步构建随机系数模型，而这种分析是GEE模型目前无法完成的。

除基于稳健Poisson模型的多水平和GEE模型外，在Log-binomial模型基础上扩展的多水平模型和GEE模型也可处理具有层次结构特征的常见结局资料。与多水平稳健Poisson模型和稳健Poisson-GEE模型的比较研究也是我们进一步关注的方向。

1．Barros AJ，Hirakata VN．Alternatives for logistic regression in cross-sectional studies:an empirical comparison of models that directly estimate the prevalence ratio．BMC Med Res Methodol，2003，3:21．

2．Petersen MR，Deddens JA．A comparison of two methods for estimating prevalence ratios．BMC Med Res Methodol，2008，8:9．

3．Lumley T，Kronmal R，Ma S:Relative risk regression in medical research:models，contrasts，estimators，and algorithms．UW Biostatistics working Paper Series．2006:293．http://www．bepress．com/uwbiostat/paper293．

4．Deddens JA，Petersen MR，Lei X．Estimation of prevalence ratios when proc genmod does not converge．In:Proceedings of the 28th Annual SAS Users Group International Conference，Paper 270 － 28．Cary，NC:SAS Institute Inc 2003．

5．Deddens JA，Petersen MR．Re:“Estimating the relative risk in cohort studies and clinical trials of common outcomes”．Am J Epidemiol，2004，159(2):213 －4;author reply 214 －5．

6．Zou G．A modified poisson regression approach to prospective studies with binary data．Am J Epidemiol，2004，159(7):702 －706．

7．Badi H，Baltagi．A companion to theoretical econometrics，Blackwell，Oxford(U．K．)，2001:331 －348．

8．Zou G，Allan D．Extension of the modified Poisson regression model to prospective studies with correlated binary data．Stat Methods Med Res，2011，11，8．［Epub ahead of print］

9．Lee J，Chia KS．Estimation of prevalence rate ratios for cross sectional data:an example in occupational epidemiology．Br J Ind Med，1993，50(9):861－862．

10．张春茂，李婵娟，蒋志伟，等．具有相关关系的二分类资料处理方法比较．中国卫生统计，2010，10(27):464 －467．

11．张华君，闵捷．广义估计方程与多水平模型在相关资料中的比较研究．中国卫生统计，2012，4(29):214 －217．

附录:

程度1:稳健Poisson模型SAS程序

proc genmod;

class ID;/*ID为个体的标识变量*/

model y=X1 X2 X3/d=poisson link=1og;

repeated subject=ID;

run;

程度2:稳健Poisson-GEE模型SAS程序

proc genmod;

class K;/*K为类体的标识变量*/

model y=X1 X2 X3/d=poisson link=1og;

repeated subject=K;

run;

程度3:两水平稳健Poisson方差分量模型SAS程序

proc glimmix;

MPIRICAL=CLASSICAL;

class K;/*K为类体的标识变量*/

model Y=X1 X2 X3/d=poisson link=1og solution;

random int/subject=K;

run;

A Comparison between Two-level and GEE Based on Robust Poisson Regression Models in the Estimation of Relative Risk or Prevalenc Rti

Zhou Shudong，Gao Yanhui，Li Lixia，et al．Department of Epidemiology and Biostatistics，School of Public Health，Guangdong Pharmaceutical University，Guangdong Key Laboratory f Mlecul Epidmig 510310，Gnzh

Objective To compare two-level and GEE based robust Poisson regression models in estimation of relative risk(RR)or prevalence ratio(PR)for common outcome data with intra-class correlation.MethodsTwo-level and GEE based robust Poisson regression models were compared by examing factors associated with life satisfaction using data from the 2010 European Social Survey．Prevalence ratios and 95%confidence intervals(95%CIs)were estimated.ResultsCompared to results from regular robust Poisson model，the GEE based robust Poisson model provided the same PR point estimates but wider 95%CIs．The two level robust Poisson model revealed lower point estimates，indicating potential confounding effects caused by random effects on the assocation of interest.ConclusionBoth two-level and GEE based methods are suitable for estimating relative risk or prevalence ratio for common outcomes with the hierarchical structure．The two-level model is superior when there are random effects，and can be easily extended for higher hierarchical structures．

Hierarchical structure;Non-rare outcome;Robust Poisson regression;Generalized estimating equations;Relative risk;Prevalence ratio

2010年广东省自然科学基金资助(10151022401000018)

book=686,ebook=333

(责任编辑:丁海龙)