为何scipy.multivariate_normal积分结果与蒙特卡洛模拟概率不符？

你代码里的数值积分部分存在核心参数错误：在定义多元正态分布的协方差矩阵时，误用了标准差（sigmaX、sigmaY）而非方差。

scipy.stats.multivariate_normal的cov参数要求传入方差-协方差矩阵，而你传入的是标准差数值。对于独立二元正态分布，协方差矩阵的对角线元素应为各自变量的方差（即标准差的平方），而非标准差本身。

修正后的数值积分代码 #

def pdf(x,y):
    # 协方差矩阵使用方差：sigmaX² 和 sigmaY²
    return multivariate_normal.pdf([x,y], mean=[0,0], cov=[[sigmaX**2, 0], [0, sigmaY**2]])

验证逻辑 #

由于两个变量独立，联合概率等于边缘概率的乘积：

• X~N(0, 0.5²)，P(|X|<0.5) = 2*Φ(0.5/0.5)-1 ≈ 0.6827
• Y~N(0, 0.8²)，P(|Y|<0.5) = 2*Φ(0.5/0.8)-1 ≈ 0.4648
• 联合概率≈0.6827*0.4648≈0.317，和蒙特卡洛模拟的31.9%（模拟误差范围内）完全匹配，修正后的数值积分结果也会与此一致。

内容的提问来源于stack exchange，提问作者feetwet