Как эффективно вычислить рабочее стандартное отклонение?

У меня есть массив списков чисел, например:

[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) или псевдокод прекрасны.

80
задан Alex Reynolds 22 July 2019 в 18:00
поделиться

8 ответов

The answer is to use Welford's algorithm, which is very clearly defined after the "naive methods" in:

It's more numerically stable than either the two-pass or online simple sum of squares collectors suggested in other responses. The stability only really matters when you have lots of values that are close to each other as they lead to what is known as "catastrophic cancellation" in the floating point literature.

You might also want to brush up on the difference between dividing by the number of samples (N) and N-1 in the variance calculation (squared deviation). Dividing by N-1 leads to an unbiased estimate of variance from the sample, whereas dividing by N on average underestimates variance (because it doesn't take into account the variance between the sample mean and the true mean).

I wrote two blog entries on the topic which go into more details, including how to delete previous values online:

You can also take a look at my Java implement; the javadoc, source, and unit tests are all online:

106
ответ дан 24 November 2019 в 09:45
поделиться

Думаю, эта проблема вам поможет. Стандартное отклонение

2
ответ дан 24 November 2019 в 09:45
поделиться

Вы можете посмотреть статью в Википедии о стандартном отклонении , в частности раздел о методах быстрых вычислений.

Я также нашел статью, в которой используется Python, вы должен иметь возможность использовать код в нем без особых изменений: Подсознательные сообщения - текущие стандартные отклонения .

1
ответ дан 24 November 2019 в 09:45
поделиться

Насколько велик ваш массив? Если это не миллионы элементов, не беспокойтесь о повторении цикла дважды. Код прост и легко тестируется.

Я бы предпочел использовать математическое расширение массива numpy , чтобы преобразовать ваш массив массивов в numpy 2D-массив и напрямую получить стандартное отклонение:

>>> 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, продолжайте читать ...

Если ваш массив равен

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

Это далеко не так элегантно, как решение для понимания списка выше.

3
ответ дан 24 November 2019 в 09:45
поделиться

Возможно, это не то, о чем вы спрашивали, но ... Если вы используете массив 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 ]

Между прочим, в это сообщение в блоге и комментарии к однопроходным методам вычисления средств и отклонений:

26
ответ дан 24 November 2019 в 09:45
поделиться

Statistics :: Descriptive - очень приличный модуль 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
8
ответ дан 24 November 2019 в 09:45
поделиться

Основной ответ - накопить сумму как x (назовите это «сумма_x1»), так и x 2 (назовите это » sum_x2 ') по мере продвижения. Тогда значение стандартного отклонения будет:

stdev = sqrt((sum_x2 / n) - (mean * mean)) 

, где

mean = sum_x / n

Это стандартное отклонение выборки; вы получаете стандартное отклонение генеральной совокупности, используя в качестве делителя «n» вместо «n - 1».

Возможно, вам придется побеспокоиться о числовой стабильности измерения разницы между двумя большими числами, если вы имеете дело с большими выборками. Для получения дополнительной информации перейдите к внешним ссылкам в других ответах (Википедия и т. Д.).

71
ответ дан 24 November 2019 в 09:45
поделиться

Взгляните на PDL (произносится как «пиддл!»).

Это язык данных 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]


Взгляните на PDL :: Primitive для получения дополнительной информации о функции statsover . Похоже, это предполагает, что ADEV - это «стандартное отклонение».

Однако это может быть PRMS (что показано в примере Sinan Statistics :: Descriptive) или RMS (что показано в примере Ars NumPy). Думаю, один из этих трех должен быть прав; -)

Для получения дополнительной информации о PDL посмотрите:

8
ответ дан 24 November 2019 в 09:45
поделиться
Другие вопросы по тегам:

Похожие вопросы: