如何有效地计算运行标准偏差？
我有一系列的数字列表，例如：
[0] (0.01, 0.01, 0.02, 0.04, 0.03) [1] (0.00, 0.02, 0.02, 0.03, 0.02) [2] (0.01, 0.02, 0.02, 0.03, 0.02) ... [n] (0.01, 0.00, 0.01, 0.05, 0.03) 我想要做的是有效地计算列表中每个索引处的所有数组元素的均值和标准差。
为了做到这一点，我一直循环遍历数组，并在列表的给定索引处求和值。 最后，我把我的“平均值列表”中的每个值除以n 。
为了做标准差，我再次循环，现在我已经计算了平均值。
我想避免两次，一次是平均值，然后是一次为SD（在我有意思之后）。
有没有一个有效的方法来计算两个值，只有通过数组一次？ 解释型语言（例如Perl或Python）或伪代码中的任何代码都可以。
答案是使用Welford的algorithm，这个algorithm在“天真的方法”之后非常明确的定义在：
- 维基百科： 计算方差的algorithm
它比其他答案中所提出的两次通过或在线简单的平方收集器build议的更稳定。 当浮点数文献中存在大量相互接近的值时，稳定性才真正重要，因为它们导致了所谓的“ 灾难性消除 ”。
您可能还想要在方差计算（平方偏差）中对采样数（N）和N-1除以之间的差异进行细化。 除以N-1导致样本方差的无偏估计，而平均N除以样本平均值和真实平均值之间的差异，则将N除以平均值低估方差。
我写了两个关于这个主题的博客文章，包括如何删除以前的在线值：
- 一次在线计算样本均值和方差
- 删除在线均值和方差的Welfordalgorithm中的值
你也可以看看我的Java工具; javadoc，源代码和unit testing都在线：
- Javadoc：
stats.OnlineNormalEstimator - 来源：
stats.OnlineNormalEstimator.java - JUnit来源：
test.unit.stats.OnlineNormalEstimatorTest.java - LingPipe主页
基本的答案是在你去的时候累计x （称为'sum_x1'）和x 2 （称为'sum_x2'）的总和。 那么标准差的值就是：
stdev = sqrt((sum_x2 / n) - (mean * mean))
哪里
mean = sum_x / n
这是样本标准偏差; 您使用“n”而不是“n – 1”作为除数得到人口标准差。
如果您处理的是大样本，则可能需要担心两个大数之间差异的数值稳定性。 转到其他答案（维基百科等）的外部参考了解更多信息。
也许不是你所要求的，但是…如果你使用一个numpy数组，它将为你有效地工作：
from numpy import array nums = array(((0.01, 0.01, 0.02, 0.04, 0.03), (0.00, 0.02, 0.02, 0.03, 0.02), (0.01, 0.02, 0.02, 0.03, 0.02), (0.01, 0.00, 0.01, 0.05, 0.03))) print nums.std(axis=1) # [ 0.0116619 0.00979796 0.00632456 0.01788854] print nums.mean(axis=1) # [ 0.022 0.018 0.02 0.02 ]
顺便说一下，在这篇博文中有一些有趣的讨论，并且提供了关于计算方法和差异的单程方法的评论：
这里是从http://www.johndcook.com/standard_deviation.html的Welford的algorithm实现的纯粹的Python翻译：;
https://github.com/liyanage/python-modules/blob/master/running_stats.py
class RunningStats: def __init__(self): self.n = 0 self.old_m = 0 self.new_m = 0 self.old_s = 0 self.new_s = 0 def clear(self): self.n = 0 def push(self, x): self.n += 1 if self.n == 1: self.old_m = self.new_m = x self.old_s = 0 else: self.new_m = self.old_m + (x - self.old_m) / self.n self.new_s = self.old_s + (x - self.old_m) * (x - self.new_m) self.old_m = self.new_m self.old_s = self.new_s def mean(self): return self.new_m if self.n else 0.0 def variance(self): return self.new_s / (self.n - 1) if self.n > 1 else 0.0 def standard_deviation(self): return math.sqrt(self.variance())
用法：
rs = RunningStats() rs.push(17.0); rs.push(19.0); rs.push(24.0); mean = rs.mean(); variance = rs.variance(); stdev = rs.standard_deviation();
看看PDL （发音为“piddle！”）。
这是用于高精度math和科学计算的Perl数据语言。
这是一个使用你的数字的例子….
use strict; use warnings; use PDL; my $figs = pdl [ [0.01, 0.01, 0.02, 0.04, 0.03], [0.00, 0.02, 0.02, 0.03, 0.02], [0.01, 0.02, 0.02, 0.03, 0.02], [0.01, 0.00, 0.01, 0.05, 0.03], ]; my ( $mean, $prms, $median, $min, $max, $adev, $rms ) = statsover( $figs ); say "Mean scores: ", $mean; say "Std dev? (adev): ", $adev; say "Std dev? (prms): ", $prms; say "Std dev? (rms): ", $rms;
其中产生：
Mean scores: [0.022 0.018 0.02 0.02] Std dev? (adev): [0.0104 0.0072 0.004 0.016] Std dev? (prms): [0.013038405 0.010954451 0.0070710678 0.02] Std dev? (rms): [0.011661904 0.009797959 0.0063245553 0.017888544]
有关statsover函数的更多信息，请参阅 PDL :: Primitive 。 这似乎表明，ADEV是“标准偏差”。
然而，它可能是PRMS（Sinan的Statistics :: Descriptive示例显示）或RMS（其中的NumPy示例显示）。 我猜这三个中的一个必须是正确的;-)
更多PDL信息请看：
- pdl.perl.org （官方PDL页面）。
- PerlMonks上的PDL快速参考指南
- Dobb博士关于PDL的文章
- PDL Wiki
- Wikipedia条目为PDL
- PDL的Sourceforge项目页面
Statistics :: Descriptive是用于这些types计算的非常体面的Perl模块：
#!/usr/bin/perl use strict; use warnings; use Statistics::Descriptive qw( :all ); my $data = [ [ 0.01, 0.01, 0.02, 0.04, 0.03 ], [ 0.00, 0.02, 0.02, 0.03, 0.02 ], [ 0.01, 0.02, 0.02, 0.03, 0.02 ], [ 0.01, 0.00, 0.01, 0.05, 0.03 ], ]; my $stat = Statistics::Descriptive::Full->new; # You also have the option of using sparse data structures for my $ref ( @$data ) { $stat->add_data( @$ref ); printf "Running mean: %f\n", $stat->mean; printf "Running stdev: %f\n", $stat->standard_deviation; } __END__
输出：
C:\Temp> g Running mean: 0.022000 Running stdev: 0.013038 Running mean: 0.020000 Running stdev: 0.011547 Running mean: 0.020000 Running stdev: 0.010000 Running mean: 0.020000 Running stdev: 0.012566
Python的runstats模块就是为了这样的事情。 从PyPI 安装runstats ：
pip install runstats
Runstats汇总可以产生一次数据中的均值，方差，标准偏差，偏度和峰度。 我们可以使用这个来创build你的“正在运行”的版本。
from runstats import Statistics stats = [Statistics() for num in range(len(data[0]))] for row in data: for index, val in enumerate(row): stats[index].push(val) for index, stat in enumerate(stats): print 'Index', index, 'mean:', stat.mean() print 'Index', index, 'standard deviation:', stat.stddev()
统计摘要是基于计算标准偏差的Knuth和Welford方法，如“计算机程序devise艺术”第2卷第1页所述。 232，第3版。 这样做的好处是数值稳定和准确的结果。
免责声明：我是Python runstats模块的作者。
你的arrays有多大？ 除非它是元素的无数长，不要担心两次循环。 该代码简单，易于testing。
我的首选是使用numpy数组math扩展将您的数组数组转换成numpy二维数组，并直接得到标准偏差：
>>> x = [ [ 1, 2, 4, 3, 4, 5 ], [ 3, 4, 5, 6, 7, 8 ] ] * 10 >>> import numpy >>> a = numpy.array(x) >>> a.std(axis=0) array([ 1. , 1. , 0.5, 1.5, 1.5, 1.5]) >>> a.mean(axis=0) array([ 2. , 3. , 4.5, 4.5, 5.5, 6.5])
如果这不是一个选项，你需要一个纯粹的Python解决scheme，请继续阅读…
如果你的数组是
x = [ [ 1, 2, 4, 3, 4, 5 ], [ 3, 4, 5, 6, 7, 8 ], .... ]
那么标准差是：
d = len(x[0]) n = len(x) sum_x = [ sum(v[i] for v in x) for i in range(d) ] sum_x2 = [ sum(v[i]**2 for v in x) for i in range(d) ] std_dev = [ sqrt((sx2 - sx**2)/N) for sx, sx2 in zip(sum_x, sum_x2) ]
如果您确定只循环一次数组，则可以合并运行总和。
sum_x = [ 0 ] * d sum_x2 = [ 0 ] * d for v in x: for i, t in enumerate(v): sum_x[i] += t sum_x2[i] += t**2
这不像上面的列表理解解决scheme那么优雅。
你可以看标准偏差的维基百科文章，特别是有关快速计算方法的章节。
还有一篇文章，我发现使用Python，你应该能够使用它的代码没有太多的变化： 潜意识信息 – 运行标准偏差 。
我认为这个问题会帮助你。 标准偏差
n=int(raw_input("Enter no. of terms:")) L=[] for i in range (1,n+1): x=float(raw_input("Enter term:")) L.append(x) sum=0 for i in range(n): sum=sum+L[i] avg=sum/n sumdev=0 for j in range(n): sumdev=sumdev+(L[j]-avg)**2 dev=(sumdev/n)**0.5 print "Standard deviation is", dev
正如下面的答案所描述的： pandas / scipy / numpy是否提供了一个累积的标准偏差函数？ Pythonpandas模块包含一个计算运行或累积标准偏差的方法 。 为此，您必须将数据转换为pandas数据框（如果是1D，则为一系列数据），但也有相应的function。
这是一个“单线”，分布在多个函数式编程风格的多行：
def variance(data, opt=0): return (lambda (m2, i, _): m2 / (opt + i - 1))( reduce( lambda (m2, i, avg), x: ( m2 + (x - avg) ** 2 * i / (i + 1), i + 1, avg + (x - avg) / (i + 1) ), data, (0, 0, 0)))
