scipy.stats.

pearsonr#

scipy.stats.pearsonr(x, y, *, alternative='two-sided', method=None, axis=0)[源码]#

皮尔逊相关系数和用于检验不相关性的 p 值。

皮尔逊相关系数 [1] 衡量两个数据集之间的线性关系。与其他相关系数一样,该系数在 -1 和 +1 之间变化,0 表示没有相关性。-1 或 +1 的相关性表示精确的线性关系。正相关表示随着 x 的增加,y 也随之增加。负相关表示随着 x 的增加,y 随之减少。

此函数还执行一项检验,即检验样本所依据的分布不相关且呈正态分布的零假设。(有关输入数据非正态性对相关系数分布的影响的讨论,请参阅 Kowalski [3]。)p 值粗略地表示一个不相关的系统产生至少与从这些数据集中计算出的皮尔逊相关性一样极端的 Pearson 相关性的数据集的概率。

参数:
xarray_like

输入数组。

yarray_like

输入数组。

axisint 或 None,默认值

进行计算的轴。默认值为 0。如果为 None,则在执行计算之前展平两个数组。

版本 1.14.0 中新增。

alternative{‘two-sided’, ‘greater’, ‘less’}, optional

定义备择假设。默认为“双侧”。以下选项可用:

  • ‘two-sided’:相关性不为零

  • ‘less’:相关性为负(小于零)

  • ‘greater’:相关性为正(大于零)

版本 1.9.0 中新增。

methodResamplingMethod, optional

定义用于计算 p 值的方法。如果 methodPermutationMethod/MonteCarloMethod 的实例,则使用 scipy.stats.permutation_test/scipy.stats.monte_carlo_test 和提供的配置选项以及其他适当的设置来计算 p 值。否则,p 值将按说明中所述进行计算。

1.11.0 版新增。

返回:
resultPearsonRResult

一个具有以下属性的对象

statisticfloat

皮尔逊积矩相关系数。

pvaluefloat

与所选备择假设相关的 p 值。

该对象具有以下方法

confidence_interval(confidence_level, method)

这会计算给定置信水平下相关系数 statistic 的置信区间。置信区间以具有 lowhigh 字段的 namedtuple 返回。如果未提供 method,则使用 Fisher 变换 [1] 计算置信区间。如果 methodBootstrapMethod 的实例,则使用 scipy.stats.bootstrap 和提供的配置选项以及其他适当的设置来计算置信区间。在某些情况下,由于退化的重采样,置信限可能为 NaN,这对于非常小的样本(约 6 个观测值)来说是典型的。

引发:
ValueError

如果 xy 的长度小于 2。

警告:
ConstantInputWarning

当输入是常量数组时引发。在这种情况下,相关系数未定义,因此返回 np.nan

NearConstantInputWarning

当输入“接近”恒定时引发。如果 norm(x - mean(x)) < 1e-13 * abs(mean(x)),则认为数组 x 接近恒定。在这种情况下,计算 x - mean(x) 过程中的数值误差可能会导致 r 计算不准确。

另请参阅

spearmanr

斯皮尔曼秩相关系数。

kendalltau

肯德尔的 tau,一种用于序数数据的相关性度量。

皮尔逊相关性

扩展示例

附注

相关系数计算如下:

\[r = \frac{\sum (x - m_x) (y - m_y)} {\sqrt{\sum (x - m_x)^2 \sum (y - m_y)^2}}\]

其中 \(m_x\) 是向量 x 的均值,\(m_y\) 是向量 y 的均值。

假设 x 和 y 来自独立的正态分布(因此总体相关系数为 0),样本相关系数 r 的概率密度函数为([1], [2]

\[f(r) = \frac{{(1-r^2)}^{n/2-2}}{\mathrm{B}(\frac{1}{2},\frac{n}{2}-1)}\]

其中 n 是样本数量,B 是 beta 函数。这有时被称为 r 的精确分布。这是 pearsonrmethod 参数保留其默认值(None)时用于计算 p 值的分布。该分布是区间 [-1, 1] 上的 beta 分布,其形状参数 a = b = n/2 - 1 相等。在 SciPy 对 beta 分布的实现中,r 的分布是

dist = scipy.stats.beta(n/2 - 1, n/2 - 1, loc=-1, scale=2)

默认情况下,pearsonr 返回的 p 值是双侧 p 值。对于具有相关系数 r 的给定样本,p 值是来自具有零相关性的总体的随机样本 x' 和 y' 的 abs(r') 大于或等于 abs(r) 的概率。在上述 dist 对象方面,给定 r 和长度 n 的 p 值可以计算为

p = 2*dist.cdf(-abs(r))

当 n 为 2 时,上述连续分布未定义。当形状参数 a 和 b 趋近于 a = b = 0 时,可以将 beta 分布的极限解释为具有在 r = 1 和 r = -1 处相等概率质量的离散分布。更直接地说,可以观察到,给定数据 x = [x1, x2] 和 y = [y1, y2],并假设 x1 != x2 和 y1 != y2,r 的唯一可能值是 1 和 -1。由于长度为 2 的任何样本 x' 和 y' 的 abs(r') 都将是 1,因此长度为 2 的样本的双侧 p 值始终为 1。

为了向后兼容,返回的对象也表现得像一个长度为两个的元组,其中包含统计量和 p 值。

数组 API 标准支持

pearsonr 除了 NumPy 之外,还对 Python Array API Standard 兼容后端提供了实验性支持。请通过设置环境变量 SCIPY_ARRAY_API=1 并将 CuPy、PyTorch、JAX 或 Dask 数组作为数组参数提供来测试这些功能。支持以下后端和设备(或其他功能)的组合。

CPU

GPU

NumPy

不适用

CuPy

不适用

PyTorch

JAX

Dask

不适用

有关更多信息,请参阅 对数组 API 标准的支持

参考文献

[1] (1,2,3)

“Pearson correlation coefficient”, Wikipedia, https://en.wikipedia.org/wiki/Pearson_correlation_coefficient

[2]

Student, “Probable error of a correlation coefficient”, Biometrika, Volume 6, Issue 2-3, 1 September 1908, pp. 302-310。

[3]

C. J. Kowalski, “On the Effects of Non-Normality on the Distribution of the Sample Product-Moment Correlation Coefficient” Journal of the Royal Statistical Society. Series C (Applied Statistics), Vol. 21, No. 1 (1972), pp. 1-12。

示例

>>> import numpy as np
>>> from scipy import stats
>>> x, y = [1, 2, 3, 4, 5, 6, 7], [10, 9, 2.5, 6, 4, 3, 2]
>>> res = stats.pearsonr(x, y)
>>> res
PearsonRResult(statistic=-0.828503883588428, pvalue=0.021280260007523286)

执行测试的精确置换版本

>>> rng = np.random.default_rng()
>>> method = stats.PermutationMethod(n_resamples=np.inf, random_state=rng)
>>> stats.pearsonr(x, y, method=method)
PearsonRResult(statistic=-0.828503883588428, pvalue=0.028174603174603175)

在数据来自*均匀*分布的零假设下执行测试

>>> method = stats.MonteCarloMethod(rvs=(rng.uniform, rng.uniform))
>>> stats.pearsonr(x, y, method=method)
PearsonRResult(statistic=-0.828503883588428, pvalue=0.0188)

生成一个 90% 的渐近置信区间

>>> res.confidence_interval(confidence_level=0.9)
ConfidenceInterval(low=-0.9644331982722841, high=-0.3460237473272273)

以及用于 bootstrap 置信区间的

>>> method = stats.BootstrapMethod(method='BCa', rng=rng)
>>> res.confidence_interval(confidence_level=0.9, method=method)
ConfidenceInterval(low=-0.9983163756488651, high=-0.22771001702132443)  # may vary

如果提供了 N 维数组,则根据大多数 scipy.stats 函数的约定,在一次调用中执行多个测试。

>>> rng = np.random.default_rng()
>>> x = rng.standard_normal((8, 15))
>>> y = rng.standard_normal((8, 15))
>>> stats.pearsonr(x, y, axis=0).statistic.shape  # between corresponding columns
(15,)
>>> stats.pearsonr(x, y, axis=1).statistic.shape  # between corresponding rows
(8,)

要执行数组切片之间的所有成对比较,请使用标准的 NumPy 广播技术。例如,要计算行之间所有对的关联

>>> stats.pearsonr(x[:, np.newaxis, :], y, axis=-1).statistic.shape
(8, 8)

如果 y = a + b*x + e,其中 a、b 是常数,e 是随机误差项,假设与 x 独立,则 x 和 y 之间存在线性关系。为简单起见,假设 x 是标准正态分布,a=0,b=1,令 e 服从均值为零、标准差为 s>0 的正态分布。

>>> rng = np.random.default_rng()
>>> s = 0.5
>>> x = stats.norm.rvs(size=500, random_state=rng)
>>> e = stats.norm.rvs(scale=s, size=500, random_state=rng)
>>> y = x + e
>>> stats.pearsonr(x, y).statistic
0.9001942438244763

这应该接近于由以下公式给出的精确值

>>> 1/np.sqrt(1 + s**2)
0.8944271909999159

对于 s=0.5,我们观察到高度相关性。一般来说,噪声的方差越大,相关性越小,而当误差方差趋于零时,相关性趋于一。

重要的是要记住,不相关并不意味着独立,除非 (x, y) 是联合正态的。相关性即使在存在非常简单的依赖关系时也可能为零:如果 X 服从标准正态分布,令 y = abs(x)。请注意,x 和 y 之间的相关性为零。事实上,由于 x 的期望值为零,cov(x, y) = E[x*y]。根据定义,这等于 E[x*abs(x)],由于对称性,它等于零。以下代码行说明了这一观察结果

>>> y = np.abs(x)
>>> stats.pearsonr(x, y)
PearsonRResult(statistic=-0.05444919272687482, pvalue=0.22422294836207743)

非零相关系数可能具有误导性。例如,如果 X 具有标准正态分布,定义 y = x(如果 x < 0),否则 y = 0。简单的计算表明 corr(x, y) = sqrt(2/Pi) = 0.797…,这意味着高度相关

>>> y = np.where(x < 0, x, 0)
>>> stats.pearsonr(x, y)
PearsonRResult(statistic=0.861985781588, pvalue=4.813432002751103e-149)

这是不直观的,因为如果 x 大于零(如果我们对 x 和 y 进行抽样,这种情况大约占一半),则 x 和 y 之间不存在依赖关系。

有关更详细的示例,请参阅 皮尔逊相关性