Meta系列（二）： 在R中如何计算效应值
1 简介
在上一篇文章中，我们介绍了meta分析中有哪些效应值，以及什么是随机效应模型和混合效应模型。下面我们就继续来看看这些效应值的函数在R中如何计算和运用。
2 效应值计算
2.1 预先计算效应值（以标准化均值差为例）
示例数据来源于demetar包。对照组与试验组效应大小的标准化均值差已经事先计算好。
library(tiderverse) # need for 'glimpse'
library(demetar)
library(meta)
data(ThirdWave)
glimpse(ThirdWave)
## Rows: 18
## Columns: 8
## $ Author <chr> "Call et al.", "Cavanagh et al.", "DanitzOrsillo"…
## $ TE <dbl> 0.7091362, 0.3548641, 1.7911700, 0.1824552, 0.421…
## $ seTE <dbl> 0.2608202, 0.1963624, 0.3455692, 0.1177874, 0.144…
## $ RiskOfBias <chr> "high", "low", "high", "low", "low", "low", "high…
## $ TypeControlGroup <chr> "WLC", "WLC", "WLC", "no intervention", "informat…
## $ InterventionDuration <chr> "short", "short", "short", "short", "short", "sho…
## $ InterventionType <chr> "mindfulness", "mindfulness", "ACT", "mindfulness…
## $ ModeOfDelivery <chr> "group", "online", "group", "group", "online", "g…
可以看到数据集包含的数据列，主要的效应指标是TE，效应指标的标准误是seTE。由数据集中可以看出有不同的亚组，可以猜测这些研究间存在异质性，所以采用随机效应模型。选择极大似然估计，并且使用Knapp-Hartung法降低假阳性率，也就是tau方估计选用REML。
使用函数：meta包中的metagen
必须的参数：TE:数据集中预计算的效应值大小
seTE：数据集中预计算的效应值标准误
m.gen <- metagen(TE = TE,
seTE = seTE,
studlab = Author,#设置研究的标签
data = ThirdWave,
sm = "SMD",#设置估计的效应值
fixed = FALSE,#是否估计固定效应模型
random = TRUE,#是否估计随机效应模型
method.tau = "REML",#tau方的估计方法
hakn = TRUE,
title = "Third Wave Psychotherapies")
summary(m.gen)#查看拟合结果
## Review: Third Wave Psychotherapies
## SMD 95%-CI %W(random)
## Call et al. 0.7091 [ 0.1979; 1.2203] 5.0
## Cavanagh et al. 0.3549 [-0.0300; 0.7397] 6.3
## DanitzOrsillo 1.7912 [ 1.1139; 2.4685] 3.8
## de Vibe et al. 0.1825 [-0.0484; 0.4133] 7.9
## Frazier et al. 0.4219 [ 0.1380; 0.7057] 7.3
## Frogeli et al. 0.6300 [ 0.2458; 1.0142] 6.3
## Gallego et al. 0.7249 [ 0.2846; 1.1652] 5.7
## Hazlett-Steve… 0.5287 [ 0.1162; 0.9412] 6.0
## Hintz et al. 0.2840 [-0.0453; 0.6133] 6.9
## Kang et al. 1.2751 [ 0.6142; 1.9360] 3.9
## Kuhlmann et al. 0.1036 [-0.2781; 0.4853] 6.3
## Lever Taylor… 0.3884 [-0.0639; 0.8407] 5.6
## Phang et al. 0.5407 [ 0.0619; 1.0196] 5.3
## Rasanen et al. 0.4262 [-0.0794; 0.9317] 5.1
## Ratanasiripong 0.5154 [-0.1731; 1.2039] 3.7
## Shapiro et al. 1.4797 [ 0.8618; 2.0977] 4.2
## Song & Lindquist 0.6126 [ 0.1683; 1.0569] 5.7
## Warnecke et al. 0.6000 [ 0.1120; 1.0880] 5.2
##
## Number of studies combined: k = 18
##
## SMD 95%-CI t p-value
## Random effects model 0.5771 [0.3782; 0.7760] 6.12 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.0820 [0.0295; 0.3533]; tau = 0.2863 [0.1717; 0.5944];
## I^2 = 62.6% [37.9%; 77.5%]; H = 1.64 [1.27; 2.11]
##
## Test of heterogeneity:
## Q d.f. p-value
## 45.50 17 0.0002
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
结果已经输出，如何解读结果是现在的重要议题。
结果包含什么
第一部分：单个研究的效应值和置信区间。由于效应值已经事先计算，这里没有包含太多信息。W%是每个研究的权重。de Vibe et al的研究给的权重最高（7.9%），而Ratanasiripong 给的权重最低，查看其效应值置信区间，置信区间大，意味着标准误非常大，因此，提示此项研究的效应大小并不精确。此外，数据结果还告诉我们一共纳入了18项研究。
第二部分：这部分提供了核心结果：合并效应大小。由于我们使用了Knapp-Hartung法估计参数，所以效应值估计使用T
分布。
接下来是异质性。现在我们暂时先只关注tau方，我们可以看到tau方为0.08且置信区间不包含0，因此，研究存在异质性。
第三部分：我们可以看到meta分析的详细信息，比如我们如何估计tau方。
2.1.1 结果调整
有时，我们可能想调整一些分析的细节，updata.meta函数可能会有所帮助，他需要输入meta对象，以及我们需要修啊该的参数。比如，我们想知道采用别的tau方估计方法是否会产生实质性差异，则做出如下修改：
m.gen_update <- update.meta(m.gen,
method.tau = "PM")
# Get pooled effect
m.gen_update$TE.random #0.587
# Get tau^2 estimate
m.gen_update$tau2#0.11
可以看到，效应值大小没有太大改变，但是PM法估计了一个更大的tau方。
2.2 标准化均值差
我们可以使用两组的均值和标准差形式的原始数据进行合并，计算总的效应值。我们必须为metacont函数提供以下7个参数：
-
sm = "SMD"sm = “MD” -
n.e 实验组样本量
-
mean.e 实验组均值
-
sd.e 实验组的均值标准差
-
n.c 对照组样本量
-
mean.c 对照组均值
-
sd.c 对照组的均值标准差
-
method.smd 这个参数的设置是选用，有几个选项“Hedges”(常用，默认)，“Cohen”（未校正的标准化均值差）和“Glass"（默认与推荐，该方法采用对照组的标准差而不是合并后的标准差计算，当有多个实验组时，该方法是初步的分析，但通常不是meta分析的首选方法）
这个例子的数据来源于dmetar，数据主要研究目的是评估干预措施对自杀意念是否有效。
# Make sure meta and dmetar are already loaded
library(meta)
library(dmetar)
library(meta)
# Load dataset from dmetar (or download and open manually)
data(SuicidePrevention)
# Use metcont to pool results.
m.cont <- metacont(n.e = n.e,
mean.e = mean.e,
sd.e = sd.e,
n.c = n.c,
mean.c = mean.c,
sd.c = sd.c,
studlab = author,
data = SuicidePrevention,
sm = "SMD",
method.smd = "Hedges",
fixed = FALSE,
random = TRUE,
method.tau = "REML",
hakn = TRUE,
title = "Suicide Prevention")
看看结果是什么：
summary(m.cont)
## Review: Suicide Prevention
##
## SMD 95%-CI %W(random)
## Berry et al. -0.1428 [-0.4315; 0.1459] 15.6
## DeVries et al. -0.6077 [-0.9402; -0.2752] 12.3
## Fleming et al. -0.1112 [-0.6177; 0.3953] 5.7
## Hunt & Burke -0.1270 [-0.4725; 0.2185] 11.5
## McCarthy et al. -0.3925 [-0.7884; 0.0034] 9.0
## Meijer et al. -0.2676 [-0.5331; -0.0021] 17.9
## Rivera et al. 0.0124 [-0.3454; 0.3703] 10.8
## Watkins et al. -0.2448 [-0.6848; 0.1952] 7.4
## Zaytsev et al. -0.1265 [-0.5062; 0.2533] 9.7
##
## Number of studies combined: k = 9
## Number of observations: o = 1147
##
## SMD 95%-CI t p-value
## Random effects model -0.2304 [-0.3734; -0.0874] -3.71 0.0059
##
## Quantifying heterogeneity:
## tau^2 = 0.0044 [0.0000; 0.0924]; tau = 0.0661 [0.0000; 0.3040]
## I^2 = 7.4% [0.0%; 67.4%]; H = 1.04 [1.00; 1.75]
##
## Test of heterogeneity:
## Q d.f. p-value
## 8.64 8 0.3738
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
## - Hedges' g (bias corrected standardised mean difference)
尽管使用的函数不一样，但是输出的结果形式一致。我们可以看到，合并的效应值为-0.23，且有统计学意义。效应大小具有负号，说明与对照组相比，治疗组的自杀意念降低了。我们看到tau方不显著大于0.
2.3 二分类变量
2.3.1 风险和优势比
在R中，我们使用metabin（meta包）函数计算效应值，需要提供以下参数：
- event.e 实验组发生的事件数
- n.e 实验组样本量
- event.c 对照组发生事件数
- n.c 对照组样本量
- method 选择合并效应值的方法，可以是通用的 inverse-variance pooling或者是Mantel-Haenszel（默认和推荐), 或Peto method, 或 Bakbergenuly-sample size method，分别对应的调用参数是：“Inverse”，“MH”，`“Peto"和"SSW”
- sm 需要计算的效应值：RR或OR
- incr 为校正存在0单元格而设置的参数，通常建议是不用设置，且后续设置不使用连续性校正。
- MH.exact 是否需要连续性校正，TRUE表示不使用
本例使用的数据集来源于dmatar包，是关于抑郁症与死亡率的数据集。
library(dmetar)
library(tidyverse)
library(meta)
data(DepressionMortality)
glimpse(DepressionMortality)
## Rows: 18
## Columns: 6
## $ author <chr> "Aaroma et al., 1994", "Black et al., 1998", "Bruce et al., 19…
## $ event.e <dbl> 25, 65, 5, 26, 32, 1, 24, 15, 15, 173, 37, 41, 29, 61, 15, 21,…
## $ n.e <dbl> 215, 588, 46, 67, 407, 44, 60, 61, 29, 1015, 105, 120, 258, 38…
## $ event.c <dbl> 171, 120, 107, 1168, 269, 87, 200, 437, 227, 250, 66, 9, 24, 3…
## $ n.c <dbl> 3088, 1901, 2479, 3493, 6256, 1520, 882, 2603, 853, 3375, 409,…
## $ country <chr> "Finland", "USA", "USA", "USA", "Sweden", "USA", "Canada", "Ne…
在本例中，我们将RR作为效应值，采用随机效应模型，由于我们正在计算二分类变量的结构，我们使用Paule-Mandel法估计tau方。
通过查看数据，我们发现，样本量因研究而异，所以，采用PM法可能有些许偏差，我们后续可以尝试其他的估计方法作为敏感性分析。
数据中不包含0单元格，所以不需要连续性校正。
m.bin <- metabin(event.e = event.e,
n.e = n.e,
event.c = event.c,
n.c = n.c,
studlab = author,
data = DepressionMortality,
sm = "RR",
method = "MH",
MH.exact = TRUE,
fixed = FALSE,
random = TRUE,
method.tau = "PM",
hakn = TRUE,
title = "Depression and Mortality")
summary(m.bin)
## Review: Depression and Mortality
## RR 95%-CI %W(random)
## Aaroma et al., 1994 2.09 [1.41; 3.12] 6.0
## Black et al., 1998 1.75 [1.31; 2.33] 6.6
## Bruce et al., 1989 2.51 [1.07; 5.88] 3.7
## Bruce et al., 1994 1.16 [0.85; 1.57] 6.5
## Enzell et al., 1984 1.82 [1.28; 2.60] 6.3
## Fredman et al., 1989 0.39 [0.05; 2.78] 1.2
## Murphy et al., 1987 1.76 [1.26; 2.46] 6.4
## Penninx et al., 1999 1.46 [0.93; 2.29] 5.8
## Pulska et al., 1998 1.94 [1.34; 2.81] 6.2
## Roberts et al., 1990 2.30 [1.92; 2.75] 7.0
## Saz et al., 1999 2.18 [1.55; 3.07] 6.3
## Sharma et al., 1998 2.05 [1.07; 3.91] 4.7
## Takeida et al., 1997 6.97 [4.13; 11.79] 5.3
## Takeida et al., 1999 5.81 [3.88; 8.70] 6.0
## Thomas et al., 1992 1.33 [0.77; 2.27] 5.3
## Thomas et al., 1992 1.77 [1.10; 2.83] 5.6
## Weissman et al., 1986 1.25 [0.66; 2.33] 4.8
## Zheng et al., 1997 1.98 [1.40; 2.80] 6.3
##
## Number of studies combined: k = 18
##
## RR 95%-CI t p-value
## Random effects model 2.0217 [1.5786; 2.5892] 6.00 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.1865 [0.0739; 0.5568]; tau = 0.4319 [0.2718; 0.7462];
## I^2 = 77.2% [64.3%; 85.4%]; H = 2.09 [1.67; 2.62]
##
## Test of heterogeneity:
## Q d.f. p-value
## 74.49 17 < 0.0001
##
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Paule-Mandel estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
由结果可知，RR=2.02（P<0.001),表明抑郁会使死亡风险增加一倍。tau方=0.19，置信区间不包含0，所以研究存在较大的异质性。由于前面提过，PM法可能对tau方估计会产生较大的影响，所以，这次使用REML法。
m.bin_update <- update.meta(m.bin,
method.tau = "REML")
exp(m.bin_update$TE.random)#2.02365
m.bin_update$tau2#0.16473
虽然结果有些许差异，但是还没到推翻结论的地步。如果我们需要汇报OR，预期重写函数，不如直接在原函数上更改参数。
m.bin_or <- update.meta(m.bin, sm = "OR")
m.bin_or
## Review: Depression and Mortality
##
## [...]
##
## Number of studies combined: k = 18
##
## OR 95%-CI t p-value
## Random effects model 2.2901 [1.7512; 2.9949] 6.52 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.2032 [0.0744; 0.6314]; tau = 0.4508 [0.2728; 0.7946];
## I^2 = 72.9% [56.7%; 83.0%]; H = 1.92 [1.52; 2.43]
##
## Test of heterogeneity:
## Q d.f. p-value
## 62.73 17 < 0.0001
##
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Paule-Mandel estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
2.3.2 使用预先计算的数据计算二进制效应大小
有时，我们无法直接提取计算OR或者RR的原始数据，此时，我们采用预先计算的数据计算二进制效应大小（metagen）。
**注意：**如果除使用预先计算的数据之外我们别无选择，数据处理时需要十分注意，该函数只能使用inverse-variance method去计算效应值，无法使用别的方法。但是，在毫无他法的情况下，这种方法仍然是不错的选择。
现在模拟预计算效应大小，将上面计算的效应值提取出来：
DepressionMortality$TE <- m.bin$TE
DepressionMortality$seTE <- m.bin$seTE
现在，假设我们有一个效应，我们知道效应值和置信区间的上下限，但是不知道标准误差，这其实在研究者十分常见，我们需要在数据集增加两列表示上下限，具体模拟如下：
# Set seTE of study 7 to NA
DepressionMortality$seTE[7] <- NA
# Create empty columns 'lower' and 'upper'
DepressionMortality[,"lower"] <- NA
DepressionMortality[,"upper"] <- NA
# Fill in values for 'lower' and 'upper' in study 7
# As always, binary effect sizes need to be log-transformed
DepressionMortality$lower[7] <- log(1.26)
DepressionMortality$upper[7] <- log(2.46)
DepressionMortality[,c("author", "TE", "seTE", "lower", "upper")]#查看数据集
## author TE seTE lower upper
## 1 Aaroma et al., 1994 0.7418 0.20217 NA NA
## 2 Black et al., 1998 0.5603 0.14659 NA NA
## 3 Bruce et al., 1989 0.9235 0.43266 NA NA
## 4 Bruce et al., 1994 0.1488 0.15526 NA NA
## 5 Enzell et al., 1984 0.6035 0.17986 NA NA
## 6 Fredman et al., 1989 -0.9236 0.99403 NA NA
## 7 Murphy et al., 1987 0.5675 NA 0.2311 0.9001
## 8 Penninx et al., 1999 0.3816 0.22842 NA NA
## [...]
尽管我们不知道数据的标准误，但是函数只需提供上下限的列名即可计算效应值，调用如下
m.gen_bin <- metagen(TE = TE,
seTE = seTE,
lower = lower,
upper = upper,
studlab = author,
data = DepressionMortality,
sm = "RR",
method.tau = "PM",
fixed = FALSE,
random = TRUE,
title = "Depression Mortality (Pre-calculated)")
summary(m.gen_bin)
## Review: Depression Mortality (Pre-calculated)
##
## [...]
##
## Number of studies combined: k = 18
##
## RR 95%-CI z p-value
## Random effects model 2.0218 [1.6066; 2.5442] 6.00 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.1865 [0.0739; 0.5568]; tau = 0.4319 [0.2718; 0.7462];
## I^2 = 77.2% [64.3%; 85.4%]; H = 2.09 [1.67; 2.62]
##
## [...]
研究结果与之前的借本相同
2.3.3 发病率比
发病率比，使用函数metainc，metabin，需要输入以下参数：
- event.e: 治疗组事件数
- time.e ： 治疗组人时
- event.c ： 对照组事件数
- time.c： 对照组人时
- method: 默认使用的方法是Mantel and Haenszel（MH），我们也可以使用“Inverse”
- sm: 可以选择是发病率比(IRR)还是发病率差（IRD)
- incr: 需要为0单元格存在而添加的连续性校正
默认情况下，我们不执行连续性校正，所以一般不设置参数，只有在选用逆方差法的时候，我们执行连续性校正。
现在我们使用dmetar包中的数据进行演示：
library(dmetar)
library(tidyverse)
library(meta)
data(EatingDisorderPrevention)
glimpse(EatingDisorderPrevention)# 先查看数据
## Rows: 5
## Columns: 5
## $ Author <chr> "Stice et al., 2013", "Stice et al., 2017a", "Stice et al., 20…
## $ event.e <dbl> 6, 22, 6, 8, 22
## $ time.e <dbl> 362, 235, 394, 224, 160
## $ event.c <dbl> 16, 8, 9, 13, 29
## $ time.c <dbl> 356, 74, 215, 221, 159
我们选用发病率比作为效应值，Mantel-Haenszel法合并效应值，Paule-Mandel法估计tau。
m.inc <- metainc(event.e = event.e,
time.e = time.e,
event.c = event.c,
time.c = time.c,
studlab = Author,
data = EatingDisorderPrevention,
sm = "IRR",
method = "MH",
fixed = FALSE,
random = TRUE,
method.tau = "PM",
hakn = TRUE,
title = "Eating Disorder Prevention")
summary(m.inc)
## Review: Eating Disorder Prevention
##
## IRR 95%-CI %W(random)
## Stice et al., 2013 0.3688 [0.1443; 0.9424] 13.9
## Stice et al., 2017a 0.8660 [0.3855; 1.9450] 18.7
## Stice et al., 2017b 0.3638 [0.1295; 1.0221] 11.5
## Taylor et al., 2006 0.6071 [0.2516; 1.4648] 15.8
## Taylor et al., 2016 0.7539 [0.4332; 1.3121] 40.0
##
## Number of studies combined: k = 5
## Number of events: e = 139
##
## IRR 95%-CI t p-value
## Random effects model 0.6223 [0.3955; 0.9791] -2.91 0.0439
##
## Quantifying heterogeneity:
## tau^2 = 0 [0.0000; 1.1300]; tau = 0 [0.0000; 1.0630]
## I^2 = 0.0% [0.0%; 79.2%]; H = 1.00 [1.00; 2.19]
##
## Test of heterogeneity:
## Q d.f. p-value
## 3.34 4 0.5033
##
## Details on meta-analytical method:
## - Mantel-Haenszel method
## - Paule-Mandel estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
我们可以看到，IRR=0.62，差别有统计学意义。根据合并的效应值，我们有理由人为预防性的干预措施在一年内将饮食失调的发病率降低了38%。我们看到tau方为0。
2.4 相关性
相关性采用逆方差法合并效应值，我们采用metacor函数合并效应值。我们采用metacor函数合并相关性，需要设置以下参数:
- cor: 相关系数
- n 样本量
示例数据来源于dmetar包。
library(dmetar)
library(tidyverse)
library(meta)
data(HealthWellbeing)
glimpse(HealthWellbeing)
## Rows: 29
## Columns: 5
## $ author <chr> "An, 2008", "Angner, 2013", "Barger, 2009", "Doherty, 2013"…
## $ cor <dbl> 0.620, 0.372, 0.290, 0.333, 0.730, 0.405, 0.292, 0.388, 0.3…
## $ n <dbl> 121, 383, 350000, 1764, 42331, 112, 899, 870, 70, 67, 246, …
## $ population <chr> "general population", "chronic condition", "general populat…
## $ country <chr> "South Korea", "USA", "USA", "Ireland", "Poland", "Australi…
我们预计存在较大的异质性，所以采用随机效应模型。
m.cor <- metacor(cor = cor,
n = n,
studlab = author,
data = HealthWellbeing,
fixed = FALSE,
random = TRUE,
method.tau = "REML",
hakn = TRUE,
title = "Health and Wellbeing")
summary(m.cor)
## Review: Health and Wellbeing
## COR 95%-CI %W(random)
## An, 2008 0.6200 [0.4964; 0.7189] 2.8
## Angner, 2013 0.3720 [0.2823; 0.4552] 3.4
## Barger, 2009 0.2900 [0.2870; 0.2930] 3.8
## Doherty, 2013 0.3330 [0.2908; 0.3739] 3.7
## Dubrovina, 2012 0.7300 [0.7255; 0.7344] 3.8
## Fisher, 2010 0.4050 [0.2373; 0.5493] 2.8
## [...]
##
## Number of studies combined: k = 29
##
## COR 95%-CI t p-value
## Random effects model 0.3632 [0.3092; 0.4148] 12.81 < 0.0001
##
## Quantifying heterogeneity:
## tau^2 = 0.0241 [0.0141; 0.0436]; tau = 0.1554 [0.1186; 0.2088];
## I^2 = 99.8% [99.8%; 99.8%]; H = 24.14 [23.29; 25.03]
##
## Test of heterogeneity:
## Q d.f. p-value
## 16320.87 28 0
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
## - Fisher's z transformation of correlations
从上述结果中，我们可以看到合并的相关系数为0.36，可以认定为是一个中等大小的相关性。
2.5 均值
对于均值的合并，仍然使用对数变换，使用时，我们需要事先确定使用原始均值还是对数转换的均值。建议在处理非负变量（如高度）或者是均值接近0时使用对数变换。如果我们需要使用对数变换，需要将sm参数设置为“MRAW"或”MLN“。
我们采用metamean函数估计均值，需要设定以下参数：
- n：观察的样本量
- mean：均值
- sd：均值的标准差
- sm：见上文
示例数据来源于dmetar，是关于抑郁量表得分的汇总。
library(dmetar)
library(tidyverse)
library(meta)
data(BdiScores)
# We only need the first four columns
glimpse(BdiScores[,1:4])
## Rows: 6
## Columns: 4
## $ author <chr> "DeRubeis, 2005", "Dimidjian, 2006", "Dozois, 2009", "Lesperanc…
## $ n <dbl> 180, 145, 48, 142, 301, 104
## $ mean <dbl> 32.6, 31.9, 28.6, 30.3, 31.9, 29.8
## $ sd <dbl> 9.4, 7.4, 9.9, 9.1, 9.2, 8.6
m.mean <- metamean(n = n,
mean = mean,
sd = sd,
studlab = author,
data = BdiScores,
sm = "MRAW",
fixed = FALSE,
random = TRUE,
method.tau = "REML",
hakn = TRUE,
title = "BDI-II Scores")
summary(m.mean)
## Review: BDI-II Scores
##
## mean 95%-CI %W(random)
## DeRubeis, 2005 32.6000 [31.2268; 33.9732] 18.0
## Dimidjian, 2006 31.9000 [30.6955; 33.1045] 19.4
## Dozois, 2009 28.6000 [25.7993; 31.4007] 9.1
## Lesperance, 2007 30.3000 [28.8033; 31.7967] 17.0
## McBride, 2007 31.9000 [30.8607; 32.9393] 20.7
## Quilty, 2014 29.8000 [28.1472; 31.4528] 15.8
##
## Number of studies combined: k = 6
## Number of observations: o = 920
##
## mean 95%-CI
## Random effects model 31.1221 [29.6656; 32.5786]
##
## Quantifying heterogeneity:
## tau^2 = 1.0937 [0.0603; 12.9913]; tau = 1.0458 [0.2456; 3.6043]
## I^2 = 64.3% [13.8%; 85.2%]; H = 1.67 [1.08; 2.60]
##
## Test of heterogeneity:
## Q d.f. p-value
## 14.00 5 0.0156
##
## Details on meta-analytical method:
## - Inverse variance method
## - Restricted maximum-likelihood estimator for tau^2
## - Q-profile method for confidence interval of tau^2 and tau
## - Hartung-Knapp adjustment for random effects model
## - Untransformed (raw) means
2.6 比例
metaprop函数用于合并比例。我们建议在合并比例时采用对数变换（设置参数sm= “PLOGIT”），如果需要采用原始数据也可（sm= “PRAW’),但是一遍不建议这么做。
如果采用对数变换的比例，我们采用广义线性混合模型（GLMM)计算合并的效应值。线性混合效应模型是固定效应模型和随机效应模型同时存在。GLMM不仅可以用于比例的合并，对于二分类的或是计数数据的效应值也可以计算（如比值比或发病率比）。
使用GLMM有三方面的注意事项：1）每个研究的权重不会展示；2）只需要设置tau方的估计方法；3）tau方不会输出置信区间，如果有需要则采用逆方差法。
需要提供以下参数的设定：
- event：观察数
- n：样本数
- method：method = “GLMM”``method = “Inverse”
- incr：是否需要连续性校正，逆方差法时需要
- sm：合并效应值的计算方法，默认是对数转换后的数据（m = “PLOGIT”）
示例数据来源于dmetar包，是一个关于青少年某药物滥用的流行情况的数据。
library(dmetar)
library(meta)
library(tidyverse)
data(OpioidMisuse)
glimpse(OpioidMisuse)
## Rows: 15
## Columns: 3
## $ author <chr> "Becker, 2008", "Boyd, 2009", "Boyd, 2007", "Cerda, 2014", "Fie…
## $ event <dbl> 2186, 91, 126, 543, 6496, 10850, 86, 668, 843, 647, 11521, 1111…
## $ n <dbl> 21826, 912, 1084, 7646, 55215, 114783, 527, 9403, 11274, 8888, …
我们使用GLMM和对数变换合并患病率数据
m.prop <- metaprop(event = event,
n = n,
studlab = author,
data = OpioidMisuse,
method = "GLMM",
sm = "PLOGIT",
fixed = FALSE,
random = TRUE,
hakn = TRUE,
title = "Opioid Misuse")
summary(m.prop)
## Review: Opioid Misuse
## proportion 95%-CI
## Becker, 2008 0.1002 [0.0962; 0.1042]
## Boyd, 2009 0.0998 [0.0811; 0.1211]
## Boyd, 2007 0.1162 [0.0978; 0.1368]
## Cerda, 2014 0.0710 [0.0654; 0.0770]
## Fiellin, 2013 0.1176 [0.1150; 0.1204]
## [...]
##
##
## Number of studies combined: k = 15
## Number of observations: o = 434385
## Number of events: e = 41364
##
## proportion 95%-CI
## Random effects model 0.0944 [0.0836; 0.1066]
##
## Quantifying heterogeneity:
## tau^2 = 0.0558; tau = 0.2362; I^2 = 98.3% [97.9%; 98.7%]; H = 7.74 [6.92; 8.66]
##
## Test of heterogeneity:
## Q d.f. p-value Test
## 838.21 14 < 0.0001 Wald-type
## 826.87 14 < 0.0001 Likelihood-Ratio
##
## Details on meta-analytical method:
## - Random intercept logistic regression model
## - Maximum-likelihood estimator for tau^2
## - Hartung-Knapp adjustment for random effects model
## - Logit transformation
## - Clopper-Pearson confidence interval for individual studies
我们看到青少年某药物滥用的患病率是9.4%，如前文所述，输出不显示每个研究的权重，我们也得到了异质性tau方的估计量，但没有置信区间。
3 总结
- 真实效应大小的方差τ^2，也称为研究间异质性方差，必须在随机效应的meta分析中估计。有几种方法可以做到这一点，哪一种效果最好取决于研究者的实际情况
- 计算合并效应大小的最常见方法是通过逆方差方法。然而，对于二元结局数据，其他方法（如Mantel-Haenszel方法）可能更可取
- 在meta包中，有一个函数可以对预先计算的效应值进行meta分析，以及一套可用于不同类型的“原始”数据的函数
4401
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
