88

是否有一种算法可以估计一组值的中值、众数、偏度和/或峰度，但这不需要一次将所有值存储在内存中？

我想计算基本统计数据：

  • 平均值：算术平均值
  • 方差：与平均值的平方偏差的平均值
  • 标准差：方差的平方根
  • 中位数：将较大一半的数字与较小的一半分开的值
  • 模式：在集合中找到的最频繁的值
  • 偏度：tl；博士
  • 峰度：tl；博士

计算任何这些的基本公式是小学算术，我知道它们。也有许多实现它们的统计库。

我的问题是我正在处理的集合中有大量（数十亿）值：在 Python 中工作，我不能只用数十亿个元素制作一个列表或散列。即使我用 C 语言编写了这个，十亿元素的数组也不太实用。

数据未排序。它是由其他过程随机、即时生成的。每组的大小是高度可变的，并且大小不会提前知道。

我已经想出了如何很好地处理均值和方差，以任何顺序遍历集合中的每个值。（实际上，在我的例子中，我按照它们生成的顺序来使用它们。）这是我正在使用的算法，由http://en.wikipedia.org/wiki/Algorithms_for_calculating_variance#On-line_algorithm提供：

  • 初始化三个变量：count、sum、sum_of_squares
  • 对于每个值：
    • 递增计数。
    • 将值添加到总和。
    • 将值的平方添加到 sum_of_squares。
  • 将总和除以计数，存储为变量均值。
  • 将 sum_of_squares 除以计数，存储为变量 mean_of_squares。
  • 平方均值，存储为 square_of_mean。
  • 从 mean_of_squares 中减去 square_of_mean，存储为方差。
  • 输出均值和方差。

这种“在线”算法有弱点（例如，准确性问题，因为 sum_of_squares 迅速增长到大于整数范围或浮点精度），但它基本上可以满足我的需求，而不必存储每个集合中的每个值。

但我不知道是否存在用于估计附加统计数据（中位数、众数、偏度、峰度）的类似技术。只要处理 N 个值所需的内存大大小于 O(N)，我就可以使用有偏差的估计器，甚至可以使用在一定程度上损害准确性的方法。

如果该库具有“在线”计算这些操作中的一个或多个的功能，那么将我指向现有的统计库也会有所帮助。

4

14 回答 14

57

我使用这些增量/递归均值和中值估计器，它们都使用常量存储：

mean += eta * (sample - mean)
median += eta * sgn(sample - median)

其中eta是一个小的学习率参数（例如 0.001），sgn () 是返回 {-1, 0, 1} 之一的符号函数。（如果数据是非平稳的，并且您想跟踪随时间的变化，请使用常数eta ；否则，对于固定源，您可以使用eta =1/n 之类的东西作为平均估计量，其中 n 是看到的样本数远......不幸的是，这似乎不适用于中位数估计器。）

这种类型的增量均值估计器似乎到处都在使用，例如在无监督神经网络学习规则中，但中值版本似乎不太常见，尽管它有好处（对异常值的鲁棒性）。在许多应用中，中值版本似乎可以替代均值估计量。

我很想看到一个类似形式的增量模式估计器......

更新 (2011-09-19)

我刚刚修改了增量中位数估计器来估计任意分位数。通常，分位数函数告诉您将数据分成两个部分的值：p 和 1-p。下面逐步估计这个值：

quantile += eta * (sgn(sample - quantile) + 2.0 * p - 1.0)

p 值应在 [0,1] 范围内。这实质上将sgn () 函数的对称输出 {-1,0,1} 向一侧倾斜，将数据样本划分为两个大小不等的 bin（数据的部分 p 和 1-p 小于/大于分位数估计，分别）。请注意，对于 p=0.5，这减少到中值估计量。

更新（2021-11-19）

有关此处描述的中值估计器的更多详细信息，我想强调以下评论中链接的这篇论文：Bylander & Rosen，1997，A Perceptron-Like Online Algorithm for Tracking the Median。这是作者网站的附言版本。

于 2010-01-27T05:24:00.633 回答
54

偏度和峰度

对于 Skewness 和 Kurtosis 的在线算法（沿着方差的线），请参阅此处的同一 wiki 页面中的更高矩统计的并行算法。

中位数

没有排序的数据，中位数很难。如果你知道你有多少数据点，理论上你只需要部分排序，例如通过使用选择算法。但是，这对数十亿的价值并没有太大帮助。我建议使用频率计数，请参阅下一节。

具有频率计数的中值和众数

如果它是整数，我会计算 频率，可能会在我确信它不再相关的某个值之外切断最高和最低值。对于浮点数（或太多整数），我可能会创建桶/间隔，然后使用与整数相同的方法。基于频率表，（近似）模式和中位数计算变得容易。

正态分布随机变量

如果它是正态分布的，我将使用总体样本均值、方差、偏度和峰度作为小子集的最大似然估计量。计算这些的（在线）算法，你现在已经知道了。例如，读取数十万或数百万个数据点，直到您的估计误差变得足够小。只需确保您从集合中随机选择（例如，您不会通过选择前 100'000 个值来引入偏差）。同样的方法也可用于估计正常情况的众数和中位数（因为样本均值都是估计量）。

进一步的评论

如果有帮助，上述所有算法都可以并行运行（包括许多排序和选择算法，例如 QuickSort 和 QuickSelect）。

我一直假设（除了关于正态分布的部分）我们谈论的是样本矩、中位数和众数，而不是给定已知分布的理论矩的估计量。

一般来说，在给定数据量的情况下，对数据进行采样（即只查看一个子集）应该非常成功，只要所有观察结果都是相同随机变量（具有相同分布）以及矩、模式和这个分布实际上存在中位数。最后一个警告并非无害。例如，不存在柯西分布的均值（和所有更高的矩）。在这种情况下，“小”子集的样本均值可能与整个样本的样本均值大相径庭。

于 2009-06-29T16:14:04.713 回答
13

我在我编写的一个名为LiveStats的简洁 Python 模块中实现了用于动态计算分位数和直方图的 P 方算法，而无需存储观测值。它应该可以非常有效地解决您的问题。该库支持您提到的所有统计信息，但模式除外。我还没有找到一个令人满意的模式估计解决方案。

于 2013-03-09T16:34:37.297 回答
7

瑞安，恐怕你没有做正确的平均值和方差......这是几周前出现的。在线版本（实际上以 Welford 方法的名称命名）的优点之一是它特别准确和稳定，请参阅此处的讨论。优点之一是您不需要存储总和或总平方和......

我想不出任何关于众数和中位数的在线方法，这似乎需要一次考虑整个列表。但很可能，与方差和均值类似的方法也适用于偏度和峰度......

于 2009-06-29T17:55:05.500 回答
3

问题中引用的维基百科文章包含在线计算偏度和峰度的公式。

对于模式 - 我相信 - 没有办法在线进行。为什么？假设您输入的所有值都不同，除了与前一个重复的最后一个值。在这种情况下，您必须记住在输入中已经看到的所有值，以检测最后一个值是否与之前看到的值重复并使其成为最常见的值。

对于中位数，它几乎是相同的 - 直到最后一个输入，如果所有输入值都不同，您不知道哪个值将成为中位数，因为它可能在当前中位数之前或之后。如果您知道输入的长度，您可以在不将所有值存储在内存中的情况下找到中值，但您仍然需要存储其中的许多（我猜大约是一半），因为错误的输入序列可能会在后半部分可能使前半部分的任何值成为中位数。

（请注意，我仅指精确计算。）

于 2009-06-29T16:03:37.820 回答
2

如果您有数十亿个数据点，那么您不太可能需要准确的答案，而不是接近的答案。通常，如果您有数十亿个数据点，则生成它们的基础过程可能会遵循某种统计平稳性/遍历性/混合属性。此外，您是否期望分布合理连续也可能很重要。

在这些情况下，如果您不需要确切的答案，则存在用于在线、低内存、分位数估计（中位数是 0.5 分位数的特殊情况）以及众数的算法。这是一个活跃的统计领域。

分位数估计示例： http: //www.computer.org/portal/web/csdl/doi/10.1109/WSC.2006.323014

模式估计示例：Bickel DR。连续数据的模式和偏度的稳健估计。计算统计和数据分析。2002；39：153–163。doi: 10.1016/S0167-9473(01)00057-3。

这些是计算统计的活跃领域。您正在进入的领域没有任何单一的最佳精确算法，而是它们的多样性（实际上是统计估计器），它们具有不同的属性、假设和性能。是实验数学。可能有成百上千篇关于这个主题的论文。

最后一个问题是您是否真的需要偏度和峰度本身，或者更可能需要一些其他参数，这些参数在表征概率分布时可能更可靠（假设您有一个概率分布！）。你期待高斯吗？

您是否有清理/预处理数据以使其主要是高斯的方法？（例如，金融交易金额在取对数后通常有点高斯）。你期望有限的标准偏差吗？你期待肥尾巴吗？您关心的数量是尾部还是散装？

于 2009-10-18T23:08:13.717 回答
2

每个人都在说您不能以在线方式执行该模式，但事实并非如此。这是一篇文章，描述了耶鲁大学的 Michael E. Fischer 和 Steven L. Salzberg 于 1982 年发明的解决这个问题的算法。来自文章：

多数查找算法使用它的一个寄存器来临时存储流中的单个项目；该项目是多数元素的当前候选者。第二个寄存器是一个初始化为 0 的计数器。对于流的每个元素，我们要求算法执行以下例程。如果计数器读数为 0，则将当前流元素安装为新的多数候选（替换可能已经在寄存器中的任何其他元素）。然后，如果当前元素与多数候选匹配，则递增计数器；否则，递减计数器。在循环的这一点上，如果到目前为止看到的部分流具有多数元素，则该元素在候选寄存器中，并且计数器保持大于 0 的值。如果没有多数派怎么办？如果不对数据进行第二次传递——这在流环境中是不可能的——算法在这种情况下不能总是给出明确的答案。它只是承诺正确识别多数元素（如果有的话）。

它也可以扩展以找到具有更多内存的前 N ​​，但这应该可以解决该模式。

于 2012-02-19T17:22:20.867 回答
1

最终，如果您对分布没有先验参数知识，我认为您必须存储所有值。

也就是说，除非您正在处理某种病理情况，否则补救措施（Rousseuw and Bassett 1990）可能足以满足您的目的。

非常简单，它涉及计算批次中位数的中位数。

于 2009-07-27T21:14:47.043 回答
0

仅使用可用的恒定空间无法在线计算中位数和众数。但是，由于中位数和众数无论如何都比“定量”更具“描述性”，因此您可以通过对数据集进行采样来估计它们。

如果数据从长远来看是正态分布的，那么你可以用你的平均值来估计中位数。

您还可以使用以下技术估计中值：为数据流中的每个（例如）1,000,000 个条目建立中值估计 M[i]，以便 M[0] 是前一百万个条目的中值，M[1]第二个一百万条目的中位数等。然后使用 M[0]...M[k] 的中位数作为中位数估计量。这当然可以节省空间，您可以通过“调整”参数 1,000,000 来控制要使用多少空间。这也可以递归地推广。

于 2009-06-29T16:18:00.017 回答
0

好的伙计试试这些：

对于 C++：

double skew(double* v, unsigned long n){
    double sigma = pow(svar(v, n), 0.5);
    double mu = avg(v, n);

    double* t;
    t = new double[n];

    for(unsigned long i = 0; i < n; ++i){
        t[i] = pow((v[i] - mu)/sigma, 3);
    }

    double ret = avg(t, n);

    delete [] t;
    return ret;
}

double kurt(double* v, double n){
    double sigma = pow(svar(v, n), 0.5);
    double mu = avg(v, n);

    double* t;
    t = new double[n];

    for(unsigned long i = 0; i < n; ++i){
        t[i] = pow( ((v[i] - mu[i]) / sigma) , 4) - 3;
    }

    double ret = avg(t, n);

    delete [] t;
    return ret;
}

您说您已经可以计算样本方差（svar）和平均值（avg），您将它们指向您的函数以执行此操作。

另外，看看皮尔逊的近似值。在如此大的数据集上，它会非常相似。3（平均值 - 中位数）/标准差您的中位数为最大值 - 最小值/2

对于浮动模式没有意义。人们通常会将它们放在尺寸很大的垃圾箱中（例如 1/100 * (max - min)）。

于 2013-01-29T13:37:50.900 回答
0

Pebay 等人解决了这个问题：

https://prod-ng.sandia.gov/techlib-noauth/access-control.cgi/2008/086212.pdf

于 2019-03-18T23:56:10.080 回答
0

中位数

可以在此处找到两个最近的百分位近似算法及其 python 实现：

t-文摘

DDS素描

两种算法都存储数据。由于 T-Digest 在尾部附近使用较小的 bin，因此极端情况下的准确性更好（接近中位数的情况较弱）。DDSketch 还提供相对错误保证。

于 2021-10-06T17:23:17.917 回答
-1

我倾向于使用桶，这可能是自适应的。桶大小应该是您需要的精度。然后，当每个数据点进入时，您将一个添加到相关存储桶的计数中。这些应该为您提供中位数和峰度的简单近似值，方法是将每个桶计算为其计数加权的值。

一个问题可能是数十亿次运算后浮点分辨率的损失，即加一不会再改变值！为了解决这个问题，如果最大存储桶大小超过某个限制，您可以从所有计数中取出一个很大的数字。

于 2010-10-28T23:38:11.960 回答
-1 
for j in range (1,M):
    y=np.zeros(M) # build the vector y
    y[0]=y0

    #generate the white noise
    eps=npr.randn(M-1)*np.sqrt(var)

    #increment the y vector
    for k in range(1,T):
        y[k]=corr*y[k-1]+eps[k-1]

    yy[j]=y

list.append(y)
于 2016-01-18T14:27:41.890 回答