x. Range Order Statistics and Wavelets

Queries like “$k$-th smallest in $a[l..r)$”, “how many values in $[x,y)$ lie in $a[l..r)$”, predecessor/successor within a range. Static and online: wavelet matrix or persistent segment tree, $O(\log n)$. Static and offline: sweep with a Fenwick tree. With point updates: Fenwick of dynamic segment trees, or offline parallel binary search.

Methods

Method	Query	Memory	Updates	Online
sort the range (naive)	$O(n\log n)$	$O(1)$	yes	yes
wavelet matrix	$O(\log \sigma )$	$n\log \sigma$ bits	no	yes
persistent segment tree over values	$O(\log n)$	$O(n\log n)$ ints	no	yes (PST)
merge sort tree	count $O({\log }^{2}n)$, $k$-th $O({\log }^{3}n)$	$O(n\log n)$	no	yes
offline sweep by value + Fenwick	count $O(\log n)$	$O(n)$	yes	no
overall binary search (parallel) + Fenwick	$k$-th $O({\log }^{2}n)$ amortized	$O(n)$	yes	no (PBS)
Fenwick of dynamic value trees (“tree in tree”)	$O({\log }^{2}n)$	$O(n{\log }^{2}n)$	yes	yes
sqrt decomposition with sorted blocks	$O(\sqrt{n\log n})$ per count	$O(n)$	yes	yes
Mo’s algorithm + a value-indexed Fenwick	$O(\sqrt{n}\log n)$ per count	$O(n)$	no	no (offline)

Wavelet matrix

Compress values to $[0,\sigma )$. At level $h$ (from the top bit down), record bit $h$ of every value in the current order, then stably move zeros before ones for the next level. A range $[l,r)$ at level $h$ maps to $[{\mathrm{rank}}_0(l),{\mathrm{rank}}_0(r))$ among zeros and to $[Z_h+{\mathrm{rank}}_1(l),Z_h+{\mathrm{rank}}_1(r))$ among ones, where $Z_h$ is the number of zeros on that level.

struct BitVec {
    vector<ull> b; vector<int> pre;
    BitVec(int n = 0) : b(n / 64 + 1), pre(n / 64 + 1) {}
    void set(int i) { b[i >> 6] |= 1ULL << (i & 63); }
    void build() { for (int j = 1; j < sz(b); j++) pre[j] = pre[j - 1] + __builtin_popcountll(b[j - 1]); }
    int rank1(int i) const { return pre[i >> 6] + __builtin_popcountll(b[i >> 6] & ((1ULL << (i & 63)) - 1)); }
    int rank0(int i) const { return i - rank1(i); }
};
template<class T> struct WaveletMatrix {
    int n, lg; vector<T> vals; vector<BitVec> lv; vector<int> Z;
    WaveletMatrix(const vector<T>& a) : n(sz(a)), vals(a) {
        sort(all(vals)); vals.erase(unique(all(vals)), vals.end());
        lg = __lg(max(sz(vals), 1)) + 1;
        vector<int> cur(n), nxt(n);
        for (int i = 0; i < n; i++) cur[i] = lower_bound(all(vals), a[i]) - vals.begin();
        lv.assign(lg, BitVec(n)); Z.assign(lg, 0);
        for (int h = lg - 1; h >= 0; h--) {
            int z = 0;
            for (int i = 0; i < n; i++) if (cur[i] >> h & 1) lv[h].set(i); else nxt[z++] = cur[i];
            lv[h].build(); Z[h] = z;
            for (int i = 0; i < n; i++) if (cur[i] >> h & 1) nxt[z++] = cur[i];
            swap(cur, nxt);
        }
    }
    T kth(int l, int r, int k) const {                 // 0-indexed k-th smallest of [l, r)
        int res = 0;
        for (int h = lg - 1; h >= 0; h--) {
            int l0 = lv[h].rank0(l), r0 = lv[h].rank0(r);
            if (k < r0 - l0) l = l0, r = r0;
            else k -= r0 - l0, res |= 1 << h, l += Z[h] - l0, r += Z[h] - r0;
        }
        return vals[res];
    }
    int countLess(int l, int r, T x) const {           // # i in [l, r) with a[i] < x
        int u = lower_bound(all(vals), x) - vals.begin(), res = 0;
        for (int h = lg - 1; h >= 0; h--) {
            int l0 = lv[h].rank0(l), r0 = lv[h].rank0(r);
            if (u >> h & 1) res += r0 - l0, l += Z[h] - l0, r += Z[h] - r0;
            else l = l0, r = r0;
        }
        return res;
    }
};

Predecessor of $x$ in the range: c = countLess(l, r, x), then kth(l, r, c - 1) if $c>0$. Successor: kth(l, r, c) if $c<r-l$. Storing prefix sums of original values per level also gives “sum of the $k$ smallest”.

No update support: the bit vectors would need re-partitioning on every change, which defeats the structure’s whole point. Use an offline sweep or the Fenwick-of-value-trees variant below instead when updates are needed. Beyond contest range queries, this is the structure behind succinct text indexing: the FM-index (a compressed full-text index built on the Burrows-Wheeler transform, used by DNA aligners like BWA and Bowtie), compressed suffix arrays, and document retrieval (“which documents contain this pattern”) all use a wavelet tree/matrix for their rank queries, precisely because the bits-not-integers memory saving matters at that scale.

Proof: the range mapping is correct

By induction on levels. Suppose that before level $h$ the current order lists the elements sorted stably by their bits above $h$, and the original range $[l,r)$ corresponds to a contiguous block $[l,r)$ in this order consisting of exactly the elements whose higher bits match the path chosen so far. The stable partition by bit $h$ places all zeros (in order) before all ones. The elements of the block with bit $h=0$ are exactly those among the first $r$ positions minus those among the first $l$ that have bit 0, and they appear contiguously at positions $[{\mathrm{rank}}_0(l),{\mathrm{rank}}_0(r))$ of the zeros part; the ones occupy $[Z_h+{\mathrm{rank}}_1(l),Z_h+{\mathrm{rank}}_1(r))$. For $k$-th smallest, choosing the zero side when $k<$ (number of zeros in the block) keeps the invariant “the answer is the $k$-th smallest of the current block”, since every value with bit $h=0$ is smaller than every value with bit $h=1$ given equal higher bits.

Merge sort tree

A segment tree whose node stores the sorted values of its range (built by merging children). “Count $\le x$ in $[l,r)$” binary searches in $O(\log n)$ canonical nodes: $O({\log }^{2}n)$. Fractional cascading (store, for each element of a node, its position in each child’s list) removes the inner binary searches for $O(\log n)$ queries, at the cost of significant bookkeeping for a modest constant — rarely worth writing in a contest, but the standard example of the technique. It’s simple and hard to get wrong, but the wavelet matrix dominates it on time and memory.

Updates: replacing an element needs re-sorting every node on its root path, $O(n)$ worst case — the merge sort tree is a static structure by design. If updates are genuinely needed, use the offline BIT sweep above, a sqrt decomposition with sorted blocks (rebuild one block in $O(B)$), or replace each node’s sorted vector with a multiset/BIT ($O({\log }^{2}n)$ updates, worse constant, $O(n\log n)$ node overhead).

With updates: Fenwick of value trees

A Fenwick tree over positions whose nodes are dynamic segment trees over values. A point update touches $O(\log n)$ Fenwick nodes, each $O(\log V)$. A query collects the $O(\log n)$ roots for prefix $r$ (sign $+$) and prefix $l$ (sign $-$) and descends all of them at once.

int kth(int l, int r, int k) {                         // positions [l, r), 0-indexed k
    vector<int> P, M;
    for (int i = r; i > 0; i &= i - 1) P.push_back(root[i - 1]);
    for (int i = l; i > 0; i &= i - 1) M.push_back(root[i - 1]);
    int lo = 0, hi = V;
    while (hi - lo > 1) {
        int mid = (lo + hi) / 2, left = 0;
        for (int p : P) left += cnt[lc[p]];
        for (int p : M) left -= cnt[lc[p]];
        if (k < left) { for (int& p : P) p = lc[p]; for (int& p : M) p = lc[p]; hi = mid; }
        else { k -= left; for (int& p : P) p = rc[p]; for (int& p : M) p = rc[p]; lo = mid; }
    }
    return lo;
}

$O(\log n\log V)$ per operation and $O((n+q)\log n\log V)$ memory. The offline alternative (overall binary search with a single Fenwick) uses $O(n)$ memory and is usually preferable.

On trees

$k$-th smallest on the path $u$–$v$: persistent value trees along root paths, then descend the combination $T_u+T_v-T_{\text{lca}}-T_{\text{parent(lca)}}$. In subtrees: Euler tour order turns it into the array problem.

Related queries

Query	How
median of $a[l..r)$	$k$-th with $k=\lfloor (r-l-1)/2\rfloor$
number of distinct values in a range	not an order statistic: offline sweep with last occurrences, or persistent tree on positions (below)
$k$-th distinct value in a range	the distinct-count structure below, plus a descent
sum of $|a_i-x|$ over a range	count and sum of values below $x$ (wavelet with sums, or PST with sums)
range mode	sqrt decomposition (Sqrt Decomposition)
$k$-th smallest of the whole array	nth_element, or a Fenwick descent if dynamic

Distinct values in a range

CSES “Distinct Values Queries”: count distinct values in $a[l,r)$. Not an order statistic (no comparison structure needed), but the same toolbox solves it.

Key idea: a value is “new” at position $i$ iff $i$ is its last occurrence so far looking backward, i.e. prev[i] (the previous equal element, or $-1$) is $<l$. So the count over $[l,r)$ is the number of $i\in [l,r)$ with prev[i] < l.

vector<int> lastOccurrence(const vector<int>& a) {     // prev[i] = largest j < i with a[j] == a[i], or -1
    int n = sz(a); vector<int> order(n), prev(n, -1);
    iota(all(order), 0);
    stable_sort(all(order), [&](int i, int j) { return a[i] < a[j]; });
    for (int k = 1; k < n; k++) if (a[order[k - 1]] == a[order[k]]) prev[order[k]] = order[k - 1];
    return prev;
}
// Offline: sweep r, keep a 1 at i iff i is currently the last occurrence of its value among [0, r).
vector<int> distinctInRanges(const vector<int>& a, const vector<pair<int, int>>& qs) {   // [l, r)
    int n = sz(a), q = sz(qs);
    vector<int> prev = lastOccurrence(a), res(q);
    vector<vector<int>> byRight(n + 1);
    for (int i = 0; i < q; i++) byRight[qs[i].second].push_back(i);
    Fenwick<int> fw(n);
    for (int r = 0; r <= n; r++) {
        for (int i : byRight[r]) res[i] = fw.sum(qs[i].first, r);
        if (r == n) break;
        if (prev[r] >= 0) fw.add(prev[r], -1);
        fw.add(r, 1);
    }
    return res;
}

$O((n+q)\log n)$ offline. For online queries, build a persistent segment tree over positions: version $r$ is the indicator array “is $i$ the last occurrence of its value among $[0,r)$”, and query(l, r) descends version $r$ summing over $[l,n)$… equivalently sums the indicator over $[l,r)$ directly.

struct DistinctOnline {                                 // query(l, r): # distinct values in [l, r)
    int n; vector<int> lc{0}, rc{0}, cnt{0}, roots{0};
    DistinctOnline(const vector<int>& a) : n(sz(a)) {
        vector<int> prev = lastOccurrence(a);
        for (int r = 0; r < n; r++) {
            int root = roots.back();
            if (prev[r] >= 0) root = add(root, 0, n, prev[r], -1);  // the old last-occurrence stops counting
            roots.push_back(add(root, 0, n, r, 1));                 // r becomes the new last occurrence
        }
    }
    int add(int p, int l, int r, int pos, int v) {
        int q = sz(lc); lc.push_back(lc[p]); rc.push_back(rc[p]); cnt.push_back(cnt[p] + v);
        if (r - l > 1) {
            int m = (l + r) / 2;
            if (pos < m) lc[q] = add(lc[p], l, m, pos, v); else rc[q] = add(rc[p], m, r, pos, v);
        }
        return q;
    }
    int sumFrom(int p, int l, int r, int x) const {      // sum of the indicator over [x, n), in version p
        if (!p || r <= x) return 0;
        if (x <= l) return cnt[p];
        int m = (l + r) / 2;
        return sumFrom(lc[p], l, m, x) + sumFrom(rc[p], m, r, x);
    }
    int query(int l, int r) const { return l < r ? sumFrom(roots[r], 0, n, l) : 0; }  // version r, positions >= l
};

$O(\log n)$ per query, $O(n\log n)$ build and memory.

Proof: the last-occurrence indicator counts distinct values correctly

For a fixed right endpoint $r$, exactly one position $i<r$ carries the indicator $1$ per distinct value present in $a[0,r)$: among all occurrences of that value before $r$, the indicator is set only at the last one (prev[i]<l is not the condition here — the sweep version marks $i$ with $+1$ and demotes the previous occurrence of the same value to $-1$ the moment a later one appears, so at any snapshot exactly the most recent occurrence of each value carries $+1$). So ${\sum }_{i<r}\text{indicator}(i)=$ number of distinct values in $a[0,r)$, and restricting the sum to $i\ge l$ counts a value iff its last occurrence before $r$ is $\ge l$ — which is exactly “this value appears somewhere in $[l,r)$” (if its last occurrence before $r$ is $\ge l$, that occurrence itself witnesses membership in $[l,r)$; if its last occurrence before $r$ is $<l$, every occurrence of it is $<l$, so it doesn’t appear in $[l,r)$ at all). The offline version maintains this indicator incrementally in a Fenwick tree as $r$ sweeps forward (removing the previous occurrence’s $+1$, adding one at $r$) and answers each query the moment its right endpoint is reached. The online version stores one persistent-tree version per value of $r$, so query(l, r) reads off exactly the same sum from the frozen snapshot at that $r$. $■$

Tested against a brute-force set-based distinct count on random arrays (values from a small alphabet to force repeats) and random ranges, for both the offline and online versions — all match.