基于分布函数生成随机数

Question

如果我们想使随机数在一个区间 [a,b] 内均匀生成，那么很容易实现：

A=rand()*(b-a) + a

rand()是一个函数，它可以生成介于0和1之间的均匀随机数。所以A是[a,b]中的随机数。

对于基于分布函数（例如y=x-x^2）生成随机数时，我遇到了问题。

我想使用这里提到的方法。但我不想使用Python函数inverse_cdf(np.random.uniform())。

我可以通过在0和X上进行积分来计算函数“y”的CDF，并将其称为“f”。但是，当我将rand()函数（介于0和1之间的数字）放入f的反函数中时，会得到一个复数！ 也就是说：A=f^(-1) (rand())返回一个复数。

这是基于分布函数生成随机数的正确方法吗？

我使用了这个网站来计算f=x^2/2 - x^3/3的反函数，下面的代码是计算的一部分，显示tmp1始终为负数。

for i=1:10
rnd1=rand;
tmp1  = 2*sqrt(6)*sqrt(6*rnd1^2-rnd1)-12*rnd1+1
cTmp1 = tmp1^(1/3)
end

- Denis

---

我很难理解这个问题，你能否添加一些示例代码并标记使用的编程语言？ - Sam Mason

@AnderBiguri，我在问题中添加了更多信息。 - Denis

具体来说，x-x^2 仅在 0 ≤ x ≤ 1 时为正。在该范围内进行积分可得面积为 1/6。因此，您可以通过指定 f(x) = 6(x - x^2) for 0 ≤ x ≤ 1; 0 otherwise 来使其成为有效密度。但是，您还可以添加一个常数来将方程向上或向下移动，这将影响范围和结果面积。那么您想要的实际密度是什么？ - pjs

@pjs 我计算了从0到z的有效f(x)的积分，得到了CDF为F(z)=3z^2-2z^3。问题在于反演这个函数，再次给我sqrt(X^2-X)。因此，如果我将一个随机变量放在0和1之间，我会得到一个复数。 - Denis

当我添加约束条件时，Wolfram无法给出闭合形式的解答，因此只能通过绘图来解决。由于即使Wolfram也无法得出干净的代数解，因此您应该能够使用牛顿法在数值上解决它。 - pjs

Answer 1

问题是：“基于分布函数生成随机数是否正确？”因此，基本上您要么有连续的概率密度函数（PDF），要么有离散的概率质量函数（PMF），记为f（x），您正在寻找一种找到随机变量x的方法。

至少有两种方法可以做到这一点：

1.使用反演变换分布。 2.使用拒绝方法 使用反演变换： 如果我们知道概率分布函数的函数，则对于某些累积分布函数（CDF），我们可以找到随机变量的闭合形式。假设您的概率函数为f（x），CDF为F（x），则假定您可以获得反函数，则可以获得随机变量。

x=inverse F(U)

其中U是指随机均匀分布。

使用拒绝抽样方法：如果概率分布函数（CDF）没有封闭形式的反函数，那么您可以始终使用拒绝抽样方法。拒绝法的主要思想是生成一个二维随机点（一对随机数）：(r1, r2)，然后该点要么位于概率密度函数(PDF)曲线下方，要么位于其上方。如果该点在曲线下方，则将其作为我们的随机数值；否则，我们将重新采样另一个点。
假设概率密度函数 f(x) 的上界为 M。

r1 is generated within interval [a, b] of horizontal axis
r2 is generated within interval [0, M] of vertical axis
If r2 <f (r1) then accept r1 and output r1
   Else reject r1 and generate new random point

逆变换法比拒绝法更优，如果你可以找到CDF的反函数，因为你可以得到闭合形式。例如，速率为1的指数分布的CDF是F(x) = 1-exp(-x)。逆变换将会是：

x= inverse F(U) = -ln(1-U) = -ln(U)

因为 1-U2=U2

- Kardi Teknomo