我正在尝试用Python(scipy)替换Python脚本中的一些rpy2代码。在这种情况下,我需要用(Python:kruskal.test()
)替换Kruskal-Wallis测试(R:scipy.stats.kruskal
)。
仅比较整数/浮点数时,scipy.stats.kruskal
返回类似的H统计量和P值。但是,在应用以字符串表示的组时遇到一些困难。
以下是数据的子样本:
y = [4.33917022422, 2.96541899883, 6.70475220836, 9.19889096119, 2.14087398016,
5.39520023918, 1.58443224287, 3.59625224078, 4.01998599966, 2.58058624352]
x = ['High_O2', 'High_O2', 'High_O2', 'High_O2', 'Low_O2',
'Low_O2', 'Low_O2', 'Low_O2', 'Mid_O2', 'Mid_O2']
在R中,只需键入:
kruskal.test(y,as.factor(x))
使用scipy(0.17)在Python(2.7)中做同样的事情:
from scipy import stats
stats.kruskal(y,x)
但是,使用scipy时,我会得到非常低的p值
(p<e-07)
和相当高的H统计量(26)
,这是错误的。我尝试将x
列表替换为{0,1,2}
没有任何改善。如何在排名过程中告诉scipy将
x
视为组? 最佳答案
传递给scipy.stats.kruskal
的每个非关键字参数都被视为一组单独的y值。通过将x
作为参数之一,kruskal
尝试将标签字符串视为第二组y值。字符串将被强制转换为NaN(应该引起RuntimeWarning
)。
相反,您需要按标签对y
值进行分组,然后将它们作为单独的输入数组传递给kruskal
。例如:
# convert `y` to a numpy array for more convenient indexing
y = np.array(y)
# find unique group labels and their corresponding indices
label, idx = np.unique(x, return_inverse=True)
# make a list of arrays containing the y-values corresponding to each unique label
groups = [y[idx == i] for i, l in enumerate(label)]
# use `*` to unpack the list as a sequence of arguments to `stats.kruskal`
H, p = stats.kruskal(*groups)
print(H, p)
# 2.94545454545 0.22929927