第3章混杂因素的处理 3.1混杂因素处理的基本问题 3.1.1混杂因素与混杂偏倚 在实际的医学问题研究中,为了探讨某处理因素(如某种药物、治疗方法)与结果(如生存时间、智力恢复)的关系,需要设立处理组和对照组进行比较,比较的前提是二者具有可比性,也就是说二者除了具备所研究的因素之外,其他因素应该尽可能在两组间是一样的。如果被研究人群中存在一个或多个既与研究结果有关,又与处理因素有关的因素,可能会掩盖或夸大所研究的处理因素与结果之间的联系(Cochran,1973),这些因素称为混杂因素(confounding factors)。这种由混杂因素造成的偏差称为混杂偏倚(confounding bias)。 例如,在对非随机分组的观察性研究中,研究对象被分配到各组的机会往往取决于研究对象的基线特征(如年龄、性别、并发症、病情严重程度以及分级等),而这些基线特征又会对治疗的结果产生影响。这时,利用各种统计检验直接比较各组间的治疗结局(如治愈率)的差异是不恰当的,因为结果不同可能是由于基线特征不同而导致的,如病程、年龄等。 混杂偏倚的本质是既与所研究的处理因素有关,又与研究结局有关的混杂因素在处理组和对照组中分布不均或不平衡造成的(耿直,2004)。 3.1.2混杂因素处理的常用方法 对混杂因素处理的目的就是要控制混杂偏倚,传统控制混杂偏倚的方法包括在研究设计阶段进行匹配,限制一定条件的研究对象进入; 在数据分析阶段使用标准化法,或按照混杂因素分层,以及采用多因素数学模型进行调整等。但这些方法都有一定的局限性: 如匹配设计、分层分析需要考虑的混杂因素都不能太多,否则由于匹配的混杂因素太多会导致找不到合适的匹配对象,分层因素太多会导致所分层数太多而使每个层内的分析样本量太少而无法分析,多因素数学模型较为常用,但往往需要注意数学模型的适用条件。 而倾向评分法则不受以上限制,它可以在分析和设计阶段有效平衡非随机对照研究中的混杂偏倚,使研究结果接近随机对照研究的结果。 相较于上述传统方法的局限性,2000年后兴起的倾向评分法则具有以下优点: (1) 适用于混杂因素很多,而结局变量发生率很低的情况; (2) 通过倾向值调整组间的混杂因素,使临床观察性数据可以成为循证医学的诊疗证据,而这些数据获取成本低且量大,更能够反映医疗实践中实际存在的疾病谱; (3) 在无法实现随机化的药物临床试验以及医疗器械临床试验中,可以通过倾向评分方法,平衡组间的混杂因素; 另外,在意向性治疗(intention to treat,ITT)分析中,综合考虑脱落病例的基线水平与结局发生情况,采用倾向评分方法对其完成临床试验的条件概率进行估计并纳入ITT分析,与传统分析中对于脱落病例只采用末次观察推进法(last observation carried forward)进行数据接转完成ITT分析相比,具有更强的外推性。 倾向评分法是一种可以在分析和设计阶段有效平衡非随机对照研究中的混杂偏倚,使研究结果接近随机对照研究结果的混杂因素控制方法。 3.1.3混杂因素处理的倾向评分法 倾向评分(propensity score)是Rosenbaum和Rubin于20世纪80年代提出的一种方法。它将考虑到的混杂因素综合为一个变量(倾向评分或称倾向值),通过平衡两组的倾向评分有效地均衡各个混杂因素的分布,达到一种类似随机化的状态,实现控制混杂偏倚(Rosenbaum,1983)。它可以在分析和设计阶段有效平衡非随机对照研究中的混杂偏倚,使研究结果接近随机对照研究的效果(Paul,1984)。2000年之后, 倾向评分法日益受到人们的关注。国际上越来越多的研究者将倾向评分法应用到流行病学、健康服务研究、经济学以及社会科学等许多领域。 1. 倾向评分法的思路 倾向评分法是计算每个研究对象的倾向值,通过匹配或其他一些方法使得处理组和对照组的倾向值同质(严格相等实际上是很难做到的),最后基于匹配样本进行统计分析; 也可以不匹配,使用倾向值作为权重进行多元分析,或者使用倾向评分进行回归调整分析(郭申阳,2012)。 倾向评分需要确定的内容: 结局变量、混杂因素、处理变量的不同及其意义。混杂变量又可称为协变量,代表着与预测变量和结局变量均有关联的外来变量。混杂变量的存在往往混淆预测变量和结局变量的联系,而且混杂因素在暴露组与对照组的分布往往是不均衡的。对混杂因素的处理,是倾向评分的重点之一应尽可能将与预测变量和结局变量都有关的所有混杂因素均包括在模型中。 2. 倾向评分法适用的资料类型 倾向评分法适合于所有非随机化研究的资料,或者说存在混杂偏倚的研究资料的处理。主要包括下面一些资料类型。 1) 观察性研究资料 观察性研究资料包括现况研究、病例对照研究以及队列研究等。在观察性研究中,处理变量不是人为给予的,而是自然存在的,每个个体暴露于处理变量的机会往往不是随机的,因此,不可避免地存在混杂偏倚。但观察性研究包括的研究对象反映实际人群特征和暴露状况,因此其结果的外推性强。随着信息技术的发展,观察性数据无论是数量还是准确性都在不断增加,因此,倾向评分法在观察性研究中具有广阔的应用前景。 2) 非随机干预研究资料 非随机干预研究也是流行病学研究中常用的研究方法,是介于队列研究和随机对照研究中间的一种研究方法。与观察性研究不同的是,在非随机干预研究中处理变量是人为给予的,但每个个体接受处理变量与否不是随机的,因此也存在混杂偏倚。 3) 随机对照方案失败的研究资料 某些随机对照研究,如药物临床试验,由于研究对象未按试验方案接受处理而导致随机化分组方案失败。对于上述非随机对照研究或随机化分组方案失败的随机对照研究,与结局变量和处理变量相关联的混杂因素的分布在比较的两组间往往是不均衡的,而忽视这种不均衡直接对研究因素的效果进行估计将可能得到有偏倚的结果。 3. 倾向评分法的基本步骤 倾向评分法是利用估计的倾向分值,使两组之间在混杂因素上达到平衡,以控制混杂对结局造成的偏倚,更真实地评价处理因素对结局的影响。 运用倾向评分方法控制混杂,从处理过程看有两大类。一类是估计倾向值后进行匹配再分析,主要有三步: 倾向值估计、匹配和匹配后分析; 另一类是估计倾向值后不采用匹配,直接利用倾向值进行分析,即在倾向值估计后进行再抽样,尽量使得处理组和对照组在倾向值上相似,达到平衡,控制选择偏倚; 或者利用倾向值作为权重进行多元分析。运用倾向评分控制选择偏倚,有不同的方法适用,实际应用中可以根据自己研究的需要选择。 倾向值的估计在实际研究中很有意义,如在观察性研究中,特别关注在使用了某中药注射剂后患者的不良反应,期望用该中药处理后相对于未用该中药产生的不良反应会更少。结局值y0对每个接受处理的个体而言是无法观测到的,所以E(y0|z=1)必须从对照组的数据中估计得到。Rosenbaum和Rubin(1983)曾讨论,如果多个协变量在处理组和对照组之间存在差异,则此估计值是有偏的,从而ATE1的估计也会有偏,需要考虑利用倾向评分调节这些差异。 3.2倾向值的估计 利用倾向评分平衡混杂因素,首先需要估计倾向值,其目的是寻找影响选择偏差的观测变量,也就是要找到使倾向评分估计最佳的一组条件变量。倾向值估计方法有很多种,包括logistic回归、probit回归以及判别分析等。可以根据分组变量和协变量的不同类型选用不同的函数,如分组变量为二分类变量时通常选用logistic回归模型、probit回归模型或者判别分析,其中如果协变量均为正态分布的变量,可以选用判别分析估计各个研究对象的倾向值; 如果协变量中包含有二分类变量,选用logistic回归方法; 若为多分类无序时选用多分类logistic模型。 3.2.1倾向评分的基本原理 1. 倾向评分值 Rosenbaum和Rubin对倾向评分值亦称倾向值的定义为(Rosenbaum,1983): 倾向评分值是在给定某些协变量的条件下,研究对象进入某一特定干预如处理组的条件概率,如式(3.1)所示。 e(xi)=P(Wi=1|Xi=xi)(3.1) 其中,e(xi)表示研究对象i的倾向值; Wi=1表示i进入处理组; Wi=0表示i进入对照组; Xi=xi表示控制了除i处理因素以外的所有已知的混杂因素。 Rosenbaum和Rubin推导并证明了一系列反映倾向值性质的原理。 倾向评分是指在给定一组观测到的特征变量(observed characteristics variable)条件下,一个研究个体被分配到处理组而非对照组的概率。Rosenbaum和Rubin 在1983年证明,以倾向评分为条件,所有观测到的协变量(covariates variable)独立于组别分配,且在大样本的情况下,两组的协变量分布几乎相同且不会混淆估计的处理效应。 估计倾向评分需要首先确定混杂因素,寻找合适的可能会导致研究结果产生偏倚的混杂因素,将这些混杂因素以协变量的形式放到模型中估计出倾向值。这一阶段的主要难点在于确定影响研究结果的混杂因素并进一步为倾向值模型中的变量设定函数形式。一般来说,混杂因素需具备以下三个条件: (1) 必须是所研究结局的独立危险因素,且在两比较组间分布不均衡。 (2) 必须与研究因素有关,但不是这一研究因素的结局。 (3) 一定不是研究因素与所研究结局因果链上的中间变量。 2. 平均处理效应 假定总体人群中每个个体都有两个潜在的结局值(potential values for any outcome): 一个是个体被分配或接受处理条件时的结局值y1,一个是个体被分配或接受对照条件时的结局值y0。对每个个体而言,这两个值仅有一个被观察到,另一个虚拟结局是不可能被观察到的。总体人群的处理效应定义为E(y1)-E(y0)。然而,通常感兴趣的效应是接受了某个特殊过程或特殊类型处理的处理效应,即所谓的处理组平均处理效应(average treatment effect on the treated),记为ATE1。令z为处理分配指标,如果个体接受处理,则z=1,否则z=0。E(y1|z=1)是处理组个体接受处理条件的平均结局值,E(y0|z=1)是处理组个体接受对照条件(comparison condition)的平均结局值。那么,处理组平均处理效应如式(3.2)。 ATE1=E(y1|z=1)-E(y0|z=1)(3.2) 倾向评分定义为P(X)=P(Z=1|X),其中X表示对分组会产生影响的混杂变量。它是评价两组间混杂变量或特征变量均衡性的近似函数。 为估计放入处理组的条件下,对照组中个体的结局指标的平均值E(Y0|Z=1),需要给对照组中个体i赋予权重如式(3.3): wi=P(Xi)1-P(Xi)(3.3) 其中Xi表示个体i对应混杂变量的具体取值。如果个体i来自处理组,可以观察到Yi=y1i; 如果个体i来自对照组,则可以观察到Yi=y0i。从而,对照组中个体被观察到结局指标的平均值的估计值如式(3.4): E^(Y0|Z=1)=∑i∈Cwiyi∑i∈Cwi(3.4) 其中,i∈C表示对照组中第i个观测值。 令NT表示处理组的个体数,i∈T表示属于此组的个体i。对这些个体样本均值的估计值为式(3.5): E^(Y1|Z=1)=∑i∈TyiNT(3.5) 估计的处理效应为式(3.6): EATE=E^(Y1|Z=1)-E^0(Y0|Z=1)(3.6) 3.2.2二分类logistic回归估计倾向值 1. 估计方法 当存在两种分组状态(即处理和对照)时,接受处理的条件概率通过二分类logistic回归进行估计。接受处理的条件概率表达如式(3.7): P(Wi|Xi=xi)=E(Wi)=exiβi1+exiβi=11+e-xiβi(3.7) 其中,Wi是第i个对象的二分类状态,即如果研究对象处于处理组,Wi=1; 若处于对照组,Wi=0,Xi代表各协变量,βi是各协变量相应的参数。 式(3.7)经过logit变换后可写为式(3.8): lnP1-P=xiβi(3.8) 式(3.8)中,P代表公式(3.7)中的P(Wi)。式(3.8)可以采用最大似然法进行估计,找到βi的估计值。估计得到的βi是使样本观测再现的可能性最大化时的logistic回归系数。将这些回归系数代入式(3.7)中,得到每一研究对象接受处理的预测概率,即估计的倾向值。 2. 拟合检验 利用数据建模的时候,需要评估所建模型对数据的拟合情况。目前已有很多统计量评估模型的拟合优度。 (1) 皮尔逊卡方拟合优度检验(Person chisquare goodnessoffit test) 该检验检测对logistic反应函数的偏离程度。当统计量的值较大,即相应的P值较小时,表明该logistic反应函数是不恰当的。该检验对较小的偏离并不敏感。 (2) 所有系数的χ2检验(Chisquare test of all coefficients) 该检验是一个似然比检验,使用对数似然比构造卡方统计量进行检验。 H0: 增加模型参数是多余的 H1: 增加模型参数不是多余的 检验统计量 LR=-2logL(β^R,σ^2R)L(β^UR,σ^2UR)~χ2(m)(3.9) 其中,R是带限制条件的模型,即仅含截距项的模型; UR是不带限制条件的模型,即增加参数后的模型; m是自由度,即限制参数的个数。 如果计算得到的LR值大于χ2α(m),拒绝原假设,即拒绝增加的所有系数都等于0的假设,增加模型的参数不是多余的。似然比检验不适宜小样本情况。 (3) HosmerLemeshow拟合优度检验(HosmerLemeshow goodnessoffit test) 该检验先将样本分为较小的组,如g个组,然后计算由2g个观测频数和估计的期望频数所组成的列联表的皮尔逊χ2检验统计量。如果统计量小于χ21-α(g-2)意味着模型拟合效果好。其中,α是给定的显著性水平,通常取0.05,g-2是自由度。该检验对样本量很敏感,在分组简化数据的过程中,可能会因为一小部分个体数据点造成对拟合的重大偏离。在判断模型拟合情况之前,最好对个体残差和有关诊断统计量进行分析。 (4) 虚拟R2(Pseudo R2) 由于logistic回归是非线性估计,无法得到类似线性回归中的决定系数即拟合优度R2,但可以有类似的虚拟R2。这些虚拟R2包括调整R2、计数R2、调整的计数R2。一般来说,虚拟R2取值较高表明拟合效果较好。但需注意: 虚拟R2不能用于比较不同数据间的拟合效果,只能用于比较同一数据的同一结果的多个模型拟合效果。 3.2.3GBM法 Logistic回归估计的倾向值是否正确,在很大程度上依赖于所选入的协变量是否以正确的函数形式纳入模型。若所选的协变量未以正确形式纳入模型,估计得到倾向值的正确性值得怀疑。函数形式的设定通常是主观的,这使得logistic回归法应用受到限制。McCaffrey等(2004)发展出一种程序,这种程序使用一般化加速建模(generalized boosted modeling, GBM)寻找两个组在协变量上的最佳平衡。 1. GBM法的基本思路 GBM倾向评分法是基于广义的Boosted回归(generalized boosted regression)估计倾向评分的方法。Boosting是一种自动的、数据自适应算法,可以被用来估计所关注的处理变量和大量协变量之间的非线性关系。Boosting用现有的方法建立一些精度不高的弱分类器或回归树,将这些预测效果粗糙的模型组合起来,组合成一个整体分类系统,达到改善整体模型性能的效果。GBM的Boosting思想主要体现在把许多简单的函数组合在一起去估计关于大量协变量的光滑函数。虽然每个简单函数缺乏平滑性且它们只能很差地近似目标函数,但是如果把这些简单函数加在一起则可以近似一个光滑函数(smooth function)。这就好比把一连串的线段连接在一起就可以去近似光滑曲线一样(Friedman,2000)。GBM实现过程中,每个简单的函数就是一个有限深度的回归树(regression tree with limited depth)。GBM也被称为一般化加速回归,是通过回归树的方式拟合多个模型,将每个模型预测值加以合并的方法。它不需要事先设定协变量的函数形式,也就不需要提供那些如β^0、β^1、β^2等的估计值。 Friedman 在2000年(Friedman,2000)证明在预测误差方面,Boosting方法优于其他方法。许多Boosting算法的变种已经出现在机器学习和统计计算文献中,比如AdaBoost算法、GBM算法、LogitBoost算法以及Gradient Boosting Machine算法等。特别是当模型中存在大量协变量,且它们与处理变量之间线性、非线性或交互效应等函数形式无法确定时,此方法最具优势。 2. GBM法的特点 回归树方法不需要设定预测变量的函数形式,回归树的结果不会因为自变量的一对一转换而变,例如,无论使用年龄、年龄的对数还是年龄的平方作为研究对象的特征,都会获得完全相同的倾向值。GBM法不产生估计的回归系数,但会给出影响力(influence),代表每一个输入变量所解释的对数似然函数的百分比,所有需要控制变量的影响力的总和为100%。例如,假设有3个变量: 年龄、病程以及处理前的风险因素,GBM法的输出结果可能显示年龄的影响力是20%,病程的影响力是30%,风险因素的影响力是50%。这说明,处理前风险因素对估计的对数似然函数的贡献最大,即该因素在两组间分布最为不平衡。 如用倾向评分法探究服用某药物的患者与未服用该药物患者治疗“肝硬化、病毒性肝炎”的疗效差异。 暴露组: 选取某药物数据库中的患者,并且用药前7天有谷丙转氨酶检查,且检查提示异常,停药后7天内有谷丙转氨酶检查的患肝硬化或病毒性肝炎的人群,同时,需满足双环醇片用药天数15天以上、住院天数30天以内的要求。最后选取基线ALT为40~200的患者,暴露组共251人。 非暴露组: 选取数据库中肝硬化、病毒性肝炎的患者,住院天数30天以内,15天以上,未使用该药物,且住院期间有两次及以上的谷丙转氨酶检查,第一次检查提示谷丙转氨酶异常者,最后选取基线ALT为40~200的患者,非暴露组共5988人。 使用GBM法得到估计的倾向评分值,并根据各个协变量对模型对数似然函数的贡献,对它们在处理分配上的重要程度进行测量和排序。图31选取了相对影响程度前十位的协变量进行展示,表31列出全部协变量的相对影响程度。 图31相对影响程度前十位的协变量 表31混杂因素对处理分配的影响程度表(全部协变量) 协变量 重 要 程 度 出院科室 60.25404 住院天数 20.06045 复方茵陈注射液 8.298267 年龄 6.361182续表 协变量 重 要 程 度 人血白蛋白 1.564307 病危天数 1.251326 总费用 0.869639 胰岛素 0.638881 入院方式 0.404482 电解质代谢紊乱 0.231666 职业 0.065757 病重天数 0 入院科室 0 婚姻 0 费别 0 入院病情 0 性别 0 阿德福韦酯 0 奥美拉唑 0 多烯磷脂酰胆碱 0 复方氨基酸注射液 0 还原型谷胱甘肽 0 螺内酯片 0 乳果糖口服溶液 0 维生素K1 0 胸腺素 0 呋塞米 0 腹腔感染 0 腹腔积液 0 肝良性肿瘤 0 乙肝肝硬化 0 原发性肝癌 0 3.3匹配和匹配后分析 3.3.1倾向评分匹配的原理 获得倾向值后,可以使用这些值匹配处理组个体和对照组个体。传统的匹配只能针对较少的协变量进行一对一匹配,当存在高维数据时,并不适用。倾向评分匹配可以综合多个变量的影响,克服传统匹配的缺点。通过计算对照组、处理组中