拒绝抽样
您可以使用rejection sampling。下面的方法通过在比原始空间小 1 维的空间中进行采样来实现这一点。
- 第 1 步：通过从均匀分布中抽样每个 x(i) 来随机抽样 x(1)、x(2)、...、x(n-1)
- 步骤 2：如果总和 S = x(1) + x(2) + ... + x(n-1) 小于 0 或大于 2，则拒绝并在步骤 1 中重新开始。
- 第 3 步：计算第 n 个变量为 x(n) = 1-S
直觉
您可以在笛卡尔坐标为 ±1, ±1 的 n 维立方体内部查看向量 x(1), x(2), ..., x(n-1), x(n) ，...，±1。这样您就可以遵循约束 -1
坐标总和必须等于 1 的附加约束将坐标约束到比超立方体更小的空间，并且将是维度为 n-1 的hyperplane。
如果您进行常规拒绝采样，从均匀分布中采样所有坐标，那么您将永远不会命中约束。采样点永远不会在超平面中。因此，您考虑 n-1 个坐标的子空间。现在您可以使用拒绝抽样。
视觉上
假设您有维度 4，那么您可以从 4 中绘制坐标中的 3。该图（均匀地）填充了一个多面体。下面通过在切片中绘制多面体来说明这一点。每个切片对应不同的总和 S = x(1) + x(2) + ... + x(n-1) 和不同的 x(n) 值。
图像：3 个坐标的域。每个彩色表面与第 4 个坐标的不同值相关。
边际分布
对于大维度，拒绝抽样将变得不那么有效，因为拒绝的比例随着维度的数量而增加。
“解决”这个问题的一种方法是从边缘分布中抽样。但是，计算这些边际分布有点乏味。比较：对于从狄利克雷分布生成样本，存在similar algorithm，但在这种情况下，边缘分布相对容易。 （然而，推导出这些分布并非不可能，见下文“与 Irwin Hall 分布的关系”）
在上面的示例中，x(4) 坐标的边缘分布对应于切口的表面积。因此，对于 4 维，您可能能够根据该图计算出计算（您需要计算那些不规则多边形的面积），但对于更大的维度，它开始变得更加复杂。
与 Irwin Hall 分布的关系
要获得边缘分布，您可以使用截断的Irwin Hall distributions。欧文霍尔分布是均匀分布变量之和的分布，将遵循一些分段多项式形状。下面以一个示例进行演示。
代码
由于我的 python 生锈了，我将主要添加 R 代码。该算法非常基本，因此我想任何 Python 编码器都可以轻松地将其改编成 Python 代码。在我看来，这个问题的难点在于算法而不是如何在 Python 中编码（虽然我不是 Python 编码器，所以我把它留给其他人）。
图像：采样输出。 4 条黑色曲线是四个坐标的边缘分布。红色曲线是基于 Irwin Hall 分布的计算。这可以通过直接计算而不是拒绝抽样来扩展到抽样方法。
python中的拒绝采样
import numpy as np
def sampler(size):
reject = 1
while reject:
x = np.random.rand(size - 1) # step 1
S = np.sum(x)
reject = (S<0) or (S>2) # step 2
x = np.append(x,1-S) # step 3
return[x]
y = sampler(5)
print(y, np.sum(y))
R 中的更多代码，包括与 Irwin Hall 分布的比较。此分布可用于计算边际分布，并可用于设计一种比拒绝抽样更有效的算法。
### function to do rejection sample
samp <- function(n) {
S <- -1
## a while loop that performs step 1 (sample) and 2 (compare sum)
while((S<0) || (S>2) ) {
x <- runif(n-1,-1,1)
S <- sum(x)
}
x <- c(x,1-S) ## step 3 (generate n-th coordinate)
x
}
### compute 10^5 samples
y <- replicate(10^5,samp(4))
### plot histograms
h1 <- hist(y[1,], breaks = seq(-1,1,0.05))
h2 <- hist(y[2,], breaks = seq(-1,1,0.05))
h3 <- hist(y[3,], breaks = seq(-1,1,0.05))
h4 <- hist(y[4,], breaks = seq(-1,1,0.05))
### histograms together in a line plot
plot(h1$mids,h1$density, type = 'l', ylim = c(0,1),
xlab = "x[i]", ylab = "frequency", main = "marginal distributions")
lines(h2$mids,h2$density)
lines(h3$mids,h3$density)
lines(h4$mids,h4$density)
### add distribution based on Irwin Hall distribution
### Irwin Hall PDF
dih <- function(x,n=3) {
k <- 0:(floor(x))
terms <- (-1)^k * choose(n,k) *(x-k)^(n-1)
sum(terms)/prod(1:(n-1))
}
dih <- Vectorize(dih)
### Irwin Hall CDF
pih <- function(x,n=3) {
k <- 0:(floor(x))
terms <- (-1)^k * choose(n,k) *(x-k)^n
sum(terms)/prod(1:(n))
}
pih <- Vectorize(pih)
### adding the line
### (note we need to scale the variable for the Erwin Hall distribution)
xn <- seq(-1,1,0.001)
range <- c(-1,1)
cum <- pih(1.5+(1-range)/2,3)
scale <- 0.5/(cum[1]-cum[2]) ### renormalize
### (the factor 0.5 is due to the scale difference)
lines(xn,scale*dih(1.5+(1-xn)/2,3),col = 2)