scipy 自制的pearson相关性实现在传递两组相同的数据时返回0.999...2

pftdvrlh  于 2023-01-26  发布在  其他
关注(0)|答案(2)|浏览(145)

我厌倦了scipy和numpy,决定继续工作,并在另一个实现的基础上,在某处的SO答案。

from statistics import pstdev, mean

def pearson(x, y):
    sx = []
    sy = []

    mx = mean(x)
    my = mean(y)

    stdx = pstdev(x)
    stdy = pstdev(y)

    for i in x:
        sx.append((i - mx) / stdx)

    for j in y:
        sy.append((j - my) / stdy)

    return sum([i * j for i, j in zip(sx, sy)]) / len(x)

我向它传递了一些数字,看看它是否给出了与scipy.stats.pearsonr相同的结果,它似乎很好。接近末尾的一个数字左右是不同的,但没有什么突破性的东西...
直到我尝试传递相同的数据集给它作为xy,当我这样做,我得到返回0.9999999999999992,当scipy和numpy都说它是1.0
这个实现有什么问题吗?我用的是population stdev而不是sample,据我所知,numpy和scipy都用了这个。我想知道为什么没有返回1.0。可能是python本身的浮点数问题吗?我在Py 2和Py 3中都试过了,都得到了0.999...
如果需要,我传递给它的数据集是:
[7, 1, 5, 1, 8, 5, 9, 8, 5, 10, 5, 8, 1, 8, 8, 8, 10, 4, 8, 9, 9, 6, 8, 7, 8, 5, 10, 5, 6, 8, 8, 7, 9, 4, 6, 10, 7, 10, 4, 5, 4, 7, 4, 8, 9, 10, 9, 8, 7, 8, 6, 8, 6, 6, 5, 7, 7, 7, 7, 3, 7, 8, 6, 8, 5, 7, 8, 7, 8, 6, 8, 6, 9, 6, 6, 6, 8, 9, 5, 7, 9, 2, 9, 6, 7, 6, 7, 7, 5, 5, 7, 7, 8, 6, 9, 1, 3, 6, 7, 9, 7, 7, 6, 9, 9, 4, 9, 9, 7, 9, 6, 2, 2, 8, 4, 7, 7, 6, 3, 7, 3, 5, 10, 9, 8, 10, 8, 7, 4, 7, 8, 9, 8, 4, 7, 9, 7, 7, 6, 8, 8, 9, 9, 7, 4, 4, 7, 3, 9, 3, 1, 8, 3, 9, 4, 8, 3, 9, 8, 8, 7, 9, 9, 8, 10, 8, 3, 10, 4, 7, 7, 10, 8, 7, 8, 7, 1, 8, 9, 5, 7, 5, 5, 3, 5, 7, 7, 7, 2, 4, 1, 6, 9, 9, 7, 7, 10, 9, 2, 9, 8, 2, 5, 1, 2, 5, 9, 1, 4, 8, 9, 6, 4, 4, 7, 3, 7, 9, 4, 3, 7, 8, 7, 6, 8, 8, 7]

sbtkgmzw

sbtkgmzw1#

你对浮点行为的期望太过乐观了。根据经验,你不会对结果不精确为1.0感到惊讶。例如,试试下面这个小得多的输入:

[7, 1, 5]

在我的盒子上,你的函数返回1.0000000000000002。“接近”1.0,但不完全是1.0。这是你能期望的最好结果。
要了解原因,请思考一下这个“应该”计算什么:

math.sqrt(x)**2 == x

数学上(在无限精度下),它应该总是返回True,但在浮点数中(无论使用多少精度,只要精度是有界的),它不可能总是True,事实上,反例很容易找到;就像刚才在我盒子上写的

>>> math.sqrt(2)**2
2.0000000000000004

问题是,* 具有有限精度 *,sqrt()必然是多对一函数。它将域1..N压缩到范围1..sqrt(N)中,并且具有有限精度的域的基数大于范围的基数。因此,在Map到范围中的相同值的域中必然存在不同的xy。因此不存在精确函数逆。
您的函数比普通的sqrt更复杂,但工作原理相同。

62lalag4

62lalag42#

是的,这与浮点行为有关。你可以尝试使用decimal模块

from decimal import *
data = [7, 1, 5, 1, 8, 5, 9, 8, 5, 10, 5, 8, 1, 8, 8, 8, 10, 4, 8]
data = [Decimal(x) for x in data]
print(pearson(data, data))

请注意,您还需要使用小数计算平均值和标准差。
蒂姆·彼得斯解释说:
默认情况下,小数使用的精度比原生二进制浮点数高。每个有限精度的sqrt都必须是多对一函数。对于其他函数:Decimal(0.5)与Decimal(“0.5”)相同,因为0.5恰好可以表示为二进制浮点。您应该使用sqrt()而不是0.5,因为sqrt()保证将结果正确舍入到全精度,而不能。

相关问题