Нетривиальный алгоритм вычисления медианы скользящего окна

25

Мне нужно рассчитать бегущую медиану:

  • Ввод: , , вектор .k ( x 1 , x 2 , , x n )nk(x1,x2,,xn)

  • Вывод: vector , где - это медиана .y i ( x i , x i + 1 , , x i + k - 1 )(y1,y2,,ynk+1)yi(xi,xi+1,,xi+k1)

(Нет мошенничества с приближениями; я хотел бы иметь точные решения. Элементы являются большими целыми числами.)xi

Существует тривиальный алгоритм, который поддерживает дерево поиска размера ; общее время работы . (Здесь «дерево поиска» относится к некоторой эффективной структуре данных, которая поддерживает вставки, удаления и срединные запросы в логарифмическом времени.)O ( n log k )kO(nlogk)

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

Можем ли мы сделать что-нибудь существенно лучше?

Существуют ли нижние оценки (например, является ли тривиальный алгоритм асимптотически оптимальным для модели сравнения)?


Изменить: Дэвид Эппштейн дал хороший нижний предел для модели сравнения! Интересно, возможно ли все же сделать что-то немного более умное, чем тривиальный алгоритм?

Например, можем ли мы сделать что-то в этом духе: разделить входной вектор на части размера ; сортировать каждую часть (отслеживая исходные позиции каждого элемента); а затем использовать кусочно отсортированный вектор, чтобы эффективно найти текущие медианы без каких-либо вспомогательных структур данных? Конечно, это все равно будет , но на практике сортировка массивов, как правило, происходит намного быстрее, чем поддержка деревьев поиска.O ( n log k )kO(nlogk)


Изменить 2: Saeed хотел увидеть некоторые причины, почему я думаю, сортировка быстрее, чем операции дерева поиска. Вот очень быстрые тесты для , : n = 10 8k=107n=108

  • ≈ 8 с: сортировка векторов с элементами каждыйкn/kk
  • ≈ 10 с: сортировка вектора с элементамиn
  • ≈ 80 с: вставок и удалений в хеш-таблице размеракnk
  • ≈ 390 с: вставок и удалений в сбалансированном дереве поиска размеракnk

Хеш-таблица существует только для сравнения; это не имеет прямого использования в этом приложении.

Таким образом, мы имеем почти 50-кратное различие в производительности сортировки и сбалансированных операциях дерева поиска. И все станет намного хуже, если мы увеличим .k

(Технические детали: Данные = случайные 32-разрядные целые числа. Компьютер = типичный современный ноутбук. Тестовый код был написан на C ++ с использованием стандартных библиотечных процедур (std :: sort) и структур данных (std :: multiset, std :: unsorted_multiset). Я использовал два разных компилятора C ++ (GCC и Clang) и две разные реализации стандартной библиотеки (libstdc ++ и libc ++). Традиционно std :: multiset был реализован как высокооптимизированное красно-черное дерево.)

Юкка Суомела
источник
1
Я не думаю , что вы будете в состоянии улучшить . Причина в том, если вы посмотрите на окно х т , . , , , x t + k - 1 , вы никогда не можете исключить ни одно из чисел x t + knlogkxt,...,xt+k1из медианы будущего окна. Это означаетчто в любое время вы должны держатькрайней мерекxt+k2,...,xt+k1 целых числа в структуре данных, и, похоже, они не обновляются меньше, чем за время регистрации. k2
РБ
Мне кажется, что ваш тривиальный алгоритм не O ( n log k ) , я что-то не так понял? И я думаю, что из-за этого у вас есть проблема с большим k , иначе логарифмический фактор ничего не значит в практических приложениях, также нет большой скрытой константы в этом алгоритме. O((nk)klogk)O(nlogk)k
Саид
@Saeed: в тривиальном алгоритме вы обрабатываете элементы один за другим; на шаге вы добавляете x i в дерево поиска и (если i > k ) вы также удаляете x i - k из дерева поиска. Это n шагов, каждый из которых занимает O ( log k ) времени. ixii>kxiknO(logk)
Юкка Суомела
То есть вы имеете в виду сбалансированное дерево поиска, а не случайное дерево поиска?
Саид
1
@Saeed: Обратите внимание, что в моих тестах я даже не пытался найти медианы. Я только что сделал вставок и n удалений в дереве поиска размера k , и эти операции гарантированно займут O ( log k ) времени. Вам просто нужно признать, что операции с деревом поиска очень медленны на практике по сравнению с сортировкой. Это легко увидеть, если вы попытаетесь написать алгоритм сортировки, который работает путем добавления элементов в сбалансированное дерево поиска - он, безусловно, работает за O ( n log n ) , но на практике он будет смехотворно медленным, а также тратит много времени. памяти. nnkO(logk)O(nlogn)
Юкка Суомела

Ответы:

32

Вот нижняя граница от сортировки. Учитывая, что входной набор длины n должен быть отсортирован, создайте вход для вашей текущей задачи медианы, состоящей из n - 1 копий числа, меньшего, чем минимум S , затем самого S , затем n - 1 копий числа, большего, чем максимум S , и установите k = 2 n - 1 . Ходовые медианы этого входа такие же , как в отсортированном порядке S .Snn1SSn1Sk=2n1S

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

Дэвид Эппштейн
источник
6
Этот ответ действительно заставляет меня задуматься, верно ли и обратное: получим ли мы эффективный алгоритм сортировки, получим ли мы эффективный алгоритм работы медианы? (Например, подразумевается ли в алгоритме эффективной целочисленной сортировки эффективный алгоритм выполнения медианы для целых чисел? Или алгоритм IO-эффективной сортировки обеспечивает алгоритм медианной обработки, эффективный в IO?)
Юкка Суомела
1
Еще раз, большое спасибо за ваш ответ, это действительно поставило меня на правильный путь и дало вдохновение для алгоритма медианного фильтра на основе сортировки! В конце концов, мне удалось найти статью 1991 года, в которой приводились те же аргументы, что и здесь, а Пэт Морин дал указатель на другую соответствующую статью 2005 года; см. ссылки [6] и [9] здесь .
Юкка Суомела
9

Изменить: Этот алгоритм теперь представлен здесь: http://arxiv.org/abs/1406.1717


Да, для решения этой проблемы достаточно выполнить следующие операции:

  • Сортировать векторов, каждый из которых содержит k элементов.n/kk
  • Делать линейную обработку времени.

Очень грубо, идея заключается в следующем:

  • Рассмотрим два смежных блока ввода, и b , оба с k элементами; пусть элементы будут 1 , 2 , . , , , К и б 1 , б 2 , . , , , b k в порядке появления во входном векторе x .abka1,a2,...,akb1,b2,...,bkx
  • Сортируйте эти блоки и изучите ранг каждого элемента в блоке.
  • Дополните векторы и b указателями предшественника / преемника, чтобы, следуя цепочкам указателей, мы могли проходить элементы в возрастающем порядке. Таким образом, мы построили двусвязные списки a и b .abab
  • Один за другим, удалить все элементы из связанного списка , в обратном порядке появления б к , б к - 1 , . , , , б 1 . Всякий раз, когда мы удаляем элемент, помните, каким был его преемник и предшественник на момент удаления .bbk,bk1,...,b1
  • Теперь сохраните «срединные указатели» и q, которые указывают на списки a и b соответственно. Initialise р к средней точке ' , и Initialise ц к хвосту пустого списка Ь ' .pqabpaqb
  • Для каждого :i

    • Удаление из списка а ' (это O ( 1 ) время, просто удалите из связанного списка). Сравните a i с элементом, на который указывает p, чтобы увидеть, удалили ли мы до или после p .aiaO(1)aipp
    • Поместите обратно в список b в его исходное положение (это O ( 1 ) раз, мы запомнили предшественника и преемника b i ). Сравните b i с элементом, на который указывает q, чтобы увидеть, добавили ли мы элемент до или после q .bibO(1)bibiqq
    • Обновите указатели и q таким образом, чтобы медиана объединенного списка a b находилась либо в p, либо в q . (Это O ( 1 ) раз, просто следуйте связанным спискам один или два шага, чтобы все исправить. Мы будем отслеживать, сколько элементов находится до / после p и q в каждом списке, и мы будем поддерживать инвариант, что оба p и q указывают на элементы, которые находятся максимально близко к медиане.)pqabpqO(1)pqpq

Связанные списки - это просто массивы -элементных индексов, поэтому они легковесны (за исключением того, что локальность доступа к памяти плохая).k


Вот пример реализации и тесты:

Вот график времени работы (для ):n2106

  • Синий = сортировка + постобработка, .O(nlogk)
  • Зеленый = поддерживать две кучи, , реализация с https://github.com/craffel/median-filterO(nlogk)
  • Красный = поддерживать два дерева поиска, .O(nlogk)
  • Черный = сохранить отсортированный вектор, .O(nk)
  • k/2
  • Ось Y = время работы в секундах.
  • Данные = 32-разрядные целые и случайные 64-разрядные целые числа из различных распределений.

время работы

Юкка Суомела
источник
3

mO(nlogm+mlogn)

O(logm)O(logn)O(logn) заряд происходит только один раз за медиану.

O(nlogm+mlogk)

Джеффри Ирвинг
источник
К сожалению, это не работает так, как написано, поскольку, если вы не удалите элементы, счетчики не будут отражать новое окно. Я не уверен, что это можно исправить, но я оставлю ответ, если есть способ.
Джеффри Ирвинг
O(nlogm)
примечание стороны: Вопрос не ясен, подчиненная структура данных не определена, мы просто знаем что-то очень расплывчатое. Как вы хотите улучшить то, что вы не знаете, что это такое? как вы хотите сравнить свой подход?
Саид
1
Я извиняюсь за незавершенную работу. Я задал конкретный вопрос, необходимый для исправления этого ответа здесь: cstheory.stackexchange.com/questions/21778/… . Если вы считаете, что это уместно, я могу удалить этот ответ, пока не будет решен дополнительный вопрос.
Джеффри Ирвинг,