python：从自定义概率函数中随机抽样

Question

我有一个带概率密度函数的分段四次分布：

p(x)= c(x/a)^2 if 0?x<a; 
      c((b+a-x)^2/b)^2 if a?x?b;
      0 otherwise

假设 c、a、b 是已知的，我试图从分布中抽取 100 个随机样本。我怎样才能用 numpy/scipy 做到这一点？

Answer 1

一种标准方法是为G = F^-1累积分布函数的逆找到一个明确的公式。这在这里是可行的（尽管它自然会分段定义），然后使用[0,1] 上的G(U)where Uis uniform 来生成您的样本。

在这种情况下，我认为我已经弄清楚了细节，但是您需要检查微积分/代数。

首先，为了简化事情，引入几个新参数是有帮助的。让

f(a,b,c,d,x) = c*x**2 #if 0 <= x <= a

和

f(a,b,c,d,x) = d*(x-e)**4 #if a < x <= b

那么你p(x)的由

p(x) = f(a,b,c/a**2,c/b**2,a+b)

我整合f以找到累积分布，然后反转并得到以下结果：

def Finverse(a,b,c,d,e,x):
    if x <= (c*a**3)/3:
        return (3*x/c)**(1/3)
    else:
        return e + ((a-e)**5 - (5*c*a**3)/(3*d))**(1/5)

假设这是正确的，那么简单地：

def randX(a,b,c):
    u = random.random()
    return Finverse(a,b,c/a**2,c/b**2,a+b,u)

在这种情况下，可以计算出一个明确的公式。当您无法计算出这样的逆公式时，请考虑使用@lucianopaz 描述的蒙特卡罗方法