R实现类似NEJM格式的逻辑回归亚组分析森林图技术求助

你可以使用forestplot包对亚组分析结果进行自定义绘图，完全匹配NEJM期刊的森林图样式，完整可运行代码如下：

# 安装依赖包（未安装时运行）
# install.packages(c("data.table", "Publish", "survival", "forestplot", "dplyr"))

# 加载包
library(data.table)
library(Publish)
library(survival)
library(forestplot)
library(dplyr)

# 数据预处理和亚组分析（与你的原有代码一致）
data(traceR)
setDT(traceR)
traceR[,':='(wmi2=factor(wallMotionIndex<0.9,levels=c(TRUE,FALSE), 
                         labels=c("bad","good")),
             abd2=factor(abdominalCircumference<95, levels=c(TRUE,FALSE), 
                         labels=c("slim","fat")))]
traceR[,sex:=as.factor(sex)]
traceR[observationTime==0,observationTime:=1]
traceR=na.omit(traceR)

# 拟合回归模型+亚组分析
fit_log <- glm(dead~treatment+age,family="binomial",data=traceR)
sub_log <- subgroupAnalysis(fit_log,traceR,treatment="treatment",
                            subgroups=~smoking+sex+wmi2+abd2, factor.reference="inline")

# ---------------------- 定制NEJM风格森林图 ----------------------
# 整理绘图数据
plot_data <- sub_log %>%
  mutate(
    # 拼接亚组标签
    label = ifelse(level == "", as.character(var), paste0("  ", level)),
    # 拼接效应值文本
    or_text = sprintf("%.2f (%.2f-%.2f)", OR, lower, upper),
    p_text = ifelse(is.na(p), "", sprintf("%.3f", p))
  ) %>%
  select(label, or_text, p_text, OR, lower, upper)

# 加上表头行
header <- tibble(
  label = "亚组",
  or_text = "OR (95% CI)",
  p_text = "P值",
  OR = NA,
  lower = NA,
  upper = NA
)
plot_data <- bind_rows(header, plot_data)

# 绘制森林图
forestplot(
  labeltext = plot_data[, c("label", "or_text", "p_text")],
  mean = plot_data$OR,
  lower = plot_data$lower,
  upper = plot_data$upper,
  # 对齐NEJM样式设置
  new_page = TRUE,
  is.summary = c(TRUE, rep(FALSE, nrow(plot_data)-1)),
  xlog = FALSE,
  col = fpColors(
    box = "black",
    line = "black",
    zero = "gray50"
  ),
  xlab = "风险比 (OR)",
  zero = 1,
  lwd.zero = 1,
  lwd.ci = 1,
  boxsize = 0.3,
  graph.pos = 2,
  txt_gp = fpTxtGp(
    label = gpar(fontfamily = "sans", cex = 0.9),
    xlab = gpar(fontfamily = "sans", cex = 1),
    ticks = gpar(fontfamily = "sans", cex = 0.8)
  ),
  ci.vertices = TRUE,
  ci.vertices.height = 0.1
)

你可以根据需求调整以下参数匹配个人排版需要：

• boxsize：森林图中点的大小
• cex相关参数：各类文字的大小
• ci.vertices：是否显示置信区间两端的小竖线
• xlab：横轴标题文本

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