斯皮尔曼相关系数#
斯皮尔曼等级顺序相关系数是一种非参数度量，用于衡量两个数据集之间关系的单调性。
考虑以下来自 [1] 的数据，该研究探讨了不健康的人类肝脏中游离脯氨酸（一种氨基酸）和总胶原蛋白（一种常存在于结缔组织中的蛋白质）之间的关系。
下面的 x 和 y 数组记录了这两种化合物的测量值。这些观测值是成对的：每个游离脯氨酸测量值与同一索引处的总胶原蛋白测量值取自同一肝脏。
import numpy as np
# total collagen (mg/g dry weight of liver)
x = np.array([7.1, 7.1, 7.2, 8.3, 9.4, 10.5, 11.4])
# free proline (μ mole/g dry weight of liver)
y = np.array([2.8, 2.9, 2.8, 2.6, 3.5, 4.6, 5.0])
这些数据在 [2] 中使用斯皮尔曼相关系数进行了分析，该统计量对样本之间的单调相关性敏感，并以 scipy.stats.spearmanr 的形式实现。
from scipy import stats
res = stats.spearmanr(x, y)
res.statistic
0.7000000000000001
对于具有强正序数相关的样本，此统计量的值往往较高（接近 1），对于具有强负序数相关的样本，此统计量的值较低（接近 -1），而对于具有弱序数相关的样本，此统计量的值在幅度上较小（接近于零）。
该测试是通过将观测到的统计量值与零分布进行比较来进行的：零分布是在总胶原蛋白和游离脯氨酸测量值相互独立的零假设下导出的统计量值的分布。
对于此测试，可以转换统计量，使得大样本的零分布为自由度为 len(x) - 2 的学生 t 分布。
import matplotlib.pyplot as plt
dof = len(x)-2 # len(x) == len(y)
dist = stats.t(df=dof)
t_vals = np.linspace(-5, 5, 100)
pdf = dist.pdf(t_vals)
fig, ax = plt.subplots(figsize=(8, 5))
def plot(ax): # we'll reuse this
ax.plot(t_vals, pdf)
ax.set_title("Spearman's Rho Test Null Distribution")
ax.set_xlabel("statistic")
ax.set_ylabel("probability density")
plot(ax)
plt.show()
通过 p 值量化比较：p 值是零分布中大于或等于观测到的统计量值的极端值的比例。在统计量为正的双侧检验中，大于转换后的统计量的零分布元素和小于观测到的统计量负数的零分布元素都被认为是“更极端”。
fig, ax = plt.subplots(figsize=(8, 5))
plot(ax)
rs = res.statistic # original statistic
transformed = rs * np.sqrt(dof / ((rs+1.0)*(1.0-rs)))
pvalue = dist.cdf(-transformed) + dist.sf(transformed)
annotation = (f'p-value={pvalue:.4f}\n(shaded area)')
props = dict(facecolor='black', width=1, headwidth=5, headlength=8)
_ = ax.annotate(annotation, (2.7, 0.025), (3, 0.03), arrowprops=props)
i = t_vals >= transformed
ax.fill_between(t_vals[i], y1=0, y2=pdf[i], color='C0')
i = t_vals <= -transformed
ax.fill_between(t_vals[i], y1=0, y2=pdf[i], color='C0')
ax.set_xlim(-5, 5)
ax.set_ylim(0, 0.1)
plt.show()
res.pvalue
0.07991669030889909
如果 p 值“小”——也就是说，如果从独立分布中抽样数据产生如此极端的统计量值的概率很低——这可能被视为反对零假设的证据，支持备择假设：总胶原蛋白和游离脯氨酸的分布不是独立的。 请注意
反之则不然；也就是说，该检验不用来为零假设提供证据。
将被认为是“小”的值的阈值应在分析数据之前做出选择 [3]，并考虑假阳性（错误地拒绝零假设）和假阴性（未能拒绝错误的零假设）的风险。
小的 p 值不是大效应的证据；相反，它们只能为“显著”效应提供证据，这意味着它们不太可能在零假设下发生。
假设在进行实验之前，作者有理由预测总胶原蛋白和游离脯氨酸测量值之间存在正相关关系，并且他们选择针对单侧备择假设来评估零假设的合理性：游离脯氨酸与总胶原蛋白具有正序数相关性。在这种情况下，仅考虑零分布中大于或等于观测统计量的值被认为更极端。
res = stats.spearmanr(x, y, alternative='greater')
res.statistic
0.7000000000000001
fig, ax = plt.subplots(figsize=(8, 5))
plot(ax)
pvalue = dist.sf(transformed)
annotation = (f'p-value={pvalue:.6f}\n(shaded area)')
props = dict(facecolor='black', width=1, headwidth=5, headlength=8)
_ = ax.annotate(annotation, (3, 0.018), (3.5, 0.03), arrowprops=props)
i = t_vals >= transformed
ax.fill_between(t_vals[i], y1=0, y2=pdf[i], color='C0')
ax.set_xlim(1, 5)
ax.set_ylim(0, 0.1)
plt.show()
res.pvalue
0.03995834515444954
请注意，t 分布提供了零分布的渐近近似值；它仅在具有大量观测值的样本中是准确的。对于小样本，执行置换检验可能更合适 [4]：在总胶原蛋白和游离脯氨酸相互独立的零假设下，每个游离脯氨酸测量值都有可能与任何总胶原蛋白测量值一起观测到。 因此，我们可以通过计算 x 和 y 之间元素之间每种可能的配对下的统计量来形成精确的零分布。
def statistic(x): # explore all possible pairings by permuting `x`
rs = stats.spearmanr(x, y).statistic # ignore pvalue
transformed = rs * np.sqrt(dof / ((rs+1.0)*(1.0-rs)))
return transformed
ref = stats.permutation_test((x,), statistic, alternative='greater',
permutation_type='pairings')
fig, ax = plt.subplots(figsize=(8, 5))
plot(ax)
ax.hist(ref.null_distribution, np.linspace(-5, 5, 26),
density=True)
ax.legend(['asymptotic approximation\n(many observations)',
f'exact \n({len(ref.null_distribution)} permutations)'])
plt.show()
ref.pvalue
0.04563492063492063