拒绝采样详解（包含例子）

为什么需要拒绝采样？

当很难直接从已知的密度函数 $f(x)$ 采样时，我们就可以使用拒绝采样，通过从建议分布 $g(x)$ 生成样本并根据一定的接受-拒绝准则来保留或拒绝这些样本，从而得到符合 $f(x)$ 的样本。在特别是高维或复杂的分布，一个简单的、我们知道如何采样的分布（如均匀分布或标准正态分布）可能更容易处理。拒绝采样为我们提供了一种桥梁，使用易于采样的分布来近似我们真正感兴趣的复杂分布。

原理解释

假设你有一个单位圆，半径为1，中心在原点。你的目标是从这个圆内均匀地采样点。但你只知道如何从其包围的正方形（边长为2，从-1到1）中均匀地采样点。

拒绝采样的方法如下：

1. 均匀撒点：随机选择一个点，其x和y坐标都是从-1到1的均匀分布。
2. 接受或拒绝：检查这个点是否在单位圆内。最简单的方法是计算点到原点的距离，如果距离小于或等于1，那么这个点就在圆内；否则，它在正方形内，但不在圆内。
3. 重复：继续此过程，直到你获得了所需数量的在圆内的点。

如果你尝试这个过程，你会发现大约有四分之一的点被拒绝（因为正方形的面积是4，而圆的面积是π，所以接受率大约是π/4）。但这完全取决于你如何选择你的“建议分布”和你的“目标分布”——在这个例子中，是正方形和圆。

具体步骤

1. 选择一个建议分布 $g(x)$
2. 找到一个常数 $c$ ，使得对于所有 $t$ ， $f(t)>0$ ，满足 $f(t)/g(t)\le c$ 。
3. 重复以下步骤，直到得到所需数量的样本：
   从建议分布 $g(x)$ 生成随机值 $y$ 。
   从 Uniform(0, 1) 分布生成一个随机 $u$ 。
   如果 $u<\frac{f(y)}{cg(y)}$ ，接受 $y$ ，并代入 $x=y$ ；否则拒绝 $y$

我们可以用以下公式描述给定一个样本 $Y$ 的情况下，它被接受的概率。

$P(accept|Y)=P(U<\frac{f(Y)}{cg(Y)}|Y)=\frac{f(Y)}{cg(Y)}$

$P(U<\frac{f(Y)}{cg(Y)}|Y)$ 这表示，给定样本 $Y$ 的情况下，随机数 $U$ 小于 $\frac{f(Y)}{cg(Y)}$ 的概率。
$\frac{f(Y)}{cg(Y)}$ ：这是接受 $Y$ 的实际概率。直观地说，这比较了我们真正感兴趣的目标分布 $f(Y)$ 和调整后的提议分布 $c\times g(Y)$ 之间的关系。如果 $f(Y)$ 在某个点 $Y$ 的值很高，而 $g(Y)$ 的值相对较低，那么这个比率就会接近1，这意味着在这个点上几乎总是接受样本。

例子

使用拒绝采样从 $Beta(\alpha =2,\beta =2)$ 分布中采样的R代码：

  reject_sample <- function(n) {
    samples <- numeric(0)
    
    while (length(samples) < n) {
      y <- runif(1, 0, 1)  # 建议分布
      u <- runif(1, 0, 1) # 决定是否接受这个样本
      
      px <- dbeta(y, 2, 2)   # Density of Beta(2, 2)
      qx <- 1  # Density of U(0, 1)
      C <- 2
      
      if (u <= px / (C * qx)) {
        samples <- c(samples, y)
      }
    }
    
    return(samples)
  }
  
  # Generate 1000 samples from Beta(2, 2)
  samples <- reject_sample(1000)

绘图：

  # 绘制拒绝采样生成的样本的直方图
  hist(samples, probability = TRUE, breaks = 30, col = rgb(0.2, 0.8, 0.5, 0.5), main = "Samples from Beta(2, 2) using Rejection Sampling", xlab = "x", ylab = "Density")
  
  # 添加理论的 Beta(2, 2) 密度曲线
  curve(dbeta(x, 2, 2), add = TRUE, col = "red", lwd = 2)
  
  # 添加图例
  legend("topright", legend = c("Sampled", "Theoretical"), fill = c(rgb(0.2, 0.8, 0.5, 0.5), "red"))