W przeszłości stosowałem następujące podejście do umiarkowanie wydajnego obliczania odchylenia rozgrzeszenia (zauważ, że jest to podejście programistów, a nie statystyków, więc niewątpliwie mogą istnieć sprytne sztuczki, takie jak shabbychef, które mogą być bardziej wydajne).
OSTRZEŻENIE: To nie jest algorytm online. Wymaga O(n)pamięci. Co więcej, ma najgorszy wynik w O(n)przypadku takich zestawów danych [1, -2, 4, -8, 16, -32, ...](tj. Taki sam jak w przypadku pełnego przeliczenia). [1]
Ponieważ jednak nadal działa dobrze w wielu przypadkach użycia, warto opublikować tutaj. Na przykład, aby obliczyć absolutne odchylenie 10000 losowo liczb od -100 do 100 w miarę dostarczania każdego elementu, mój algorytm zajmuje mniej niż jedną sekundę, podczas gdy pełne ponowne obliczenie zajmuje ponad 17 sekund (na mojej maszynie będą się różnić w zależności od maszyny i zgodnie z danymi wejściowymi). Musisz jednak zachować cały wektor w pamięci, co może być ograniczeniem dla niektórych zastosowań. Zarys algorytmu jest następujący:
- Zamiast pojedynczego wektora do przechowywania poprzednich pomiarów, użyj trzech posortowanych kolejek priorytetowych (coś w rodzaju sterty min / max). Te trzy listy dzielą dane wejściowe na trzy: pozycje większe niż średnia, pozycje mniejsze niż średnia i pozycje równe średniej.
- (Prawie) za każdym razem, gdy dodajesz przedmiot, zmienia się średnia, więc musimy podzielić na części. Kluczową sprawą jest posortowana natura partycji, co oznacza, że zamiast skanować każdy element na liście do podziału, musimy tylko czytać te elementy, które przenosimy. W najgorszym przypadku będzie to nadal wymagać
O(n) operacji przenoszenia, w wielu przypadkach użycia tak nie jest.
- Używając sprytnej księgowości, możemy upewnić się, że odchylenie jest poprawnie obliczane przez cały czas, podczas podziału partycji i dodawania nowych elementów.
Przykładowy kod w pythonie znajduje się poniżej. Pamiętaj, że pozwala tylko dodawać elementy do listy, a nie usuwać. Można to łatwo dodać, ale w chwili, gdy to pisałem, nie było takiej potrzeby. Zamiast samodzielnie wdrażać kolejki priorytetowe, skorzystałem z sortowanej listy z doskonałego pakietu blist Daniela Stutzbacha , który wykorzystuje B + Tree .
Rozważ ten kod na licencji MIT . Nie został znacznie zoptymalizowany ani dopracowany, ale działał dla mnie w przeszłości. Nowe wersje będą dostępne tutaj . Daj mi znać, jeśli masz jakieś pytania lub znajdziesz jakieś błędy.
from blist import sortedlist
import operator
class deviance_list:
def __init__(self):
self.mean = 0.0
self._old_mean = 0.0
self._sum = 0L
self._n = 0 #n items
# items greater than the mean
self._toplist = sortedlist()
# items less than the mean
self._bottomlist = sortedlist(key = operator.neg)
# Since all items in the "eq list" have the same value (self.mean) we don't need
# to maintain an eq list, only a count
self._eqlistlen = 0
self._top_deviance = 0
self._bottom_deviance = 0
@property
def absolute_deviance(self):
return self._top_deviance + self._bottom_deviance
def append(self, n):
# Update summary stats
self._sum += n
self._n += 1
self._old_mean = self.mean
self.mean = self._sum / float(self._n)
# Move existing things around
going_up = self.mean > self._old_mean
self._rebalance(going_up)
# Add new item to appropriate list
if n > self.mean:
self._toplist.add(n)
self._top_deviance += n - self.mean
elif n == self.mean:
self._eqlistlen += 1
else:
self._bottomlist.add(n)
self._bottom_deviance += self.mean - n
def _move_eqs(self, going_up):
if going_up:
self._bottomlist.update([self._old_mean] * self._eqlistlen)
self._bottom_deviance += (self.mean - self._old_mean) * self._eqlistlen
self._eqlistlen = 0
else:
self._toplist.update([self._old_mean] * self._eqlistlen)
self._top_deviance += (self._old_mean - self.mean) * self._eqlistlen
self._eqlistlen = 0
def _rebalance(self, going_up):
move_count, eq_move_count = 0, 0
if going_up:
# increase the bottom deviance of the items already in the bottomlist
if self.mean != self._old_mean:
self._bottom_deviance += len(self._bottomlist) * (self.mean - self._old_mean)
self._move_eqs(going_up)
# transfer items from top to bottom (or eq) list, and change the deviances
for n in iter(self._toplist):
if n < self.mean:
self._top_deviance -= n - self._old_mean
self._bottom_deviance += (self.mean - n)
# we increment movecount and move them after the list
# has finished iterating so we don't modify the list during iteration
move_count += 1
elif n == self.mean:
self._top_deviance -= n - self._old_mean
self._eqlistlen += 1
eq_move_count += 1
else:
break
for _ in xrange(0, move_count):
self._bottomlist.add(self._toplist.pop(0))
for _ in xrange(0, eq_move_count):
self._toplist.pop(0)
# decrease the top deviance of the items remain in the toplist
self._top_deviance -= len(self._toplist) * (self.mean - self._old_mean)
else:
if self.mean != self._old_mean:
self._top_deviance += len(self._toplist) * (self._old_mean - self.mean)
self._move_eqs(going_up)
for n in iter(self._bottomlist):
if n > self.mean:
self._bottom_deviance -= self._old_mean - n
self._top_deviance += n - self.mean
move_count += 1
elif n == self.mean:
self._bottom_deviance -= self._old_mean - n
self._eqlistlen += 1
eq_move_count += 1
else:
break
for _ in xrange(0, move_count):
self._toplist.add(self._bottomlist.pop(0))
for _ in xrange(0, eq_move_count):
self._bottomlist.pop(0)
# decrease the bottom deviance of the items remain in the bottomlist
self._bottom_deviance -= len(self._bottomlist) * (self._old_mean - self.mean)
if __name__ == "__main__":
import random
dv = deviance_list()
# Test against some random data, and calculate result manually (nb. slowly) to ensure correctness
rands = [random.randint(-100, 100) for _ in range(0, 1000)]
ns = []
for n in rands:
dv.append(n)
ns.append(n)
print("added:%4d, mean:%3.2f, oldmean:%3.2f, mean ad:%3.2f" %
(n, dv.mean, dv._old_mean, dv.absolute_deviance / dv.mean))
assert sum(ns) == dv._sum, "Sums not equal!"
assert len(ns) == dv._n, "Counts not equal!"
m = sum(ns) / float(len(ns))
assert m == dv.mean, "Means not equal!"
real_abs_dev = sum([abs(m - x) for x in ns])
# Due to floating point imprecision, we check if the difference between the
# two ways of calculating the asb. dev. is small rather than checking equality
assert abs(real_abs_dev - dv.absolute_deviance) < 0.01, (
"Absolute deviances not equal. Real:%.2f, calc:%.2f" % (real_abs_dev, dv.absolute_deviance))
[1] Jeśli objawy utrzymują się, skontaktuj się z lekarzem.
