Нахождение медианы списка может казаться тривиальной задачей, но её выполнение за линейное время требует серьёзного подхода. В этом посте я расскажу об одном из самых любимых мной алгоритмов — нахождении медианы списка за детерминированное линейное время с помощью медианы медиан. Хотя доказательство того, что этот алгоритм выполняется за линейное время, довольно сложно, сам пост будет понятен и читателям с начальным уровнем знаний об анализе алгоритмов.
Нахождение медианы за O(n log n)
Самым прямолинейным способом нахождения медианы является сортировка списка и выбор медианы по её индексу. Самая быстрая сортировка сравнением выполняется за
O(n log n), поэтому от неё зависит время выполнения1, 2.def nlogn_median(l):
l = sorted(l)
if len(l) % 2 == 1:
return l[len(l) / 2]
else:
return 0.5 * (l[len(l) / 2 - 1] + l[len(l) / 2])У этого способа самый простой код, но он определённо не самый быстрый.
Нахождение медианы за среднее время O(n)
Следующим нашим шагом будет нахождение медианы в среднем за линейное время, если нам будет везти. Этот алгоритм, называемый «quickselect», разработан Тони Хоаром, который также изобрёл алгоритм сортировки с похожим названием — quicksort. Это рекурсивный алгоритм, и он может находить любой элемент (не только медиану).
- Выберем индекс списка. Способ выбора не важен, на практике вполне подходит и случайный. Элемент с этим индексом называется опорным элементом (pivot).
- Разделим список на две группы:
- Элементы меньше или равные pivot,
lesser_els - Элементы строго большие, чем pivot,
great_els
- Элементы меньше или равные pivot,
- Мы знаем, что одна из этих групп содержит медиану. Предположим, что мы ищем k-тый элемент:
- Если в
lesser_elsесть k или больше элементов, рекурсивно обходим списокlesser_elsв поисках k-того элемента. - Если в
lesser_elsменьше, чем k элементтов, рекурсивно обходим списокgreater_els. Вместо поиска k мы ищемk-len(lesser_els).
- Если в
Вот пример алгоритма, выполняемого для 11 элементов:
Возьмём представленный ниже список. Мы хотим найти медиану.
l = [9,1,0,2,3,4,6,8,7,10,5]
len(l) == 11, поэтому мы ищем шестой наименьший элемент
Сначала нам нужно выбрать опорный элемент (pivot). Мы случайным образом выбираем индекс 3.
Значение элемента с этим индексом равно 2.
Разбиваем список на группы согласно pivot:
[1,0,2], [9,3,4,6,8,7,10,5]
Нам нужен шестой элемент. 6-len(left) = 3, поэтому нам нужен
третий наименьший элемент в правом массиве
Теперь мы ищем третий наименьший элемент в следующем массиве:
[9,3,4,6,8,7,10,5]
Мы случайным образом выбираем индекс, который будет нашим pivot.
Мы выбрали индекс 2, значение в котором равно l[2]=6
Разбиваем на группы согласно pivot:
[3,4,5,6] [9,7,10]
Нам нужен третий наименьший элемент, поэтому мы знаем, что это
третий наименьший элемент в левом массиве
Теперь мы ищем третий наименьший в следующем массиве:
[3,4,5,6]
Мы случайным образом выбираем индекс, который будет нашим pivot.
Мы выбрали индекс 1, значение в котором равно l[1]=4
Разбиваем на группы согласно pivot:
[3,4] [5,6]
Нам нужен третий наименьший элемент, поэтому мы знаем, что это
наименьший элемент в правом массиве.
Теперь мы ищем наименьший элемент в следующем массиве:
[5,6]
На этом этапе у нас есть базовый вариант, выбирающий наибольший
или наименьший элемент на основании индекса.
Нам нужен наименьший элемент, то есть 5.
return 5Чтобы найти с помощью quickselect медиану, мы выделим quickselect в отдельную функцию. Наша функция
quickselect_median будет вызывать quickselect с нужными индексами.import random
def quickselect_median(l, pivot_fn=random.choice):
if len(l) % 2 == 1:
return quickselect(l, len(l) / 2, pivot_fn)
else:
return 0.5 * (quickselect(l, len(l) / 2 - 1, pivot_fn) +
quickselect(l, len(l) / 2, pivot_fn))
def quickselect(l, k, pivot_fn):
"""
Выбираем k-тый элемент в списке l (с нулевой базой)
:param l: список числовых данных
:param k: индекс
:param pivot_fn: функция выбора pivot, по умолчанию выбирает случайно
:return: k-тый элемент l
"""
if len(l) == 1:
assert k == 0
return l[0]
pivot = pivot_fn(l)
lows = [el for el in l if el < pivot]
highs = [el for el in l if el > pivot]
pivots = [el for el in l if el == pivot]
if k < len(lows):
return quickselect(lows, k, pivot_fn)
elif k < len(lows) + len(pivots):
# Нам повезло и мы угадали медиану
return pivots[0]
else:
return quickselect(highs, k - len(lows) - len(pivots), pivot_fn)В реальном мире Quickselect отлично себя проявляет: он почти не потребляет лишних ресурсов и выполняется в среднем за
O(n). Давайте докажем это.Доказательство среднего времени O(n)
В среднем pivot разбивает список на две приблизительно равных части. Поэтому каждая последующая рекурсия оперирует с 1?2 данных предыдущего шага.
Существует множество способов доказательства того, что этот ряд сходится к 2n. Вместо того, чтобы приводить их здесь, я сошлюсь на замечательную статью в Википедии, посвящённую этому бесконечному ряду.
Quickselect даёт нам линейную скорость, но только в среднем случае. Что, если нас не устраивает среднее, и мы хотим гарантированного выполнения алгоритма за линейное время?
Детерминированное O(n)
В предыдущем разделе я описал quickselect, алгоритм со средней скоростью
O(n). «Среднее» в этом контексте означает, что в среднем алгоритм будет выполняться за O(n). С технической точки зрения, нам может очень не повезти: на каждом шаге мы можем выбирать в качестве pivot наибольший элемент. На каждом этапе мы сможем избавляться от одного элемента из списка, и в результате получим скорость O(n^2), а не O(n).С учётом этого, нам нужен алгоритм для подбора опорных элементов. Нашей целью будет выбор за линейное время pivot, который в худшем случае удаляет достаточное количество элементов для обеспечения скорости
O(n) при использовании его вместе с quickselect. Этот алгоритм был разработан в 1973 году Блумом (Blum), Флойдом (Floyd), Праттом (Pratt), Ривестом (Rivest) и Тарьяном (Tarjan). Если моего объяснения вам не хватит, то можете изучить их статью 1973 года. Вместо того, чтобы описывать алгоритм, я подробно прокомментирую мою реализацию на Python:def pick_pivot(l):
"""
Выбираем хорошй pivot в списке чисел l
Этот алгоритм выполняется за время O(n).
"""
assert len(l) > 0
# Если элементов < 5, просто возвращаем медиану
if len(l) < 5:
# В этом случае мы возвращаемся к первой написанной нами функции медианы.
# Поскольку мы выполняем её только для списка из пяти или менее элементов, она не
# зависит от длины входных данных и может считаться постоянным
# временем.
return nlogn_median(l)
# Сначала разделим l на группы по 5 элементов. O(n)
chunks = chunked(l, 5)
# Для простоты мы можем отбросить все группы, которые не являются полными. O(n)
full_chunks = [chunk for chunk in chunks if len(chunk) == 5]
# Затем мы сортируем каждый фрагмент. Каждая группа имеет фиксированную длину, поэтому каждая сортировка
# занимает постоянное время. Поскольку у нас есть n/5 фрагментов, эта операция
# тоже O(n)
sorted_groups = [sorted(chunk) for chunk in full_chunks]
# Медиана каждого фрагмента имеет индекс 2
medians = [chunk[2] for chunk in sorted_groups]
# Возможно, я немного повторюсь, но я собираюсь доказать, что нахождение
# медианы списка можно произвести за доказуемое O(n).
# Мы находим медиану списка длиной n/5, поэтому эта операция также O(n)
# Мы передаём нашу текущую функцию pick_pivot в качестве создателя pivot алгоритму
# quickselect. O(n)
median_of_medians = quickselect_median(medians, pick_pivot)
return median_of_medians
def chunked(l, chunk_size):
"""Разделяем список `l` на фрагменты размером `chunk_size`."""
return [l[i:i + chunk_size] for i in range(0, len(l), chunk_size)]Давайте докажем, что медиана медиан является хорошим pivot. Нам поможет, если мы представим визуализацию нашего алгоритма выбора опорных элементов:
Красным овалом обозначены медианы фрагментов, а центральным кругом — медиана медиан. Не забывайте, мы хотим, чтобы pivot разделял список как можно ровнее. В худшем возможном случае каждый элемент в синем прямоугольнике (слева вверху) будет меньше или равен pivot. Верхний правый прямоугольник содержит 3?5 половины строк —
3/5*1/2=3/10. Поэтому на каждом этапе мы избавляемся по крайней мере от 30% строк.Но достаточно ли нам отбрасывать 30% элементов на каждом этапе? На каждом этапе наш алгоритм должен выполнять следующее:
- Выполнять работу O(n) по разбиению элементов
- Для рекурсии решать одну подзадачу размером в 7?10 от исходной
- Для вычисления медианы медиан решать одну подзадачу размером с 1?5 от исходной
В результате мы получаем следующее уравнение полного времени выполнения
T(n):Не так уж просто доказать, почему это равно
O(n). Быстрое решение заключается в том, чтобы положиться на основную теорему о рекуррентных соотношениях. Мы попадаем в третий случай теоремы, при котором работа на каждом уровне доминирует над работой подзадач. В этом случае общая работа будет просто равна работе на каждом уровне, то есть O(n).Подводим итог
У нас есть quickselect, алгоритм, который находит медиану за линейное время при условии наличия достаточно хорошей опорного элемента. У нас есть алгоритм медианы медиан, алгоритм
O(n) для выбора опорного элемента (который достаточно хорош для quickselect). Соединив их, мы получили алгоритм нахождения медианы (или n-ного элемента в списка) за линейное время!Медианы за линейное время на практике
В реальном мире почти всегда достаточно случайного выбора медианы. Хотя подход с медианой медиан всё равно выполняется за линейное время, на практике его вычисление длится слишком долго. В стандартной библиотеке
C++ используется алгоритм под названием introselect, в котором применено сочетание heapselect и quickselect; предел его выполнения O(n log n). Introselect позволяет использовать обычно быстрый алгоритм с плохим верхним пределом в сочетании с алгоритмом, который медленнее на практике, но имеет хороший верхний предел. Реализации начинают с быстрого алгоритма, но возвращаются к более медленному, если не могут выбрать эффективные опорные элементы.В завершение приведу сравнение элементов, используемых в каждой из реализаций. Это не скорость выполнения, а общее количество элементов, которые рассматривает функция quickselect. Здесь не учитывается работа по вычислению медианы медиан.
Именно этого мы и ожидали! Детерминированный опорный элемент почти всегда рассматривает при quickselect меньшее количество элементов, чем случайный. Иногда нам везёт и мы угадываем pivot с первой попытки, что проявляется как впадины на зелёной линии. Математика работает!
- Это может стать интересным применением поразрядной сортировки (radix sort), если вам нужно найти медиану в списке целых чисел, каждое из которых меньше 232.
- На самом деле в Python используется Timsort, впечатляющее сочетание теоретических пределов и практической скорости. Заметки о списках в Python.
Комментарии (32)
yorko
20.01.2018 15:11-1Кажется, что такой выбор опорного элемента (медиана из медиан) должен помогать и quicksort-у, поскольку при каждом разделении получаются более равные по длине 2 списка. Однако для quicksort доказано, что оптимальная стратегия выбора опорного элемента – случайная.
splav_asv
20.01.2018 16:54Можно ссылочку на доказательство?
yorko
20.01.2018 18:37Я это помню по курсу Тима Рафгардена. В видео "Choosing a good pivot" лектор про это говорит и доказывает, что так получится среднее время O(n logn), но формальное доказательство, что случайная стратегия – наилучшая, наверное в его книге надо искать.
splav_asv
20.01.2018 18:55Обычно везде говорится, что этастратегия, не имеет детерминированной «kill sequence» (данные при которых сваливаемся в O(n^2)). На разных данных оптимальный выбор будет отличаться. Если сортировать последовательность равномерно распределённых случайных величин(или близкую к ней), то вообще можно брать первый попавшийся pivot: среднее не пострадает, а константа лучше.
Это я всё к чему — оптимальной стратегии нет. Есть размен ухудшения константы в среднем на отсутствие детерминированного худшего случая. А дальше зависит от данных и требований.yorko
20.01.2018 19:11Пожалуй, про "оптимальную стратегию" – это громко сказано, запомнилось просто, как акцентировано это Рафгарден отмечал. В книге Сэджвика "Algorithms" про quicksort нашел только утверждение (с. 295, Proposition L), что в худшем случае quicksort делает ~ N^2/2 сравнений (очевидно), но случайное перемешивание на каждом шаге снижает число сравнений обратно до ~ n log n. Но да, это далеко не то же самое, о чем я говорил в начале.
chersanya
20.01.2018 21:16Нужно как минимум понимать, что при рассмотрении quicksort со случайным выбором элемента среднее берётся не по входным данным, а по реализациям случайных решений в ходе алгоритма. Это и значит, что не существует набора данных, дающих n^2 — среднее количество операций для любого массива будет O(n log n).
Можно даже посчитать вероятность того, что алгоритм будет выполняться например в 5 раз медленнее, чем в среднем. Эта вероятность получается весьма малой, особенно при больших n — поэтому на практике имея адекватный генератор случайных чисел нет смысла что-то оптимизировать по сравнению со случайным выбором элемента.splav_asv
20.01.2018 21:46Не вижу противоречия.
chersanya
20.01.2018 22:38Например:
Есть размен ухудшения константы в среднем на отсутствие детерминированного худшего случая.
«детерминированного худшего случая» нет ни в алгоритме со случайным выбором pivot, ни в описанном в статье.splav_asv
20.01.2018 22:46Ну так я и не говорю, что в этих вариантах есть худший случай. Случайный выбор можно сделать быстрее O(1), данный алгоритм O(n). Первый не влияет на константу по идее, второй влияет. Для первого в каждом конкретном случае есть плохой вариант. Для второго, исключая массив повторяющихся величин, нет.
В данном случае такой размен.chersanya
20.01.2018 23:18Что-то я вас вообще не понял.
Случайный выбор можно сделать быстрее O(1)
Алгоритмов быстрее О(1) не бывает.
Для первого в каждом конкретном случае есть плохой вариант.
Поясните, что имеете в виду под случаем и что под вариантом?splav_asv
21.01.2018 00:10можно сделать быстрее O(1)
Прошу прощения, пропущена запятая:
Случайный выбор можно сделать быстрее, за O(1), данный алгоритм(из статьи) O(n).
В данном контексте случай — реализация генератора псевдослучайных чисел. Вариант — набор данных, который нужно подать на вход quicksort (зная состояние ГПСЧ), чтобы свалиться в O(n^2).
MikailBag
20.01.2018 17:40Насколько я помню, quicksort используется исключительно из-за низкой константы, а при использовании "безопасного" поиска медианы он работает медленнее других алгоритмов сортировки.
chersanya
20.01.2018 18:06Причём тут случайная стратегия, если речь о детерминированных алгоритмах? Детерминированный quicksort стандартно пишется как раз с таким же выбором опорного элемента, как и в алгоритме поиска медианы. Единственное отличие quicksort от quickselect — то, что он обрабатывает обе части массива, а не только одну из них.
Соответственно случайный выбор элемента можно абсолютно так же использовать и в поиске медианы, если не ограничиваться детерминированными алгоритмами.yorko
20.01.2018 18:42Вывод по графику в статье
Детерминированный опорный элемент почти всегда рассматривает при quickselect меньшее количество элементов, чем случайный
Мой вопрос (или мысли вслух): работает ли то же самое для quicksort. Будет ли детерминированный quicksort with median pivots рассматривать меньше элементов, чем недетерминированный quicksort with random pivot.
chersanya
20.01.2018 21:12Да, это тоже верно. Однако, как отмечено и в статье, это не говорит о более быстром выполнении алгоритма. Всегда или почти всегда на практике быстрее работает случайный выбор.
nikitadanilov
20.01.2018 15:15А почему pivot выбирается именно как элемент массива с некоторым индексом? В случае quicksort это понятно, но для quickselect можно выбрать любое значение (не обязательно из массива). Например, брать среднее (O(n)), тогда массив будет гарантированно хорошо делиться.
andreybotanic
20.01.2018 15:51Вот контрпример: массив, состоящий из одного числа 100 и 99 нулей. Среднее значение равно 1. В итоге в одной части окажутся все нули, а во второй всего одно число — наша сотня. Не самой удачное разделение, не так ли?
netch80
22.01.2018 23:05Это только на один цикл вглубь. На следующем окажется, что все элементы одинаковы.
(Хотя, чтобы корректно выйти из этой ситуации, исходный алгоритм должен делить вход не на две части, а на три, или же пополамить набор элементов, равных опорному. Для quicksort некоторые методы деления отрезка специально заточены на оптимизацию такого варианта. Это можно перенести и на quickselect, но явно.)
nikitazvonarev
20.01.2018 21:36На самом деле, нас могут интересовать не только медианы, но и квартили и вообще произвольные квантили, короче говоря, k-е порядковые статистики (т.е. k-ый элемент в отсортированной последовательности). quickselect (который обрубленный quicksort) спокойно обобщается на этот случай с сохранением асимптотики среднего времени, а детерминированный подход можно ли обобщить?
В статье это никак не раскрыто.chersanya
20.01.2018 22:35Описанный в статье алгоритм это и есть детерминированный quickselect, и он абсолютно так же обобщается.
khim
21.01.2018 00:25Вы наверное какую-ту другую статью читали. Потому что то, что описано в статье сначала рассматривает нахождение медианы как нахождение k-й порядковой статистики с номером N/2, а потом, рекурсивно решает задачу нахождения k-й порядковой статистики.
Обратите внимание что один шаг алгоритма: он заключается в том, что мы отбрасываем либо 30% «слишком больших» чисел, либо 30% слишком маленьких, после чего ищем… медиану? Фиг вам: так как «перекосили» массив, то нам теперь нужна не медиана, нам нужен элемент с другим (впрочем легко высчитываемым) номером.
Куда вы это обобщать собрались?
snowSTAFF
20.01.2018 21:36Ваши алгоритмы. Все на Python 3.6.3
Что-то мне не везло. Это все на списке из 10 000 рандомных от 0 до 1000
Причем, 3-я функция работает с ошибкой.
Для 1-ой:
2.47 ms ± 39.9 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
Для 2-ой:
7.92 ms ± 159 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
Для 3-ей:
11.2 ms ± 73.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)alexeykuzmin0
20.01.2018 23:05Дык список слишком короткий, чтобы найти различие между 1 и 2. А 3 и должна работать дольше, она только сравнений меньше делает (ну и имеет теоретическую гарантию того, что не будет тупить — для варианта 2 можно подобрать контрпример).
sergio_nsk
20.01.2018 22:22В C++ уже всё сделано для вас:
#include <iostream> #include <vector> #include <algorithm> #include <functional> int main() { std::vector<int> v{5, 6, 4, 3, 2, 6, 7, 9, 3}; std::nth_element(v.begin(), v.begin() + v.size()/2, v.end()); std::cout << "The median is " << v[v.size()/2] << '\n'; }
Output:
The median is 5
Mingun
21.01.2018 11:55Вот здесь же:
# Для простоты мы можем отбросить все группы, которые не являются полными. O(n) full_chunks = [chunk for chunk in chunks if len(chunk) == 5]
можно упростить до:
last_chunk = chunks[len(chunks) - 1] full_chunks = chunks if len(last_chunk) < 5: full_chunks.pop()
и будет О(1). Разбивка ведь оригинального списка ведётся так, что меньше 5-ти элементов может быть только в последнем подсписке.
Кроме того, визуализация алгоритма неправильная же. В первом и последнем подсписках выбраны вовсе не медианы! В первом подсписке медиана 41, а не 99, а во втором — 65 вместо 116. Правда, это легко исправляется дописыванием единичек спереди к последним числам в подсписках, но всё равно странно приводить заведомо неправильную картинку в качестве визуализации алгоритма.
amakhrov
21.01.2018 22:24можно упростить
и будет О(1)Хорошее замечание. Стоит уточнить, что хотя асимптотика всего алгоритма от этого не меняется, мы уменьшаем константу.
rcoh
23.01.2018 10:41+1Привет, я автор. Извините за автоматический перевод. Спасибо, что заметили эту ошибку! Это был ужасный надзор, и я исправил его.
mikhaelkh
Рекуррентное соотношение должно быть T(n) <= T(0.2 * n) + T(0.7 * n) + O(n), нам важно, что (0.2 * n + 0.7 * n) / n = 0.9 < 1. Использовать то, что T(n) = O(n) при доказательстве этого же факта — это facepalm.
На самом деле очень просто — подберём константу c такую, что для всех n выполняется T(n) <= T(0.2 * n) + T(0.7 * n) + c * n или T(n) <= 10 * c * n , тогда T(n) <= 10 * c * n по индукции.