使用贝叶斯预测概率来进行临床试验的中期分析-wutao-blog

贝叶斯预测概率（Bayesian Predictive Probabilities）进行临床试验的中期分析（Interim analysis）

如何计算成功的预测概率（Predictive Probability of Success，PPoS）？

PPoS 的定义：根据当前的中期数据，在未来的分析中获得一个成功的结果的概率

比如先验分布是 Beta 分布，数据是二项分布，得到后验 Beta 分布，基于这个后验分布，未来再收集 N 个样本，其中至少 K 个成功的概率是多少

可以通过模拟的方法来计算这个预测概率：

在中期分析中，从基于现在数据得到的后验分布中抽样感兴趣的参数（比如从 Beta-Binomial 模型的后验分布中抽样）
利用第一步抽样得到的参数构成的预测分布中，对未来的样本进行抽样，得到未观测的成功事件数
将第二步得到的数据加入到现在的中期分析数据中，对数据进行补全，利用这个新数据重新计算后验分布，如果后验概率满足标准，则这个模拟实验成功
重复上面三步 B 次，其中成功的比例则为 PPoS

先看一个例子：现在有一个 100 个病人的单臂二期临床研究，检测的是二分类变量（是否对某个治疗有效）；目标是将实际反应率和金标准反应率 50% 做比较，看其是否有效。实际有效的病人数量 $x$ 服从二项分布：

x \sim Bin(p,N=100)

其中 $p$ 是在研究的群体中的反应率，N 是病人的数量

如果反应率大于金标准 0.5 的后验概率大于 0.95，那么就认为这个试验是成功的，也就是：

Pr(p>0.5|x)>\eta = 0.95

先验概率分布适用 Beta 分布：

p \sim Beta(\alpha=1,\beta=1)

所以这个反应率的后验分布为（依据 Beta-Binomial 共轭模型）：

p|(X=x)\sim Beta(\alpha+x,\beta + N -x)=Beta(1+x,1+100-x)

因此要大于 0.95，则需要在 100 个病人中至少 59 个病人有效：

1-pbeta(0.5,58+1,1+100-58)
[1] 0.9445408

1-pbeta(0.5,59+1,1+100-59)
[1] 0.9636215

然后我们假设分别在入组 20, 50 和 75 个病人时进行中期分析（interim analysis）：

在第一个中期分析中，20 个病人观测到 12 个病人有效，p 大于 0.5 的概率为 0.81：
```
1-pbeta(0.5,12+1,1+20-12)
[1] 0.8083448
```

那么如果后续继续入组 80 个病人，这个临床试验成功的概率是多少呢？这里面就需要计算上面提到的 PPoS 了，基于上面描述的步骤写个函数来进行模拟计算：

get_p <- function(alpha, beta, interim_N, interim_x, future_N, sample_counts, threshold, p_threshold){
  ##第一步, 从后验分布中抽样 N 个参数值
  p <- rbeta(sample_counts, interim_x + alpha, beta + interim_N - interim_x)
  ##第二步，对每个参数值构成的二项分布，从未来的数据中抽样得到相应的成功次数
  tt <- sapply(p, function(x){
    success <- rbinom(1, future_N, x)
    ##第三步，将抽样得到的数据加入现在的数据中，重新计算符合标准的后验概率
    p_success <- 1 - pbeta(threshold, 
                           alpha + success + interim_x, 
                           beta + interim_N + future_N - (success + interim_x))
    return(ifelse(p_success > p_threshold,1,0)) 
  })
  return(mean(tt == 1)) ###第四步 N次模拟中有多少比例是成功的
}

对于第一个中期分析：

##算个 10 次，看看均值
dt <- sapply(c(1:10),function(x){get_p(alpha = 1, beta = 1, 
                                       interim_N = 20, interim_x = 12,
                                       future_N = 80, sample_counts = 1000, 
                                       threshold = 0.5,
                                       p_threshold = 0.95)})
mean(dt)
0.5446

也就是说未来再入组 80 个病人，这个试验的成功概率为 0.54。

在第二次中期分析中，病人的数量增多了，变成了 50 个病人，其中 28 个有效果，p 大于 0.5 的概率为 0.799：

1-pbeta(0.5,28+1,1+50-28)
[1] 0.799469

未来再入组 50 个病人，这个试验的成功概率为 0.31：

dt <- sapply(c(1:10),function(x){get_p(alpha = 1, beta = 1, 
                                       interim_N = 50, interim_x = 28,
                                       future_N = 50, sample_counts = 1000, 
                                       threshold = 0.5,
                                       p_threshold = 0.95)})
mean(dt)
[1] 0.3053

在第三次中期分析时，病人的数量为 75，其中 41 个有效，p 大于 0.5 的概率为 0.79：

 1-pbeta(0.5,41+1,1+75-41)
[1] 0.7889044

未来再入组 25 个病人，这个试验的成功概率只有 0.09 了：

dt <- sapply(c(1:10),function(x){get_p(alpha = 1, beta = 1, 
                                       interim_N = 75, interim_x = 41,
                                       future_N = 25, sample_counts = 1000, 
                                       threshold = 0.5,
                                       p_threshold = 0.95)})
mean(dt)
[1] 0.0859

参考：The Utility of Bayesian Predictive Probabilities for Interim Monitoring of Clinical Trials

目录CONTENT

使用贝叶斯预测概率来进行临床试验的中期分析

评论区