整合 PDF 得到错误的累积分布

Question

我有一个 PDF（下面的代码），我正在尝试获取它的 cummulative distribution function (CDF)。我可以使用从 0 整合 PDF 的标准方法。（它是最小允许值）增加 x 值。这按预期工作：

PDF 为红色，"stepwise" CDF 为蓝色。但这给我留下了 table 的 (x, y) 数据，我想要描述此 CDF 的实际函数。所以我转而将 PDF 从 0. 积分到一个变量 y 来获得这个表达式。我这样做 with WolframAlpha (WA) 我发现：

我使用 upper incomplete gamma function 将这个函数写入我的代码，我得到了这个：

其中 WA 集成 CDF 为橙色。我尝试了较低的不完全伽马函数，但结果更糟。

我很确定 WA CDF 的编写没有错误，所以我不确定我在这里做错了什么。

import numpy as np
from scipy.special import gamma, gammaincc
from scipy import integrate
import matplotlib.pyplot as plt

def main():
    xx = np.linspace(0., 5., 100)
    yy = PDF(xx)
    plt.plot(xx, yy, c='r')

    # Find stepwise CDF
    cummul = []
    for x in xx:
        cummul.append([x, integrate.quad(PDF, 0., x)[0]])
    plt.plot(*np.array(cummul).T)

    # WolframAlpha's integral
    y = CDF(xx)
    plt.plot(xx, y)


def PDF(y):
    a = (343. / 15.) * np.sqrt(7. / (2. * np.pi))
    b, c = 5. / 2., -7. / 2.
    return a * (y ** b) * np.exp(c * y)


def CDF(x):
    """
    From WolframAlpha
    """
    a = (343. / 15.) * np.sqrt(7. / (2. * np.pi))
    b, c = 5. / 2., -7. / 2.
    return a * (
        x**b * (-c * x)**(-b) * (gammaincc(b + 1, -c * x) - b * gamma(b))) / c

if __name__ == '__main__':
    main()

Answer 1

scipy.special.gammaincc 是 正则化 上不完全伽马函数。如果不完全 gamma 函数是 gamma(a, x)，则重新排序的不完全 gamma 函数是 gamma(a, x)/gamma(a)。 WA 公式中的不完全伽马函数不是正则化版本。如果将 gammaincc(b + 1, -c * x) 乘以 gamma(b + 1)，您将得到预期的结果。即把return语句改成

    return a * (
        x**b * (-c * x)**(-b) * (gamma(b + 1) * gammaincc(b + 1, -c * x) - b * gamma(b))) / c

也可以这样写

    return a * (
        x**b * (-c * x)**(-b) * gamma(b + 1) * (gammaincc(b + 1, -c * x) - 1)) / c

如果我使用后一个版本，并将 Wolfram Alpha 积分图更改为

    plt.plot(xx, y, '--', linewidth=4, alpha=0.6)

我得到这个情节：

整合 PDF 得到错误的累积分布

Wrong cummulative distribution obtained from integrating the PDF

python

scipy

cdf