scipy.stats.pearsonr 的最小 p 值

Question

我在我的数据上运行 scipy.stats.pearsonr，我得到

(0.9672434106763087, 0.0)

r 值高而 p 值非常低是合理的。但是，p显然不是0，所以我想知道p=0.0是什么意思。是 p<10^-10，p<10^-100 还是极限？

Answer 1

正如@MB-F 在评论中指出的那样，它是通过分析计算的。

在版本0.19.1的代码中，您可以隔离该部分代码并根据以下形式绘制 p 值r

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import betainc
r = np.linspace(-1, 1, 1000)*(1-1e-10);

for n in [10, 100, 1000]:
    df = n - 2
    t_squared = r**2 * (df / ((1.0 - r) * (1.0 + r)))
    prob = betainc(0.5*df, 0.5, df/(df+t_squared))
    plt.semilogy(r, prob, label=f'n={n}')
plt.axvline(0.9672434106763087, ls='--', color='black', label='r value')
plt.legend()
plt.grid()

当前稳定版本1.9.3使用不同的公式

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import btdtr
r = np.linspace(-1, 1, 1000)*(1-1e-10);
for n in [10, 100, 1000]:
    ab = 0.5*n
    prob = btdtr(ab, ab, 0.5*(1-abs(r)))
    plt.semilogy(r, prob, label=f'n={n}')
plt.axvline(0.9672434106763087, ls='--', color='black', label='r value')
plt.legend()
plt.grid()

但产生相同的结果。

您可以看到，如果您有 1000 个点和您的相关性，则该p值将小于最小浮动值。

贝塔分布

Scipy 提供了概率分布的集合，其中包括beta 分布。

线路

    prob = btdtr(ab, ab, 0.5*(1-abs(r)))

可以替换为

from scipy.stats import beta
prob = beta(ab, ab).cdf(0.5*(1-abs(r)))

在那里您可以获得更多有关它的信息。