Selection in Worst-Case Linear Time: Median of Medians Algorithm Explained

Selection in Worst-Case Linear Time: Median of Medians Algorithm

I remember being genuinely surprised the first time I learned that finding the $i$-th smallest element of an unsorted array doesn’t actually require sorting the whole array first. Sorting would cost $O(n \lg n)$, but the selection problem — finding just one specific order statistic, like the median — can be solved in $O(n)$ time, even in the worst case. The algorithm that achieves this, commonly called “median of medians” or the “BFPRT algorithm” after its five inventors, is one of the cleverest pieces of algorithm design I’ve come across, precisely because of how it uses a recursive, self-referential trick to guarantee a good pivot choice every single time.

I want to walk through this algorithm carefully because it demonstrates a powerful and somewhat unusual technique: using an approximate, recursively-computed estimate (the median of medians) as a pivot, specifically engineered so that the resulting partition is guaranteed to be reasonably balanced, no matter what the input looks like.

History and Background

The median of medians algorithm was introduced in a landmark 1973 paper titled “Time Bounds for Selection,” authored by Manuel Blum, Robert W. Floyd, Vaughan Pratt, Ronald L. Rivest, and Robert Tarjan — five of the most influential names in theoretical computer science, which is why the algorithm is sometimes referred to by the acronym BFPRT, formed from their initials. This paper was groundbreaking because it was the first to demonstrate a worst-case linear-time algorithm for the general selection problem, settling an important open question about whether selection could be done as efficiently as finding the minimum or maximum (which trivially take $O(n)$ time) even for arbitrary order statistics like the median.

The algorithm builds on the earlier randomized “quickselect” approach, which Tony Hoare had introduced in 1961 (alongside quicksort) and which achieves expected $O(n)$ time but has a worst-case of $O(n^2)$ for adversarial inputs. The 1973 BFPRT paper’s key innovation was finding a deterministic way to always choose a “good enough” pivot, eliminating the dependence on randomness or luck entirely.

Problem Statement

Given an unsorted array of $n$ distinct numbers and an integer $i$ (with $1 \le i \le n$), I want to find the element that would be in position $i$ if the array were sorted — the $i$-th order statistic — using as little time as possible, and critically, I want a guarantee on the worst-case running time, not just the average case. Sorting the whole array and then indexing into position $i$ works, but costs $O(n \lg n)$; I want to do better, achieving $O(n)$ time even in the worst case.

Core Concepts

How It Works

The algorithm, SELECT(A, i), proceeds as follows:

  1. Divide the $n$ elements of the array into $\lceil n/5 \rceil$ groups, each containing 5 elements (except possibly the last group, which may have fewer).
  2. Find the median of each group by simply sorting each small group (a constant-size sort, since each group has at most 5 elements) and picking the middle element — this takes $O(1)$ time per group, so $O(n)$ time total across all $\lceil n/5 \rceil$ groups.
  3. Recursively find the median of these medians, calling SELECT recursively on this set of $\lceil n/5 \rceil$ group-medians to find their median — call this value $x$.
  4. Partition the original array around $x$ as the pivot, exactly as in quicksort, producing a position $k$ such that everything before position $k$ is less than $x$, and everything after is greater.
  5. Recurse or return: If $i == k$, return $x$. If $i < k$, recursively call SELECT on the left part of the partition looking for the $i$-th element there. If $i > k$, recursively call SELECT on the right part, looking for the $(i – k)$-th element there (since the first $k$ elements are excluded).

Working Principle

The reason this works — and works with a guaranteed worst-case linear time bound — comes down to a clever guarantee about how good a pivot the median of medians actually is. Because $x$ is chosen as the median of the group medians, I can prove that at least roughly 3/10 of all the elements in the array are guaranteed to be less than or equal to $x$, and at least roughly 3/10 are guaranteed to be greater than or equal to $x$. This means the partition step is guaranteed to be reasonably balanced — never as bad as, say, a 1-versus-(n-1) split — no matter what the input array looks like, since this bound doesn’t depend on the input being random or “nice” in any way.

This guaranteed balance is exactly what a randomized quickselect hopes for but can’t guarantee — random pivot selection could, in the worst case (even if unlikely), keep picking terrible pivots repeatedly. The median-of-medians approach removes luck from the equation entirely by paying a modest additional cost (the recursive computation of the pivot itself) in exchange for a hard worst-case guarantee.

Mathematical Foundation

Guaranteed fraction of elements bounded by the pivot. Consider the $\lceil n/5 \rceil$ group medians, and let $x$ be their median (found recursively). At least half of the $\lceil n/5 \rceil$ groups have their median $\le x$, which is:

$$ \frac{1}{2}\left\lceil \frac{n}{5} \right\rceil $$

groups. For each such group (excluding possibly the last, partial group, and the group containing $x$ itself), at least 3 of its 5 elements are $\le$ that group’s median, and therefore $\le x$. This gives a lower bound on the number of elements guaranteed to be $\le x$:

$$ 3 \left( \frac{1}{2}\left\lceil \frac{n}{5} \right\rceil – 2 \right) \ge \frac{3n}{10} – 6 $$

By a symmetric argument, at least $\frac{3n}{10} – 6$ elements are guaranteed to be $\ge x$. This means the partition step, in the worst case, calls SELECT recursively on a subproblem of size at most:

$$ \frac{7n}{10} + 6 $$

The recurrence. Combining the cost of finding group medians ($O(n)$), the recursive call to find the median of medians (on $\lceil n/5 \rceil$ elements), the partition step ($O(n)$), and the recursive call on the resulting subproblem (at most $\frac{7n}{10} + 6$ elements), I get the recurrence:

$$ T(n) \le T\left(\frac{n}{5}\right) + T\left(\frac{7n}{10} + 6\right) + O(n) $$

Solving the recurrence by substitution. I guess $T(n) \le cn$ for a suitable constant $c$ and sufficiently large $n$, and verify inductively:

$$ T(n) \le c\frac{n}{5} + c\left(\frac{7n}{10} + 6\right) + an $$

$$ = \frac{cn}{5} + \frac{7cn}{10} + 6c + an $$

$$ = \frac{9cn}{10} + 6c + an $$

$$ = cn – \left(\frac{cn}{10} – 6c – an\right) $$

For this to be $\le cn$, I need $\frac{cn}{10} – 6c – an \ge 0$, which holds for $n$ sufficiently large as long as $c$ is chosen large enough relative to $a$ (specifically, choosing $c$ large enough that $\frac{c}{10} > a$ leaves room to absorb the $6c$ term for large $n$). This confirms:

$$ T(n) = O(n) $$

The specific choice of group size 5 is what makes this work: the coefficient $\frac{1}{5} + \frac{7}{10} = \frac{9}{10} < 1$ is strictly less than 1, which is exactly what’s needed for the recurrence to resolve to a linear bound. If groups of size 3 were used instead, the analogous coefficient would be $\ge 1$, and the recurrence would fail to resolve to $O(n)$.

Diagrams

flowchart TD
    A[Array of n elements] --> B[Divide into groups of 5]
    B --> C[Sort each group, find its median: O n total]
    C --> D[Recursively find median of the ceil n/5 medians]
    D --> E[Use that value x as the pivot]
    E --> F[Partition array around x: O n]
    F --> G{Compare i to pivot position k}
    G -->|i == k| H[Return x]
    G -->|i < k| I[Recurse SELECT on left part]
    G -->|i > k| J[Recurse SELECT on right part, looking for i - k]

Pseudocode

SELECT(A, p, r, i)
    // Base case: small array, just sort and index
    if r - p + 1 <= 5
        sort A[p..r]
        return A[p + i - 1]

    // Step 1 & 2: divide into groups of 5, find each group's median
    numGroups = ceil((r - p + 1) / 5)
    let M[1..numGroups] be a new array
    for j = 1 to numGroups
        groupStart = p + (j - 1) * 5
        groupEnd = min(groupStart + 4, r)
        sort A[groupStart..groupEnd]
        M[j] = A[groupStart + floor((groupEnd - groupStart) / 2)]

    // Step 3: recursively find the median of medians
    x = SELECT(M, 1, numGroups, ceil(numGroups / 2))

    // Step 4: partition A[p..r] around x
    k = PARTITION-AROUND-VALUE(A, p, r, x)

    // Step 5: recurse on the correct side
    rank = k - p + 1
    if i == rank
        return A[k]
    elseif i < rank
        return SELECT(A, p, k - 1, i)
    else
        return SELECT(A, k + 1, r, i - rank)

Step-by-Step Example

Let me trace through finding the median (6th smallest, since $i = 6$ for $n = 11$) of the array:

A = [12, 3, 5, 7, 4, 19, 26, 21, 2, 6, 8]

Step 1 & 2: Divide into groups of 5: [12, 3, 5, 7, 4], [19, 26, 21, 2, 6], [8]. Sorting each group and finding medians:

So M = [5, 19, 8].

Step 3: Recursively find the median of M = [5, 19, 8] (looking for the 2nd smallest of 3 elements, i.e., $\lceil 3/2 \rceil = 2$). Since this is small enough for the base case, sort it: [5, 8, 19], and the 2nd element is 8. So x = 8.

Step 4: Partition the original array A around x = 8. Elements less than 8: {3, 5, 7, 4, 2, 6} (6 elements). Elements greater than 8: {12, 19, 26, 21} (4 elements). After partitioning, 8 lands at position 7 (1-indexed), meaning k‘s rank is 7.

Step 5: I want $i = 6$. Since $6 < 7$ (the rank of the pivot), I recurse into the left part, now looking for the 6th smallest among the 6 elements less than 8: {3, 5, 7, 4, 2, 6}. But there are only 6 elements there and I want the 6th (the largest of that subset), which resolves quickly: sorting {2, 3, 4, 5, 6, 7} gives the 6th element as 7.

So the overall answer is 7, meaning 7 is the 6th smallest (median) of the original 11-element array. I can double check: fully sorting A gives [2, 3, 4, 5, 6, 7, 8, 12, 19, 21, 26], and the 6th element is indeed 7. The algorithm gives the correct answer.

Time Complexity

This is in contrast to randomized quickselect, whose expected time is $O(n)$ but whose worst-case time is $O(n^2)$ for adversarially chosen (or extremely unlucky) inputs.

Space Complexity

The algorithm uses $O(n)$ additional space in typical implementations, primarily for storing the array of group medians at each level of recursion (an array of size $\lceil n/5 \rceil$), plus $O(\lg n)$ stack space for the recursion depth itself, since each recursive call operates on a fraction of the previous problem size (at most $\frac{7n}{10} + 6$), giving logarithmic recursion depth. In-place variants exist that reduce auxiliary space by reusing the original array’s storage for group medians, bringing the extra space down closer to $O(\lg n)$ (dominated by recursion stack depth) at the cost of more complex implementation.

Correctness Analysis

Correctness follows from two separate guarantees working together. First, the partition correctness is identical to quicksort’s partition logic: after partitioning around any pivot value $x$, every element before the pivot’s final position is guaranteed to be $\le x$, and every element after is guaranteed to be $\ge x$ — this is true regardless of how $x$ was chosen. Second, the recursive structure’s correctness follows from induction on array size: the base case (arrays of size $\le 5$) is trivially correct since it directly sorts and indexes. For larger arrays, assuming the recursive calls to find the median of medians and to search the appropriate partition side are both correct (by the inductive hypothesis, since both operate on strictly smaller subproblems), the top-level call correctly narrows down to and returns the true $i$-th order statistic.

The worst-case linear time bound doesn’t affect correctness — it’s a separate guarantee about efficiency — but it’s what distinguishes this algorithm from being merely correct (which quickselect also is) to being correct and asymptotically optimal in the worst case.

Advantages

Disadvantages

Applications

Implementation in C

#include <stdio.h>
#include <stdlib.h>

void insertionSort(int arr[], int left, int right) {
    for (int i = left + 1; i <= right; i++) {
        int key = arr[i];
        int j = i - 1;
        while (j >= left && arr[j] > key) {
            arr[j + 1] = arr[j];
            j--;
        }
        arr[j + 1] = key;
    }
}

void swap(int* a, int* b) {
    int temp = *a;
    *a = *b;
    *b = temp;
}

/* Partitions arr[left..right] around the given pivot value.
   Returns the final index of the pivot after partitioning. */
int partitionAroundValue(int arr[], int left, int right, int pivotValue) {
    int pivotIndex = left;
    for (int i = left; i <= right; i++) {
        if (arr[i] == pivotValue) {
            pivotIndex = i;
            break;
        }
    }
    swap(&arr[pivotIndex], &arr[right]);  /* move pivot to the end temporarily */

    int storeIndex = left;
    for (int i = left; i < right; i++) {
        if (arr[i] < pivotValue) {
            swap(&arr[i], &arr[storeIndex]);
            storeIndex++;
        }
    }
    swap(&arr[storeIndex], &arr[right]);  /* move pivot to its final position */
    return storeIndex;
}

/* Median of medians selection: finds the element that would be at
   position i (0-indexed) if arr[left..right] were fully sorted. */
int select_(int arr[], int left, int right, int i) {
    if (right - left + 1 <= 5) {
        insertionSort(arr, left, right);
        return arr[left + i];
    }

    int n = right - left + 1;
    int numGroups = (n + 4) / 5;
    int* medians = (int*)malloc(numGroups * sizeof(int));

    for (int g = 0; g < numGroups; g++) {
        int groupStart = left + g * 5;
        int groupEnd = groupStart + 4;
        if (groupEnd > right) groupEnd = right;
        insertionSort(arr, groupStart, groupEnd);
        medians[g] = arr[groupStart + (groupEnd - groupStart) / 2];
    }

    /* Recursively find the median of the group medians */
    int medianOfMedians = select_(medians, 0, numGroups - 1, (numGroups - 1) / 2);
    free(medians);

    int pivotFinalIndex = partitionAroundValue(arr, left, right, medianOfMedians);
    int rank = pivotFinalIndex - left;   /* 0-indexed rank within arr[left..right] */

    if (i == rank) {
        return arr[pivotFinalIndex];
    } else if (i < rank) {
        return select_(arr, left, pivotFinalIndex - 1, i);
    } else {
        return select_(arr, pivotFinalIndex + 1, right, i - rank - 1);
    }
}

int main(void) {
    int arr[] = {12, 3, 5, 7, 4, 19, 26, 21, 2, 6, 8};
    int n = sizeof(arr) / sizeof(arr[0]);

    int i = 5;  /* 0-indexed: 5 means the 6th smallest element (the median) */
    int result = select_(arr, 0, n - 1, i);

    printf("The %d-th smallest element (0-indexed) is: %d\n", i, result);

    return 0;
}

I chose to represent i as a 0-indexed position internally (so i = 5 corresponds to the 6th smallest element), which is a common convention in C code, though I explicitly note the distinction since the pseudocode above used 1-indexed rank for closer alignment with the textbook presentation.

Sample Input and Output

Input:
arr = [12, 3, 5, 7, 4, 19, 26, 21, 2, 6, 8]
i = 5 (0-indexed, meaning the 6th smallest element)

Output:
The 5-th smallest element (0-indexed) is: 7

This matches my hand-traced example exactly, confirming that 7 is indeed the median (6th smallest of 11 elements) of the original array.

Optimization Techniques

Common Mistakes

Further Reading

Exit mobile version