mgcv包gam cox.ph模型绘图y轴：风险率还是风险比？

问题背景 #

我正在用R的mgcv包拟合带多暴露组的非线性Cox比例风险模型（广义可加模型），代码如下：

library(mgcv)

V <- rep(1, nrow(dt))

fit <- gam(cbind(V, group_number) ~ s(time_elapsed, exposure_group, bs='fs', k=1, m=1) + covs,
  data = dt,
  family=cox.ph,
  weights=dt$outcome,
  control=gam.control(trace=TRUE, maxit=500)
)

plot.gam(fit)

我的患者水平数据集dt包含：

• exposure_group：有序因子（1、2、3、4，1为参考组），代表多暴露组
• outcome：二分类结局变量
• time_elapsed：结局发生天数
• group_number：分层变量
• covs：其他协变量

使用plot.gam()绘图时，x轴是time_elapsed很明确，但默认y轴仅显示s(time_elapsed, exposure_group)及内部数字，文档未明确说明y轴代表绝对风险率还是对数转换后的相对风险比，想了解如何解读这个y轴。

注：此问题并非《Hazard Ratio Plot from mgcv::gam cox.ph model》的重复问题，因我有多暴露组，绘图含多条曲线，但该问题的答案可能适用。

解答 #

• 在mgcv的cox.ph族GAM模型中，plot.gam()输出的y轴代表对数尺度下的风险比（相对于参考组的对数HR），更准确地说，是平滑项s(time_elapsed, exposure_group)对模型线性预测器的贡献值。
• Cox模型的线性预测器本质是对数风险的线性组合，你使用的fs（factor-smooth）平滑项会为每个非参考暴露组生成一条随time_elapsed变化的曲线，曲线数值代表该暴露组在控制协变量covs和分层后，相对于参考组的对数风险差。
• 若要将y轴转换为直观的风险比（HR），可在绘图时添加trans = exp参数，代码示例：

  plot.gam(fit, trans = exp)
  

  此时y轴显示实际风险比，数值>1代表该暴露组风险高于参考组，<1则低于参考组。
• 默认y轴的数值是中心化结果，参考组的曲线会被固定在0附近（平滑项的约束逻辑），其他组曲线偏离0的程度直接对应对数HR的大小。

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