Эффективная оценка медианы

Итак, у Вас есть какой-то поток данных. Большой такой поток. Или уже готовый набор. И хочется определить какие-то его характеристики. Алгоритм определения минимального и максимального значения могут придумать даже не программисты. Вычисление среднего уже чуть сложнее, но тоже не представляет никаких трудностей — знай подсчитывай себе сумму да инкрементируй счетчик на каждое новое значение. Среднеквадратичное отклонение — все то же самое, только числа другие. А как насчет медианы?

Для тех, кто забыл, что это такое, напоминаю — медиана (50-й перцентиль) выборки данных — это такое значение, которое делит эту выборку пополам — данные из одной половины имеют значение не меньше медианы, а из второй — не больше. Ценность её заключается в том, что её значение не зависит от величины случайных всплесков, которые могут очень сильно повлиять на среднее.

Строго говоря, из определения следует, что для вычисления точного значения медианы нам нужно хранить всю выборку, иначе нет никаких гарантий, что мы насчитали именно то, что хотели. Но для непрерывных и больших потоков данных точное значение все равно не имеет большого смысла — сейчас оно одно, а через новых 100 отсчетов — уже другое. Поэтому эффективный метод оценки медианы, который не будет требовать много памяти и ресурсов CPU, и будет давать точность порядка одного процента или лучше — как раз то что нужно.

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

Идея метода состоит в том, чтобы построить такой процесс вычисления, который будет сходиться к действительному значению медианы. Если мы уже обработали какой-то обьем данных и имеем какую-то оценку медианы, то про поступлении нового обьема (с почти такой же медианой, что важно) наша оценка должна быть улучшена. Если более точно — то оценка должна быть улучшена с большей вероятностью, чем ухудшена.

Можно использовать разного рода окна вычисления медианы, например, посчитать точную медиану последних 100 значений, и усреднить с постоянно уменьшающимся весом с предыдущей оценкой. Такой метод будет хорошо работать, но есть гораздо более легкий и практически такой же точный. А именно — просто сдвигать текущую оценку медианы к новому значению на какую-то небольшую константную дельту. В случае, если предыдущая оценка была не точной, то при обработке нового объема данных она приблизится к действительному значению, т.е. станет более точной. А если оценка уже и так точная, то на обработке нового объема данных на половине значений будет сдвиг в одну сторону, а на другой половине — в другую, в итоге оценка вернется к точному значению.

Важно, что дельта должна иметь одинаковое абсолютное значение для сдвигов как в большую, так и в меньшую сторону, иначе мы получим не 50-й перцентиль. Но теперь встает проблема подбора значения дельты — слишком маленькое даст медленную сходимость, а слишком большое — получим большую погрешность, особенно если дельта сравнима с самим значением медианы. Автоматическое вычисление дельты уже лишает её звания константы, но это и неважно, если в итоге мы получим лучший результат.

Имеет смысл делать дельту постоянно уменьшающейся, чтобы улучшить сходимость. Но не слишком быстро, иначе, при неблагоприятных условиях оценка рискует никогда не догнать действительное значение медианы. Коэфициент 1/count подходит идеально — он легко вычисляется, и ряд sum(1/n) — расходящийся, так что в конце-концов оценка достигнет действительной медианы, какой бы ни плохой она была изначально.

Привязывать значение дельты к текущей оценке медианы — рискованно. В неудачных условиях слишком заниженной оценке будет сложно расти. Лучше всего вычислять дельту исходя из другой вполне устойчивой характеристики выборки — среднего значения. С коэфициентом 1/count значение дельты будет постоянно уменьшаться (за исключением немногочисленных случаев, когда текущее среднее вырастет слишком сильно по сравнению с предыдущим). Для данных, которые могут принимать не только неотрицательные значение, такой подход уже не сработает, но мы можем легко выйти из положения, если считать два средних — положительных и отрицательных значений, и использовать сумму их абсолютных значений.

Итоговый алгоритм оценки медианы почти такой же простой, как и алгоритм вычисления среднего значения:

sum += value
count ++
delta = sum / count / count        // delta = average/count
if value < median {
  median -= delta
} else {
  median += delta
}

Погрешность и скорость сходимости, по моему мнению, у него отличная — на выборке в 100 значений погрешность обычно около 10-20%, на 1000 — уже около 1-3%, и дальше погрешность уменьшается ещё больше.

Для желающих самостоятельно оценить качество работы алгоритма — небольшой скриптик на perl-е:

#!/usr/bin/perl -s

$N = shift || 100000;
for (1 .. $N) {
  push @a, (rand(1)*rand(1))**2;
}
@t = sort { $a <=> $b } @a;
$MEDIAN = $t[$N/2-1];

median_init();
for (@a) {
  median_update($_);
}

$diff = ($MEDIAN-$median)/$MEDIAN*100;
$avg = $sum/$count;

print "Real:tt$MEDIAN
Estimated:t$median
Difference:t$diff %

Average: $avg
";

sub median_init() {
  $count = $sum = $median = 0;
}
sub median_update($) {
  $value = $_[0];
  $sum += $value;
  $count ++;
  $delta = $sum / $count / $count;
  if($median <= $value) {
    $median += $delta;
    $s2 = '+';
  } else {
    $median -= $delta;
    $s2 = '-';
  }
  if($debug) {
    $s = ($value > $MEDIAN) ? '+' : '-';
    printf "%s %.5ft%dt%.7ft  %s  %en", $s, $value, $count, $median, $s2, $delta;
  }
}

Автор: zuborg

Источник