我有一堆数据点和一个模型，现在我正在寻找可以排除模型的置信水平。
我从一个卡方值开始，伴随着一个自由度数的值。
使用 scipy，我一直在尝试的看起来像这样：
from scipy.stats import chi2
from scipy.stats import norm
chisq = 74.1
df = 21 #degrees of freedom
cdf = chi2.cdf(chisq, df,scale=1)
sigma = norm.interval(cdf)
这会产生我正在寻找的输出（5.5 sigma 置信度），但对于更高的 chi^2 值会失败。
在我看来，问题在于 scipy 使用的数据类型的精度。对于 chi^2 值 334.7（和相同数量的自由度），chi2.cdf 产生 1.0，这反过来使 norm.interval 返回 -inf/inf。这应该最好返回 16.1 sigma 置信度。
我找到了一种使用 mpmath 为 cdf 获得更好精度的方法：
import mpmath
mpmath.mp.dps = 200 # decimal digits of precision
def cdf(x,k):
x,k = mpmath.mpf(x), mpmath.mpf(k)
return mpmath.gammainc(k/2, 0, x/2, regularized=True)
但是，问题仍然存在，因为 scipy.stats.norm.interval 似乎将输入四舍五入为 1.0。
有没有办法通过更改数据类型 norm.interval 使用来规避这种情况，或者是否有其他方法/包来计算具有任意精度输入的正态分布的端点？