使用oddsratio包plot_gam函数可视化GAM模型时遇错求助

问题原因 #

oddsratio::plot_gam函数对mgcv中张量样条项（ti()）的处理逻辑存在缺陷，当模型包含ti(age, los)这类交互张量项时，函数提取主效应（如age）的预测数据时，会将其识别为原子向量而非列表，导致$x操作符调用失败。

解决方案 #

方案1：手动构造预测数据并绘图（最稳妥）

绕过plot_gam函数，直接用predict.gam生成预测值，再用ggplot可视化：

# 加载依赖包
library(mgcv)
library(ggplot2)

# 1. 固定其他协变量（这里将los设为中位数）
fixed_los <- median(dat$los, na.rm = TRUE)

# 2. 生成age的序列数据（覆盖数据中age的取值范围）
new_dat <- data.frame(
  age = seq(min(dat$age, na.rm = TRUE), max(dat$age, na.rm = TRUE), length.out = 100),
  los = fixed_los
)

# 3. 从GAM模型中预测link值及标准误
preds <- predict(fit, newdata = new_dat, type = "link", se.fit = TRUE)

# 4. 将link值转换为优势比（OR）及置信区间
new_dat$or <- exp(preds$fit)
new_dat$or_lower <- exp(preds$fit - 1.96 * preds$se.fit)
new_dat$or_upper <- exp(preds$fit + 1.96 * preds$se.fit)

# 5. 绘制主效应图
ggplot(new_dat, aes(x = age, y = or)) +
  geom_line(linewidth = 1) +
  geom_ribbon(aes(ymin = or_lower, ymax = or_upper), alpha = 0.2) +
  labs(x = "年龄", y = "优势比（OR）") +
  theme_minimal()

方案2：修正plot_gam函数源码

如果必须使用plot_gam函数，可修改其依赖的get_gam_plotdata函数中提取张量样条项数据的逻辑：

1. 查看原函数源码：

getAnywhere(get_gam_plotdata)

2. 复制源码并修改关键部分：
   找到处理预测变量的代码块，当识别到ti()类型的平滑项时，从模型的smooth列表中直接提取x值，替换原有的plot_df[[set_pred]]$x调用。例如：

# 修改后的get_gam_plotdata函数（仅展示关键修改部分）
get_gam_plotdata <- function(model, pred, int = NULL, ci = 0.95, ...) {
  # 原有代码...
  
  # 新增：判断是否为ti项
  is_ti <- any(grepl("^ti\\(", names(model$smooth)))
  if (is_ti) {
    # 找到对应pred的ti项
    smooth_idx <- which(sapply(model$smooth, function(s) grepl(paste0("ti\\(", pred, "\\)"), s$label)))
    if (length(smooth_idx) > 0) {
      plot_df[[pred]] <- list(x = model$smooth[[smooth_idx]]$x,
                              y = model$smooth[[smooth_idx]]$y)
    }
  }
  
  # 原有代码...
}

3. 重新定义修改后的函数，再调用plot_gam即可。

方案3：升级oddsratio包

先尝试升级到最新版包，看官方是否已修复该问题：

install.packages("oddsratio")

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