Wyszukiwanie mediany

WPROWADZENIE

Mediana zbioru jest środkowym elementem tego zbioru po jego uporządkowaniu. Jeżeli zbiór zawiera nieparzystą liczbę n elementów, to ma jeden element środkowy na pozycji i=(n+1)//2 i jest on medianą tego zbioru. Jeżeli zbiór zawiera parzystę liczbę n elementów, to zbiór ma dwie mediany na pozycjach i=n//2 (mediana dolna) oraz i=n//2+1 (mediana górna). W praktyce za medianę przyjmuje się na ogół medianę dolną.

Prosta metoda znajdowania mediany polega na posortowaniu zbioru i wybraniu elementu środkowego. Lepszy algorytm jest opisany w [BK1982].

Nie wiadomo dokładnie, ile porównań jest potrzebnych do znalezienia mediany. W pracy [BFPRT1973] podano dolne ograniczenie równe 3n//2-2 porównań. Bent i John [BJ1985] podali dolną granicę 2n porównań. Górną granicę 3n podali Schonhage, Paterson i Pippenger [SPP1976].

Wybrane algorytmy znajdowania mediany.

MEDIANA KROCZĄCA

Dla pewnego strumienia danych liczbowych możemy chcieć na bieżąco raportować element środkowy z ostatnich n elementów ze strumienia. Jest to mediana krocząca, którą chcemy obliczać w czasie O(log n), gdzie n jest szerokością okna danych. Zalecaną strukturą danych w tym przypadku jest indeksowana lista z przeskokami.

https://en.wikipedia.org/wiki/Skip_list

ALGORYTM HOARE'A

Algorytm jest oparty na tym samym pomyśle, co algorytm sortowania quicksort. Wersja wg Cormena.


def swap(L, left, right):
    """Zamiana miejscami elementów."""
    tmp = L[left]
    L[left] = L[right]
    L[right] = tmp

def partition(L, left, right):
    """Zwraca indeks elementu rozdzielającego."""
    # Element rozdzielający to ostatni z prawej,
    # dlatego na końcu trzeba go przerzucić do środka.
    # Będzie on na docelowej pozycji ze względu na sortowanie.
    x = L[right]   # element rozdzielający
    i = left
    for j in range(left, right):
        if L[j] <= x:
            swap(L, i, j)
            i += 1
    swap(L, i, right)
    return i

def quicksort(L, left, right):
    """Sortowanie szybkie wg Cormena str. 169."""
    if left >= right:
        return
    pivot = partition(L, left, right)
    # pivot jest na swoim miejscu
    quicksort(L, left, pivot - 1)
    quicksort(L, pivot + 1, right)

def select(L, left, right, k):
    """Selekcja na bazie Cormena."""
    if (right-left+1) < k:
        raise ValueError("k is too large")
    if left == right and k == 1:
        return L[left]
    pivot = partition(L, left, right)
    # pivot jest na swoim miejscu
    k2 = pivot - left + 1
    if k == k2:
        return L[pivot]
    elif k < k2:
        return select(L, left, pivot - 1, k)
    else:
        return select(L, pivot + 1, right, k - k2)

Algorytm selekcji Hoare'a ma pesymistyczną złożoność O(n^2), tak jak quicksort. Ale średnia liczba porównań nie przekracza 4n, co daje średnią złożoność O(n). Zauważmy, że w algorytmie selekcji tylko jedno zadanie jest wykonywane rekurencyjnie, a w quicksort są to dwa zadania.

ALGORYTM MAGICZNYCH PIĄTEK

Algorytm stosuje technikę dziel i zwyciężaj. Kluczowy jest dobry wybór punktu podziału, przy którym stała część elementów będzie poniżej i powyżej niego. Takim dobrym punktem podziału jest mediana zbioru.


def select_five(L, left, right, k):
    # 1. Jeżeli jest mało elementów, to sortuj.
    # p powinno być co najmniej 5, aby ustawiło się i > 0.
    p = 5   # może być kilkadziesiąt
    if (right-left+1) < p:
        insertsort(L, left, right)
        return left+k-1   # zwracam indeks
    # 2. Podziel listę na 5-elementowe podzbiory, najwyżej jeden 4-elementowy.
    # 3. Posortuj podzbiory.
    left2 = left
    right2 = left + 4
    i = left   # pierwszy wolny
    while right2 <= right:
        insertsort(L, left2, right2)
        # Przerzucamy mediany na początek tablicy.
        swap(L, i, left2+2)
        i += 1
        left2 += 5
        right2 += 5
    # Tu można posortować zbiory mniej niż 5-elementowe.
    if right2 == right+1 or right2 == right+2:
        insertsort(L, left2, right)
        swap(L, i, left2+1)
        i += 1
    # 5. Wyznaczamy medianę median rekurencyjnie.
    median_idx = select_five(L, left, i-1, (i-left+1) // 2)
    # 6. Dalej jak Hoare, mediana będzie punktem podziału.
    swap(L, median_idx, right)   # element podzialu na prawo
    pivot = partition(L, left, right)
    k2 = pivot - left + 1
    if k == k2:
        return pivot   # zwracam indeks
    elif k < k2:
        return select_five(L, left, pivot-1, k)
    else:
        return select_five(L, pivot+1, right, k-k2)

Niech T(n) oznacza złożoność czasową algorytmu. Wykonanie algorytmu składa się z trzech etapów.
(1) Znajdowanie median piątek w czasie O(n).
(2) Znajdowanie rekurencyjne mediany median w czasie T(n/5).
(3) Wykonania rekurencyjnego w czasie T(3n/4).
To ostatnie oszacowanie wynika z faktu, że przynajmniej 1/4 elementów jest nie większa od mediany median M (lewa górna ćwiartka), a więc 3/4 elementów jest większa lub równa M. Podobnie przynajmniej 1/4 elementów jest nie mniejsza od M (prawa dolna ćwiartka), czyli 3/4 elementów jest mniejsza lub równa M. Wygodnie jest to pokazać na rysunku.


. elementy mniejsze lub równe M (mediana median)
. [nie mniej niż n/4]]
. +-----------------------+
. |(*)  (*)  (*)  (*)  (*)| (*)  (*)  (*)  (*)
. ||/\  |/\  |/\  |/\  |/\| |/\  |/\  |/\  |/\
. |(*)  (*)  (*)  (*)  (*)| (*)  (*)  (*)  (*)
. ||/\  |/\  |/\  |/\  |/\| |/\  |/\  |/\  |/\
. |                   +---|-------------------+
. |(*)<=(*)<=(*)<=(*)<=(M)<=(*)<=(*)<=(*)<=(*)|
. +-------------------|---+                   |
.  |/\  |/\  |/\  |/\ ||/\  |/\  |/\  |/\  |/\|
.  (*)  (*)  (*)  (*) |(*)  (*)  (*)  (*)  (*)|
.  |/\  |/\  |/\  |/\ ||/\  |/\  |/\  |/\  |/\| elementy większe lub równe M
.  (*)  (*)  (*)  (*) |(*)  (*)  (*)  (*)  (*)| [nie mniej niż n/4]
.                     +-----------------------+

Całkowite oszacowanie ma postać T(n) <= O(n) + T(n/5) + T(3n/4), co daje rozwiązanie postaci T(n) = O(n). Kluczowa jest nierówność 1/5 + 3/4 = 19/20 < 1.

Zamiast podziału na zbiory 5-elementowe można zrobić podziały na zbiory 7-elementowe, co nie zmieni idei algorytmu. Istotne jest uzyskanie oszacowania T(n) <= O(n) + T(an) + T(bn), gdzie a+b < 1, co daje całkowity czas T(n) = O(n/(1-a-b)). Na koniec zauważmy, że w typowych przypadkach szybszy jest quickselect ze względu na prostszy kod.

REFERENCJE