因果推断:未测量混杂因素的敏感性分析-T249

发布时间:2021-05-19 阅读 3988

Stata连享会   主页 || 视频 || 推文 || 知乎

温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。

New! lianxh 命令发布了:
随时搜索推文、Stata 资源。安装命令如下:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh

课程详情 https://gitee.com/lianxh/Course

课程主页 https://gitee.com/lianxh/Course

⛳ Stata 系列推文:

PDF下载 - 推文合集

作者: 李胜胜 (安徽大学)
邮箱: lisheng2@foxmail.com

Source: Ariel Linden, Maya B. Mathur, Tyler J. VanderWeele, 2020, Conducting sensitivity analysis for unmeasured confounding in observational studies using E-values: The evalue package, Stata Journal, 20(1): 162–175. -PDF-, -PDF2-


目录


1. 引言

使用观察数据进行评估时的一个基本问题是,未测量的混杂因素可能会被误认为是处理效应。因此,研究人员在进行分析时,努力调整影响所有这些关联的变量。然而,在观察性研究中,不太可能获得所有潜在混杂变量的数据。因此,应进行评估后的敏感性分析,以评估未测量的混杂因素与处理分配之间,以及与结果之间的关系有多强,从而解释观察到的处理效应。

在本文中,主要介绍 evalue 软件包,该软件包使用 VanderWeele 和 Ding (2017 年,内科年鉴167:268-274) 提出的方法对观察研究中未测量的混杂因素进行敏感性分析。 evalue 报告 E 值,定义为在风险比率等级上的最小关联强度,这是未经衡量的混杂因素与处理和结果之间具有的最低关联强度。完全解释特定的处理-结果关联,具体取决于所测量的协变量。 evalue 计算几种常见结果类型的点估计值 (以及可选的置信区间) 的 E 值,包括风险比率、常见或罕见结果的胜算比、常见或罕见结果的风险比、结果标准化平均差异和风险差异。

evalue 报告的 E 值被定义为风险比 (RR) 等级的最小关联强度,该风险比等级取决于测量的协变量,而未衡量的混杂因素在处理和结果方面都需要具备解释其关联。与大多数其他敏感性分析方法 (混杂因素的指定强度是否足以解释效应估计) 相反,E 值关注的是混杂因素关联的大小,该混杂因素可能产生与观察到处理结果关联相同的混杂偏差。

E 值方法和公式适用于多种混杂因素。然后将混杂关联的大小解释为可产生的最大 RRs,该 RRs 可比较整个未测量混杂因子组的任何两个值 (以测量的协变量为条件) 。研究者不需要选择混杂变量 (或指定它们的混杂关联) ,而只是报告一个未测量的混杂因素与处理结果的关联程度,以解释估计效应。然后,读者或其他研究人员可能会评估这种程度的混杂因素是否合理。

2. 估计方法

E 值是根据 RR 等级计算得到,因此必须将除 RR 之外的统计模型的结果转换为 RR 等级。以下将介绍计算各种模型的 E 值所涉及的方法。

2.1 RR 和比率的 E 值

在 RR 等级上计算任何结果类型的 E 值 (及其最接近于零的置信区间) 的基本公式如下 (VanderWeele 和 Ding,2017) :

2.2 胜算比的 E 值

当结果相对罕见时 (例如,随访结束时的患病率 <15%) ,胜算比 (OR) 接近 RR,因此可以使用基本的 E 值公式 (见第 2.1 节) 。在一项病例对照研究中,结果只需要在基础人群中罕见,而不需要在研究样本中罕见 (当随访结束时患病率大约 >85% 时,也可以做相同的考虑,因为变量编码可以简单地反转,即反过来编码) 。如果结果并不罕见 (随访结束时患病率在 15% 至 85% 之间) ,可用 OR 的平方根代替 RR,可获得近似的 E 值 (VanderWeele 2017,即第 2.1 节中给出 E 值公式中的 RROR。请注意,当结果很罕见时,OR 变换得到的近似值很差,因此应使用“罕见”结果假设下的计算。然而,当结果的概率在 15% 到 85% 之间时,OR 近似值则非常有效 (Ding 和 VanderWeele,2016) 。

2.3 危险比的 E 值

当上述结果相对罕见时,应使用基本 E 值公式 (见第2.1节) 。当结果普遍时,可通过应用第 2.1 节 E 值公式中的近似值 RR(10.5HR)/(10.51/HR) 获得近似的 E 值 (VanderWeele,2017) 。

2.4 标准化均值差异的 E 值

使用标准化效应大小 d (结果变量的平均值除以结果混合标准偏差[SD]) 和 SD 的标准误差,通过在 E 值公式中应用近似值 RRe[0.91×d] 可以获得近似的 E 值 (Lipsey 和 Wilson, 2001; VanderWeele, 2017; Linden, 2019) 。类似地,可以通过使用近似值(e[0.91×d1.78×SD],e[0.91×d+1.78×SD])可以获得 RR 的近似置信区间 (CI) 。但此方法依赖于额外的假设和近似值,面对这种情况,其他敏感性分析技术已被开发出来。但它们还是通常需要别的假设,变量也不一定有相应的 E 值。

2.5 风险差异的 E 值

如果已处理和未处理的调整后风险为 p1 和 p0,则可以通过在 E 值公式中将 RR 替换为p1/p0来获得 E 值。当然,风险差异 (RD) 等级上 CI 的 E 值很复杂,需要计算几个衡量值,并使用网格搜索来找到相应的偏差因子。当转换为 RR 等级时,该偏差因子会引起置信下限的 E值发生改变。或者,如果结果概率 p1 和 p0 较大时 (例如,如果它们在 0.20 和 0.80 之间) ,则可以使用 2.4 节中给出的连续结果差异的近似方法。

2.6 非零假设的 E 值

到目前为止,已经描述了如何计算 E 值来评估未测量的混杂因素与处理和结果之间的最小关联强度,从而将点估计值或 CI 的一个极限移到零。然而,类似的程序可用于评估将估计值移到 RR 其他值所需的两个混杂相关的最小幅度。如果有一个观察到的 RR 的 RR,并且想要评估将估计值转移到其他值 RRT 所需的两个关联的最小强度,那么首先要取两个值的比率 RR/RRT,然后将第 2.1 节中给出的 E 值公式应用于该比值。

3. evalue 命令

3.1 命令安装

ssc install evalue,replace

3.2 语法结构

  • E-value for risk ratio (RR) and rate ratio:
    evalue rr point_estimate [, lcl(#) ucl(#) true(#) figure[(twoway_options)]]
  • E-value for odds ratio (OR):
    evalue or point_estimate [, lcl(#) ucl(#) true(#) common figure[(twoway_options)]]
  • E-value for hazard ratio (HR):
    evalue hr point_estimate [, lcl(#) ucl(#) true(#) common figure[(twoway_options)]]
  • E-value for standardized mean difference (SMD):
    evalue smd point_estimate [, se(#) true(#) figure[(twoway_options)]]
  • E-value for risk difference (RD):
    evalue rd #a #b #c #d [, true(#) level(#) grid(#) figure[(twoway_options)]]

evalue rd 的语法中,#a 是曝光的患病个体数 (E=1,D=1) ;#b 是曝光的非患病个体数 (E=1,D=0) ;#c 是未曝光的患病个体数 (E=0,D=1) ;#D 是未曝光的非患病个体数 (E=0,D=0) 。如果观察到的 RD 为负,则应首先反转曝光的编码 (反过来编码) 以产生正的 RD。

[options]选项

  • lcl(#):点估计周围的置信区间 (CI) 的下限,可用于 RR,OR 和 HR 模型
  • ucl(#):点估计周围的置信区间 (CI) 的上限,可用于 RR,OR 和 HR 模型
  • se(#):SMD 点估计的标准误差;适用于SMD模型
  • true(#):将观察点估计值 (非空值) 转移到处理效应值,对于RD和SMD模型,默认值为 true(0);对于所有其他模型,默认值为 true(1)
  • common:指定 OR 和 HR 模型随访结束时的结果患病率在 15% 和 85% 之间。,默认值为 rare
  • level(#):为 RD 模型设置置信水平,默认值为 level(95)
  • grid(#):用于 RD 模型 E 值的网格搜索间距 (公差) ,默认值为网格 (0.0001)
  • figure[(twoway_options)]:产生一条曲线,描绘出用于点估计的联合混杂因素关联的范围,以及适用时的 CI;所有twoway_options中的选项都可用

4. Stata 范例

evalue 的设计类似于 immediate 命令 (参见[U]19 immediate commands) ,因为它获取点估计和作为参数类型化的 CI,而不是从存储在内存中的数据中获取。这使得研究人员能够对文献中公布的结果 (通常仅包括点估计和 CI) 进行敏感性分析,并作为使用个体水平数据后估计命令。在下面的示例中,我们将说明这两种情况。

4.1 RR 的 E 值

在一项以人群为基础的病例对照研究中,Victora 等人 (1987) 研究了母乳喂养与婴儿呼吸道感染死亡之间的关系。在对年龄、出生体重、社会地位、母亲受教育程度和家庭收入进行控制后,研究人员发现,仅用配方奶粉喂养的婴儿死于呼吸道感染的可能性是纯母乳喂养婴儿的 3.9 倍 (95% CI,1.8-8.7) 。文献中,研究人员控制了社会经济地位的指标,但没有控制吸烟,吸烟可能会减少母乳喂养和增加呼吸系统死亡的风险。

要计算此相对风险的 E 值,我们输入点估计值 (3.9) 以及置信区间的上限和下限 (分别为 1.8 和 8.7) 。我们还应用了 figure 选项。

. evalue rr 3.9, lcl(1.8) ucl(8.7) figure
   E-value (point estimate): 7.263
   E-value (CI): 3.000

如代码输出结果和图所示,点估计的 E 值为7.26。该 E 值可解释为:观察到的风险比 3.9 可以通过与处理和结果相关的未测量的混杂因素来解释,其风险比为 7.2 倍,高于和超过测量的混杂因素,但较弱的混杂因素不能做到这一点 (VanderWeele 和 Ding,2017) 。同样,置信下限 (即,最接近于零的置信界限) 的 E 值为 3.0,这可以解释为:与呼吸系统死亡和母乳喂养相关的未经测量的混杂因素 (风险比均为3.0倍) 可以解释较低的置信度,但较弱的混杂因素无法解释。因此,这些 E 值的因果关系证据看起来相当有力,因为需要大量不可估的混杂因素来将观察到的关联或其 CI 降低为零 (VanderWeele 和 Ding,2017) 。

4.2 OR 的 E 值

在此示例中,我们通过不指定 common 选项来针对极少数的结果率 (即小于15%的情况) 执行敏感性分析。我们使用 Moorman (2008) 等人的一项研究得出的估算值,表明,在母乳喂养 6 至 12 个月的绝经前妇女中,发生卵巢癌的几率比未母乳喂养的妇女低 0.5 倍 (95% CI:[0.3,0.8]) 。

. evalue or 0.5, lcl(0.3) ucl(0.8) 
   E-value (point estimate): 3.414
   E-value (CI): 1.809

如结果所示,点估计的 E 值为 3.41,CI为 1.81。点估计似乎是适度稳健的,但与母乳喂养和卵巢癌这种程度的混杂关联可能会使 CI 变为零。

4.3 HR 的 E 值

在此示例中,我们使用 Stata 官方命令 stcox (请参阅 [ST] stcox) 通过 Cox 回归估算 HR 和 CI,然后将这些估算值传递给 evalue hr 。数据由 stcox 提供,针对癌症药物试验的 48 名参与者,在此期间 64.6% 的患者死亡 (常见结果) 。在这 48 位患者中,有 28 位接受了治疗,有 20 位接受了安慰剂。参与者的年龄从 47 岁到 67 岁不等。我们拟合了一个模型来评估治疗效果,并根据患者年龄进行了调整。

. webuse drugtr, clear
(Patient Survival in Drug Trial)

. stcox drug age

         failure _d:  died
   analysis time _t:  studytime

Iteration 0:   log likelihood = -99.911448
Iteration 1:   log likelihood = -83.551879
Iteration 2:   log likelihood = -83.324009
Iteration 3:   log likelihood = -83.323546
Refining estimates:
Iteration 0:   log likelihood = -83.323546

Cox regression -- Breslow method for ties

No. of subjects =           48                  Number of obs    =          48
No. of failures =           31
Time at risk    =          744
                                                LR chi2(2)       =       33.18
Log likelihood  =   -83.323546                  Prob > chi2      =      0.0000

------------------------------------------------------------------------------
          _t | Haz. Ratio   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
        drug |   .1048772   .0477017    -4.96   0.000     .0430057    .2557622
         age |   1.120325   .0417711     3.05   0.002     1.041375     1.20526
------------------------------------------------------------------------------

. evalue hr .1048772, lcl(.0430057) ucl(.2557622) common
   E-value (point estimate): 8.245
   E-value (CI): 4.483

stcox 的输出结果中可以看到,与安慰剂相比,该药物的危害性更低 (HR = 0.105; 95%CI:[0.04,0.26]) ,因此存活时间更长 (P <0.001) 。HR 的 E 值和置信上限分别为 8.25 和 4.48,这表明即使考虑到现有的控制混杂因素较差,因果关系的证据也相当强。

4.4 SMD 的 E 值

在这个例子中,我们说明了如何将从线性回归模型中得出的治疗效果估计值与 SMD 的二元曝光进行拟合,然后如何传递该估计值 (及其标准误差) 来评估 SMD。我们首先使用 Cattaneo (2010) 的数据子集,进行 regress (见 [R] reg ) 来评估母亲怀孕期间吸烟状况对婴儿出生体重的影响。使用母亲的年龄、教育水平、婚姻状况、以及这个婴儿是否是母亲的第一个孩子作为协变量。

. webuse cattaneo2, clear
(Excerpt from Cattaneo (2010) Journal of Econometrics 155: 138-154)

. regress bweight mbsmoke mage medu mmarried fbaby

      Source |       SS           df       MS      Number of obs   =     4,642
-------------+----------------------------------   F(5, 4636)      =     56.14
       Model |  88775811.1         5  17755162.2   Prob > F        =    0.0000
    Residual |  1.4661e+09     4,636  316244.268   R-squared       =    0.0571
-------------+----------------------------------   Adj R-squared   =    0.0561
       Total |  1.5549e+09     4,641  335032.156   Root MSE        =    562.36

------------------------------------------------------------------------------
     bweight |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
     mbsmoke |   -224.422   22.07908   -10.16   0.000    -267.7075   -181.1365
        mage |   .4146478   1.821294     0.23   0.820    -3.155956    3.985252
        medu |   7.914262   3.753915     2.11   0.035     .5548011    15.27372
    mmarried |   159.1408   20.92324     7.61   0.000     118.1213    200.1603
       fbaby |  -52.95197    17.8149    -2.97   0.003    -87.87765    -18.0263
       _cons |   3203.872   53.88544    59.46   0.000     3098.231    3309.513
------------------------------------------------------------------------------

我们发现吸烟母亲所生婴儿的平均出生体重比不吸烟母亲所生婴儿的平均出生体重少 224 克。为了将该点估计值转换为 SMD,我们使用了社区贡献的命令 esizeregi (Linden 2019) 。为了计算 SMD,还需要检索因变量 (bweight) 的 SD,并获得每组的观测值 (mbsmoke) 。

summarize bweight
tabulate mbsmoke
esizeregi -224.422, sdy(578.8196) n1(864) n2(3778)
--------------------------------------------------------------------
        Effect Size |  Estimate   Std. Err.     [95% Conf. Interval]
--------------------+-----------------------------------------------
          Cohen's d | -0.383382   0.037920     -0.457704   -0.309060
--------------------------------------------------------------------

接下来,将估计值和标准误插入 evalue smd 中以计算 E 值。

. evalue smd -0.383382, se(0.037920)
   E-value (point estimate): 2.187
   E-value (CI): 1.981

4.5 RD 的 E 值

在此示例中,我们说明了如何为二进制结果的 RD 计算 E 值。与其他 evalue 子命令不同,用户输入点估计和 CI 后进行模型估计, evalue rd 的用户输入 2 个 2×2 表的四个值 (请参见 epitab cscsi) 。其中,#a 是曝光的患病个体数 (E=1,D=1) ;#b 是曝光的非患病个体数 (E=1,D=0) ;#c是未曝光的患病个体数 (E=0,D=1) ;#D 是未曝光的非患病个体数 (E=0,D=0) 。

Hammond 和 Horn (1958a,b) 报告了一项针对 187,783 名男性的队列研究,结果表明吸烟与肺癌死亡之间存在关联,其中 42% 被分类为有定期吸烟 (曝光) 与 58% 其他人 (不吸烟或仅偶尔吸烟) 。我们可以使用 csi 计算 RD 和 CI,从而得出 RD 估计值和 CI 为 0.00456(95%CI:[0.00405,0.00507])。

csi 397 51 78557 108778

使用 evalue rd ,输入这些数据以计算 E 值:

. evalue rd 397 78557 51 108778
   E-value (point estimate): 20.947
   E-value (CI): 15.957

这些结果可以解释如下:在观察到的风险差异 RD = 0.00456 的情况下,一个与经常吸烟和肺癌死亡相关的未测量的混杂因素,其风险比为 20.95 倍,高于和超过测量的混杂因素,这可以解释估计值,但较弱的混杂因素不能解释;移动 CI 包括空值,这是一种未测量的混杂因素,即与经常吸烟和肺癌死亡相关的风险比为 15.96 倍,但较弱的混杂因素不能做到这一点 (VanderWeele 和 Ding,2017) 。

4.6 非零假设的 E 值

在这一点上,我们描述了如何计算 E 值来评估未测量的混杂因素与治疗 (处理) 和结果之间的最小关联强度,从而将点估计值或 CI 的一个极限移到零。然而,我们可以使用类似的程序来评估将估计值移到 RR 其他值所需的两个混杂相关的最小值。在 evalue 中,我们通过在 true() 选项中指定所需的值来实现这一点。

例如,卫生保健研究和质量机构的一项研究 (Ip 等,2007) 中报告母乳喂养与儿童白血病之间的 RR 为 0.80 (95% 置信信区间:[0.71,0.91]) 。通过计算零值效应的 E 值,得到点估计值和 CI 分别为 1.81和 1.43。

. evalue rr 0.80, lcl(0.71) ucl(0.91)
   E-value (point estimate): 1.809
   E-value (CI): 1.429

假设我们有兴趣评估两个未测量的混杂关联需要多大才能将 RR 估计值从 0.80 转移到 0.90。我们只需指定 true(0.90)

. evalue rr 0.80, lcl(0.71) ucl(0.91) true(0.90)
   E-value (point estimate): 1.500
   E-value (CI): 1.000

如结果所示,我们得到了 1.50 点估计值的 E 值,它描述了一个未测量的混杂因素与母乳喂养和儿童白血病之间的关联程度,以将观察到的 RR 从 0.80 移动到 0.90。这个非零 E 值的解释是“对于一个未测量的混杂因素,将观察到的 RR=0.80 的估计值转换为 RRT=0.90 的估计值,这是与母乳喂养和儿童白血病相关的未测混杂因素,风险比为 1.5 倍可以做到这一点,但更弱的混杂却无法做到” (VanderWeele 和 Ding,2017) 。因为 CI 已经包含了 0.90 的值,所以不需要额外的未测量的混杂来包含该值,因此 CI 包含 0.90 的 E 值仅为 E-value = 1.0。

我们还可以在原假设的另一侧为 RR 的值计算 E 值。因此,如果我们要评估将 RR 估计值从 0.80 移至 RR 估计值 1.20 所需的两个混杂关联的最小强度,则只需指定 true(1.20)

. evalue rr 0.80, lcl(0.71) u(0.91) true(1.20)
   E-value (point estimate): 2.366
   E-value (CI): 1.967

如结果所示,为了将估计值转换为 RR 为1.20,我们得到点的 E 值为 2.37,CI 上限的 E 值为 1.97。这些非空 E 值的解释是,对于一个未测量的混杂因素,将观察到的 RR 估计值 0.80 转换为 RR为 1.20,一个未测量的混杂因素与母乳喂养和儿童白血病相关,二者的 RR 为2.37 倍可以做到,但较弱的混杂因素不能做到这一点。同样,为了将置信上限从 0.91 变为 1.20,与母乳喂养和白血病相关的未测量的混杂因素 (RR 分别为 1.97倍) 可以做到这一点,但较弱的混杂因素不能做到这一点 (VanderWeele 和 Ding,2017) 。

5. 结语

与其他方法相比,该方法的主要优势在于,根据常见的处理效果分析,该方法易于使用。研究人员以通常的方式 (例如,使用协变量进行回归) 拟合调整后的模型,然后将 evalue 应用于处理变量及其 CI 的系数。同样, evalue 可以很容易地计算出发表研究中通常报告中的处理效果估计值的 E 值,从而使读者能够评估这种程度的混杂因素关联是否合理。

总之, evalue 工具包,可以方便用于在观察性研究中评估处理效果后进行敏感性分析。我们主张在此类研究中报告 E 值,以协助研究人员和其他研究人员权衡证据以证明混杂的稳健性,并最终得出因果关系。

6. 参考资料

  • Ariel Linden, Maya B. Mathur, Tyler J. VanderWeele. 2020. Conducting sensitivity analysis for unmeasured confounding in observational studies using E-values: The evalue package, Stata Journal, 20(1): 162–175. -PDF-, -PDF2-
  • VanderWeele, T., & Ding, P. .2017. Sensitivity Analysis in Observational Research: Introducing the E-Value. Annals of Internal Medicine, 167, 268-274.-PDF-,-Link-

7. 相关推文

Note:产生如下推文列表的 Stata 命令为:
lianxh 敏感性 因果推断
安装最新版 lianxh 命令:
ssc install lianxh, replace

相关课程

免费公开课

最新课程-直播课

专题 嘉宾 直播/回看视频
最新专题 文本分析、机器学习、效率专题、生存分析等
研究设计 连玉君 我的特斯拉-实证研究设计-幻灯片-
面板模型 连玉君 动态面板模型-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课,2小时]
  • Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。

课程主页

课程主页

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 连享会-主页知乎专栏,400+ 推文,实证分析不再抓狂。直播间 有很多视频课程,可以随时观看。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会-常见问题解答:
https://gitee.com/lianxh/Course/wikis

New! lianxh 命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh