R中如何从stanfit对象获取类brmsfit的后验条件效应估计值

brms的conditional_effects()没有特殊逻辑，核心流程是基于参数的全后验抽样，逐点计算目标统计量的后验分布，再做分位数/均值汇总，这套方法适配任意结构的自定义Stan模型，不需要手动推导预测的解析形式。

具体实现步骤 #

• 第一步：构造预测网格
  确定你要展示效应的自变量取值序列，其余协变量固定为典型值（连续变量取均值/中位数、分类变量取众数，和brms默认规则一致）。比如要展示x的效应，就生成x从样本最小值到最大值的等距取值点：

  # 单自变量示例，多变量直接在数据框里补充固定取值的协变量列即可
  pred_grid <- data.frame(
    x = seq(min(original_data$x), max(original_data$x), length.out = 100)
  )
  

• 第二步：提取全量后验抽样
  不要用参数的汇总统计值，要提取每一次MCMC迭代对应的参数抽样值，这是适配任意复杂模型的核心，用rstan自带的extract()函数即可：

  # permuted=TRUE时返回按参数名存储的抽样列表，每个元素是向量/矩阵，长度/行数等于有效迭代次数
  post_samples <- rstan::extract(your_stanfit_object)
  

• 第三步：逐点计算目标量的后验分布
  完全复现你Stan模型中目标统计量的计算逻辑即可：如果要展示期望响应的效应，就照搬Stan模型中线性预测器mu的计算代码；如果要展示考虑观测噪声的预测区间，就照搬生成观测值y的分布抽样逻辑。
  以简单线性回归（模型中mu = alpha + beta * x，观测噪声为sigma）为例：

  n_iter <- length(post_samples$alpha) # 有效后验迭代次数
  n_point <- nrow(pred_grid) # 预测网格点数
  # 初始化矩阵：行对应迭代次数，列对应预测点，存储每个点在每次迭代下的计算值
  post_pred <- matrix(NA, nrow = n_iter, ncol = n_point)
  
  for (i in 1:n_point) {
    # 这里替换成你自己Stan模型里的目标量计算逻辑即可，不管多复杂都直接照搬
    post_pred[, i] <- post_samples$alpha + post_samples$beta * pred_grid$x[i]
    
    # 如果要算带观测噪声的预测区间，把上面一行替换成下面的形式即可
    # post_pred[, i] <- rnorm(n_iter, 
    #                        mean = post_samples$alpha + post_samples$beta * pred_grid$x[i],
    #                        sd = post_samples$sigma)
  }
  

  不管模型是非线性结构、含多层随机效应、加了特殊参数化形式，只要和你Stan代码里的计算逻辑保持一致，结果就不会出错。

• 第四步：逐点汇总得到标准格式结果
  对每个预测点对应的后验向量，计算对应统计量，和brms输出的列名完全对齐：

  result <- pred_grid
  result$estimate__ <- colMeans(post_pred) # 后验均值
  result$se__ <- apply(post_pred, 2, sd) # 后验标准误
  result$lower__ <- apply(post_pred, 2, quantile, probs = 0.025) # 95%可信区间下界
  result$upper__ <- apply(post_pred, 2, quantile, probs = 0.975) # 95%可信区间上界
  

  得到的result可以直接用于绘制效应曲线，格式和conditional_effects()返回的对应效应表完全一致。如果需要调整可信区间宽度，修改quantile()里的概率值即可。

注意事项 #

• 构造多变量模型的预测网格时，非焦点协变量不要取随机值，固定为样本典型值才能得到和brms一致的调整后条件效应
• 如果需要计算连续变量的瞬时边际效应，只需要把第三步中计算的目标量替换为“x取x0和x0+极小值时的预测值差/极小值”即可，其余流程不变
• 后验抽样计算时尽量避免用参数的点估计（比如后验均值）代入计算，必须用全量迭代的抽样值，否则得到的可信区间会忽略参数的协方差结构，结果存在偏差

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