如何绘制经验cdf (ecdf)

Question

如何在Python的matplotlib中绘制一组数字的经验CDF？我正在寻找与pylab的"hist“函数类似的cdf。

我能想到的一件事是：

from scipy.stats import cumfreq
a = array([...]) # my array of numbers
num_bins =  20
b = cumfreq(a, num_bins)
plt.plot(b)

Answer 1

这看起来(几乎)就是你想要的。两件事：

首先，结果是一个包含四个项目的元组。第三个是垃圾箱的大小。第二个是最小bin的起点。第一个是每个bin中或其下的点的数量。(最后一个是超出限制的点数，但由于您没有设置任何点数，所有点数都将被二进制。)

其次，您需要重新缩放结果，以便最终值为1，以遵循CDF的通常约定，但除此之外，它是正确的。

下面是它在幕后的作用：

def cumfreq(a, numbins=10, defaultreallimits=None):
    # docstring omitted
    h,l,b,e = histogram(a,numbins,defaultreallimits)
    cumhist = np.cumsum(h*1, axis=0)
    return cumhist,l,b,e

它执行组织图，然后生成每个bin中计数的累积和。因此，结果的第i个值是小于或等于第i个bin的最大值的数组值的数量。所以，最终的值就是初始数组的大小。

最后，要绘制它，需要使用bin的初始值和bin大小来确定需要的x轴值。

另一种选择是使用numpy.histogram，它可以进行归一化并返回bin边缘。您将需要自己计算结果计数的累积和。

a = array([...]) # your array of numbers
num_bins = 20
counts, bin_edges = numpy.histogram(a, bins=num_bins, normed=True)
cdf = numpy.cumsum(counts)
pylab.plot(bin_edges[1:], cdf)

(bin_edges[1:]是每个存储箱的上边缘。)

Answer 2

如果你喜欢linspace并且更喜欢一行程序，你可以这样做：

plt.plot(np.sort(a), np.linspace(0, 1, len(a), endpoint=False))

考虑到我的品味，我几乎总是这样做：

# a is the data array
x = np.sort(a)
y = np.arange(len(x))/float(len(x))
plt.plot(x, y)

即使有>O(1e6)数据值，它也适用于我。如果你真的需要下采样，我会设置

x = np.sort(a)[::down_sampling_step]

编辑以响应有关我使用endpoint=False或上面定义的y的原因的注释/编辑。以下是一些技术细节。

经验CDF通常被正式定义为

CDF(x) = "number of samples <= x"/"number of samples"

为了精确地匹配这个正式的定义，你需要使用y = np.arange(1,len(x)+1)/float(len(x))，这样我们才能得到y = [1/N, 2/N ... 1]。这个估计器是一个无偏的估计器，它将在无限样本Wikipedia ref.的限制下收敛到真正的CDF。

我倾向于使用y = [0, 1/N, 2/N ... (N-1)/N]，因为(a)它更容易编码/更地道，(b)但仍然是正式合理的，因为人们总是可以在收敛证明中将CDF(x)与1-CDF(x)交换，以及(c)使用上述(简单)下采样方法。

在某些特定情况下，定义

y = (arange(len(x))+0.5)/len(x)

它是这两个约定之间的中间部分。这实际上是说，“有一个1/(2N)概率小于我在我的样本中见过的最低值，而一个1/(2N)概率值大于我到目前为止见过的最大值。

请注意，如果将where显示为分段常量函数更有用，则此约定的选择将与plt.step中使用的CDF参数进行交互。为了与上面提到的正式定义完全匹配，人们需要使用建议的y=[0,1/N..., 1-1/N]约定，或者将where=post与y=[1/N, 2/N ... 1]约定一起使用，但不是反过来。

但是，对于大样本和合理分布，正文中给出的约定答案很容易写，是真实CDF的无偏估计器，并与下采样方法一起工作。

Answer 3

您可以使用scikits.statsmodels库中的ECDF函数：

import numpy as np
import scikits.statsmodels as sm
import matplotlib.pyplot as plt

sample = np.random.uniform(0, 1, 50)
ecdf = sm.tools.ECDF(sample)

x = np.linspace(min(sample), max(sample))
y = ecdf(x)
plt.step(x, y)

在0.4版中，scicits.statsmodels被重命名为statsmodels。ECDF现在位于distributions模块中(而statsmodels.tools.tools.ECDF已折旧)。

import numpy as np
import statsmodels.api as sm # recommended import according to the docs
import matplotlib.pyplot as plt

sample = np.random.uniform(0, 1, 50)
ecdf = sm.distributions.ECDF(sample)

x = np.linspace(min(sample), max(sample))
y = ecdf(x)
plt.step(x, y)
plt.show()