scipy.integrate.quad 在大范围内给出错误的结果答案

【问题标题】：scipy.integrate.quad gives wrong result on large rangesscipy.integrate.quad 在大范围内给出错误的结果
【发布时间】：2015-09-03 23:48:14
【问题描述】：

我正在尝试对两个“半”正态分布的总和进行积分。 scipy.integrate.quad 在我尝试在小范围内积分时工作正常，但在我为大范围积分时返回 0。代码如下：

mu1 = 0
mu2 = 0
std1 = 1
std2 = 1

def integral_fun(x):
    nor1 = 0.5 * ((1 / (np.sqrt(2 * np.pi) * std1)) * (np.e ** ((-(x-mu1) ** 2) / (2 * std1 **2))))
    nor2 = 0.5 * ((1 / (np.sqrt(2 * np.pi) * std2)) * (np.e ** ((-(x-mu2) ** 2) / (2 * std2 **2))))
    return nor1 + nor2


integrate.quad(integral_fun, -5, 5)
Out[54]: (0.9999994266968564, 8.668320228277793e-10)

integrate.quad(integral_fun, -10, 10)
Out[55]: (1.0000000000000002, 8.671029607900576e-10)

integrate.quad(integral_fun, -100000, 100000)
Out[56]: (0.0, 0.0)

为什么会这样？

【问题讨论】：

数值积分通过函数值的加权和来近似。如果我不得不猜测，我会说函数值都是 0，因为该范围内几乎所有点的密度（几乎）为 0。
如果范围是从-inf 到inf，它可以工作
@ThePredator 这是因为在这种情况下，scipy 会将积分转换为(0,1] 范围内的一个整数，转换为x = (1-t)/t。这样就可以很好地对原点附近的点进行采样。

标签： python python-3.x scipy integration quad

【解决方案1】：

这里的原因是，您的函数仅在积分区域的一个非常小的区域中达到非常强烈的峰值，并且在其他任何地方实际上为零，quad 从未找到此峰值，因此只能看到被积函数为零。

由于在这种情况下您知道峰在哪里，因此拆分积分的限制以便单独考虑峰周围的区域是合理的。

为此，您可以使用points 参数，以一种稍微有点混蛋的方式强制quad 分别考虑峰值。

In [3]: integrate.quad(integral_fun, -100000, 100000, points=[-10,10])
Out[3]: (1.0000000000000002, 8.671029607900576e-10)

【讨论】：