如何在metafor中绘制同研究效应量同色的漏斗图（funnel plot）？

我正尝试在metafor包中绘制漏斗图（funnel plot）以检验发表偏倚，期望同一研究的效应量在图中显示为相同颜色。近期我看到一篇文献采用如下代码实现该需求：

mod1 <- rma(yi, vi, data = dat3)
pdf(NULL)
dev.control(displaylist="enable")
par(mar=c(4,4,0.1,0) + 0.1)
funnel(mod1, col =inferno(25)[pos],
       ylim = c(0, 0.26), xlim = c(-0.6, 0.8),
       ylab = "Standard error (SE)", xlab = "Effect size (lnRR)")

想请教col参数中的inferno(25)[pos]具体含义，以及如何实现这类同一研究效应量同色的漏斗图。

1. inferno(25)[pos]的含义 #

• inferno()是viridis包中的配色函数，生成的是一组色盲友好、高对比度的渐变色系，inferno(25)表示生成包含25种颜色的向量。
• pos是一个与效应量数量等长的数值向量，同一研究的所有效应量对应相同的数值。比如某研究包含3个效应量，这3个效应量对应的pos值完全一致，这样inferno(25)[pos]就会把同一种颜色分配给该研究的所有效应量。

2. 实现同一研究效应量同色的漏斗图步骤 #

步骤1：标记研究分组 #

首先确保你的数据集dat3有一个区分不同研究的变量（比如study_id），将其转换为数值型索引向量pos：

# 假设dat3中存在study_id列，将其转为分组索引
dat3$pos <- as.integer(factor(dat3$study_id))

同一研究的所有行，pos值会保持一致。

步骤2：安装并加载viridis包 #

如果未安装viridis包，先执行安装：

install.packages("viridis")
library(viridis)

除了inferno()，你也可以选择viridis()、magma()等其他同系列配色函数，效果类似。

步骤3：绘制漏斗图 #

在funnel()函数的col参数中使用inferno(n)[pos]，其中n只需大于等于你的研究总数即可：

mod1 <- rma(yi, vi, data = dat3)
# 调整绘图边距
par(mar=c(4,4,0.1,0) + 0.1)
# 绘制漏斗图
funnel(mod1, 
       col = inferno(25)[dat3$pos],  # 用pos匹配同一研究的颜色
       ylim = c(0, 0.26), 
       xlim = c(-0.6, 0.8),
       ylab = "Standard error (SE)", 
       xlab = "Effect size (lnRR)")

可选：添加图例（辅助解读） #

如果需要让读者明确颜色对应的研究，可以添加图例：

# 获取唯一研究ID和对应颜色
unique_studies <- unique(dat3$study_id)
unique_colors <- inferno(25)[unique(dat3$pos)]
# 添加图例，位置可按需调整
legend("topright", 
       legend = unique_studies, 
       fill = unique_colors, 
       cex = 0.7, 
       title = "Study ID")

内容的提问来源于stack exchange，提问作者madina_b