为什么我的峰度函数不能产生与scipy.stats.kurtosis相同的输出？

Question

我有一个家庭作业的问题，我应该为峰度写一个函数，如下所述：

​

​

分母中的θ是标准差(方差的平方根)，分子中的x-with-bar是x的平均值。

我已经实现了如下函数：

import numpy as np
from scipy.stats import kurtosis

testdata = np.array([1, 2, 3, 4, 5])

def mean(obs):
    return (1. / len(obs)) * np.sum(obs)

def variance(obs):
    return (1. / len(obs)) * np.sum((obs - mean(obs)) ** 2)

def kurt(obs):
    num = np.sqrt((1. / len(obs)) * np.sum((obs - mean(obs)) ** 4))
    denom = variance(obs) ** 2  # avoid losing precision with np.sqrt call
    return num / denom

前两个函数mean和variance分别与numpy.mean和numpy.var成功交叉验证。

我尝试使用以下语句交叉验证kurt：

>>> kurtosis(testdata) == kurt(testdata)
False

下面是两个峰度函数的输出：

>>> kurtosis(testdata)  # scipy.stats
-1.3

>>> kurt(testdata)  # my crappy attempt
0.65192024052026476

我哪里错了？scipy.stats.kurtosis是不是做了一些比我得到的等式中的内容更奇特的事情？

Answer 1

缺省情况下，scipy.stats.kurtosis()

1. 计算超额峰度(即，从统计偏差的result).
2. Corrects中减去3(这会影响一些denominators).

这两种行为都可以通过scipy.stats.kurtosis()的可选参数进行配置。

最后，方法中的np.sqrt()调用是不必要的，因为公式中没有平方根。一旦我删除它，你的函数的输出就会与我从kurtosis(testdata, False, False)得到的结果相匹配。

I尝试使用以下语句交叉验证kurt

你不应该为了精确相等而比较浮点数。即使数学公式是相同的，它们如何翻译成计算机代码的微小差异也会影响计算结果。

最后，如果您打算编写数字代码，我强烈建议您阅读What Every Computer Scientist Should Know About Floating-Point Arithmetic。

附注:这是我用过的函数：

In [51]: def kurt(obs):
   ....:     num = np.sum((obs - mean(obs)) ** 4)/ len(obs)
   ....:     denom = variance(obs) ** 2  # avoid losing precision with np.sqrt call
   ....:     return num / denom