Order Statistics & Quantiles

Figure 1 plots realizations of $F_n$ for iid standard-Normal samples at $n = 20$, $200$, and $1000$. At $n = 20$, $F_n$ is a jagged step-function whose staircase shape is very different from the smooth $\Phi$. At $n = 200$, the steps compress into what visually reads as $\Phi$ with small wobbles. At $n = 1000$, the wobbles are invisible at this resolution — the ECDF is a plausible stand-in for $\Phi$. The annotated supremum distance $D_n = \sup_x |F_n(x) - \Phi(x)|$ shrinks from 0.22 at $n = 20$ to 0.06 at $n = 200$ to 0.03 at $n = 1000$, following the $1/\sqrt n$ rate that §29.5’s DKW inequality proves non-asymptotically.

Nonparametric statistics has two distinct meanings in the ML literature. The one Topic 29 cares about is distribution-free inference: the validity of a procedure does not depend on the form of $F$. A bootstrap CI, a KS test, a Wilcoxon test, and a quantile CI built from order-statistic pairs are all distribution-free. The other meaning — “flexible modeling with infinite-dimensional parameters” — covers Bayesian nonparametrics (Dirichlet processes, Gaussian processes) and kernel methods; those are in formalml.

Distribution-freeness is what makes KS useful for real ML workflows: detecting distribution shift between training and deployment samples, validating that a model’s calibration stays the same across data slices, or checking that a simulator reproduces the empirical distribution of a sensor. None of those applications specify $F$ up front, so a t-test is the wrong tool.

Classical statistics assumes a parametric model; nonparametric statistics lets the data specify the model via the empirical distribution, and the order statistics are the lens through which we study that empirical distribution. §29.10’s forward-map spells out how Topics 30 (KDE), 31 (Bootstrap), and 32 (Empirical Processes) each sharpen one aspect of this sentence.

For iid $X_1, \dots, X_n \sim F$ with $F$ continuous, the order statistics are the sorted sample, written with parenthesized indices:

The parentheses distinguish $X_{(i)}$ (the $i$-th smallest value in the sample) from $X_i$ (the $i$-th draw, in the order collected). When $F$ is continuous, ties occur with probability zero, so the inequalities are strict almost surely.

For iid $X_1, \dots, X_n$ with continuous CDF $F$ and density $f$, the joint density of $(X_{(1)}, \dots, X_{(n)})$ on the ordered simplex $\{x_1 \le x_2 \le \cdots \le x_n\}$ is

and zero outside that simplex. The marginal density of the $i$-th order statistic is

(DAV2003 §2.1–§2.2.)

For $i = 1$, the marginal density reduces to $n\,(1 - F(x))^{n - 1}\,f(x)$ — the minimum’s density is the density at $x$ times the probability that the other $n - 1$ draws are all above. For $i = n$, the marginal is $n\,F(x)^{n - 1}\,f(x)$ — the maximum’s density is $f(x)$ times the probability the rest are all below. In both cases, as $n$ grows, the density concentrates near the lower or upper support endpoint respectively. The rate is governed by $F(x)^{n-1}$ (upper tail) or $(1 - F(x))^{n-1}$ (lower tail), which decays exponentially in $n$ away from the boundary — a topic fully developed at formalml extreme-value-theory.

The joint density of $(X_{(i)}, X_{(j)})$ for $i < j$ comes from marginalizing Thm 1:

valid for $x_i < x_j$. The combinatorial factor counts placements: $i - 1$ below $x_i$, one at $x_i$, $j - i - 1$ between, one at $x_j$, $n - j$ above. §29.7 uses the special case $(i, j)$ to build distribution-free quantile CIs — pick $r, s$ so that $\Pr\bigl(X_{(r)} \le \xi_p \le X_{(s)}\bigr) \ge 1 - \alpha$.

For non-iid data (record values, concomitants, ranked-set sampling, progressive censoring) the joint density takes a different combinatorial form — the $n!$ symmetry is broken. DAV2003 Chapters 5–6 develop the non-iid theory in full; Topic 29 restricts to iid throughout. The one use case outside iid order statistics that recurs in ML is censoring — survival analysis estimates $F$ from partially observed data — which we flag as an extension but do not develop.

Topic 16 §16 Thm 1 established that for any iid model, the order-statistic tuple $T(X) = (X_{(1)}, \dots, X_{(n)})$ is sufficient for $\theta$: the conditional distribution of $X$ given $T$ does not depend on $\theta$. The same result holds nonparametrically — for any $F$, the order statistic captures every feature of the sample that $F$ influences. This is the structural justification for Track 8: the objects §29 develops are the raw features of every nonparametric procedure, because no statistic that respects the iid assumption can ever discard information in the order statistic. Independence notation $\perp\!\!\!\perp$ (Topic 16 §16.9) carries forward: $X_{(1)}, \dots, X_{(n)}$ are not independent (ordering is the whole point), but for iid $X_i$, the vector of ranks $\{R_i\}$ satisfies $(R_1, \dots, R_n) \perp\!\!\!\perp (X_{(1)}, \dots, X_{(n)})$ — a fact exploited by rank tests (forthcoming on formalml).

For iid $U_1, \dots, U_n \sim \mathrm{Uniform}(0, 1)$, the $i$-th order statistic $U_{(i)}$ has density

i.e. $U_{(i)} \sim \mathrm{Beta}(i,\, n - i + 1)$. This is Thm 1 specialized to $F(u) = u,\, f(u) = 1$ on $(0,1)$; the Beta-density normalization follows from $B(i, n-i+1) = (i-1)!(n-i)!/n!$.

Corollary 1 (general $F$). If $X_{(i)}$ is the $i$-th order statistic of an iid $F$-sample and $F$ is continuous, then $F(X_{(i)}) \sim \mathrm{Beta}(i, n - i + 1)$, since $F(X_i) \sim \mathrm{Uniform}(0, 1)$ under the probability-integral transform, so sorting the $F$-values gives the Uniform order statistics.

(DAV2003 §2.2.)

Let $Y_1, \dots, Y_n$ be iid $\mathrm{Exp}(1)$. Define the normalized spacings $W_k = Y_{(k)} - Y_{(k - 1)}$ for $k = 1, \dots, n$ (with $Y_{(0)} = 0$). Then the rescaled spacings $E_k := (n - k + 1)\, W_k$ are iid $\mathrm{Exp}(1)$, and

where $E_1, \dots, E_n$ on the right-hand side are iid $\mathrm{Exp}(1)$.

From Theorem 3,

where $H_m = 1 + \tfrac{1}{2} + \cdots + \tfrac{1}{m}$ is the harmonic number. The minimum’s expectation is $H_n - H_{n - 1} = 1/n$ — the smallest of $n$ iid Exp(1) variables is itself Exp($n$), mean $1/n$. The maximum’s expectation is $H_n - H_0 = H_n \approx \log n + \gamma$ — the largest of $n$ iid Exp(1) grows logarithmically. Both are lessons that recur in EVT (forthcoming on formalml): the maximum of iid Exp grows slowly, which is why the Gumbel distribution is the right asymptotic limit for Gumbel-domain tails.

The factor $n - k + 1$ in Theorem 3 counts the “remaining” iid samples at step $k$: among $n$ iid Exp(1) draws, the minimum is Exp($n$); conditional on the minimum, the second smallest is the minimum of the remaining $n - 1$ — which by memorylessness is Exp($n - 1$); and so on. The spacing structure is memorylessness in disguise. This is the exponential-family property that makes §29.6’s Bahadur residual have a clean form — the tail of the ECDF at quantile $p$ decomposes into iid contributions via the probability-integral transform onto Uniform, and then onto Exp via $U \mapsto -\log U$.

Rényi’s representation is the order-statistic reading of a Poisson-process fact: if we sprinkle $n$ iid $\mathrm{Exp}(\lambda)$-distributed inter-arrival times on $[0, \infty)$, we get a Poisson process of rate $\lambda$. The sorted arrival times are $Y_{(1)} < Y_{(2)} < \cdots$, and the conditional distribution given $Y_{(n+1)} = T$ is uniform order statistics on $[0, T]$. Theorem 3 is a direct restatement of the inter-arrival independence in the forward direction. Point-process machinery (formalml point-processes, forthcoming) makes this connection structural.

Topic 16 §16.14 Ex 13 proved an unexpected consequence of Basu’s theorem: for iid $\mathrm{Uniform}(0, \theta)$, the ratio $X_{(n - 1)} / X_{(n)}$ is independent of the maximum $X_{(n)}$. The uniform-order-statistics structure Theorem 2 and Rényi’s representation exploit is the reason: once you change variables to the Uniform-order-statistics simplex, the joint density factorizes in ways that let Basu’s theorem detect independence cleanly. Basu in Topic 16 was an early hint of the uniform-order-statistics structure that §29.3 develops in full.

The population quantile at level $p \in (0, 1)$ is

For continuous strictly-increasing $F$, $\xi_p = F^{-1}(p)$ in the usual inverse-function sense. The left-continuous-inverse definition handles atoms and flat regions of $F$ cleanly.

Let $X_{(1)} \le \cdots \le X_{(n)}$ be the order statistics of an iid sample. The Type-7 sample quantile at level $p \in [0, 1]$ is

Equivalently: $\hat\xi_p$ linearly interpolates between the two order statistics bracketing the position $(n - 1)\, p + 1$ in the sorted sample. For $p = 0$, $\hat\xi_p = X_{(1)}$; for $p = 1$, $\hat\xi_p = X_{(n)}$; for $p = 0.5$ with odd $n$, $\hat\xi_p = X_{((n + 1)/2)}$, the sample median.

(HYN1996, default convention in R quantile(), NumPy quantile, SciPy scipy.stats.mstats.mquantiles.)

For the sample $\{1, 2, 3, 4, 5\}$ with $n = 5$: at $p = 0.5$, $h = (5 - 1)(0.5) = 2$, so $\hat\xi_{0.5} = X_{(3)} = 3$ exactly. At $p = 0.25$, $h = 1$, so $\hat\xi_{0.25} = X_{(2)} = 2$. At $p = 0.9$, $h = 3.6$, so $\hat\xi_{0.9} = X_{(4)} + 0.6 (X_{(5)} - X_{(4)}) = 4 + 0.6 = 4.6$. These match the numpy.quantile([1,2,3,4,5], q) defaults exactly — a sanity check that the Type-7 convention is what the reader’s tooling uses.

The nine Hyndman–Fan conventions differ in how they interpolate when $p$ is not at an integer order-statistic position. Type 1 (inverse-CDF) picks $X_{(\lceil n p \rceil)}$ with no interpolation — it’s discontinuous in $p$. Type 5 and Type 6 are older conventions still used in some engineering packages. Type 7 is the unique member of the family that (i) is a continuous piecewise-linear function of $p$ and (ii) hits the boundary exactly ($\hat\xi_0 = X_{(1)}$ and $\hat\xi_1 = X_{(n)}$). Because R, NumPy, and SciPy all adopt Type 7 by default, a reader writing Python code to reproduce §29’s figures gets Type-7 numbers without setting a flag. Any other convention would introduce a silent discrepancy.

Exact sample quantiles require sorting the sample — $O(n \log n)$ and $O(n)$ memory. In streaming or high-cardinality settings, approximate quantile sketches (Jain–Chlamtac 1985’s P² algorithm; Greenwald–Khanna 2001; Dunning’s t-digest 2021) trade small guaranteed error for bounded memory. These are widely used in ML production systems for real-time latency monitoring. Topic 29 does not develop them — streaming quantiles are an engineering sub-topic of its own. The relevant entry point is Dunning 2021 (t-digest) for practical use.

For iid $X_1, \dots, X_n \sim F$ with $F$ any CDF,

The supremum distance $D_n := \sup_x |F_n(x) - F(x)|$ is the one-sample Kolmogorov–Smirnov statistic. Glivenko–Cantelli says $D_n \to 0$ almost surely. (Topic 10 §10.7 Thm 5; see there for the ε-net argument.)

For iid $X_1, \dots, X_n \sim F$ and every $\varepsilon > 0$,

The constant $2$ is optimal (Massart 1990) when $\varepsilon \ge \sqrt{\ln 2 / (2n)}$; DVO1956 originally proved the bound with a larger, unspecified constant. Inverting the bound at confidence level $1 - \alpha$ gives a distribution-free finite-sample confidence envelope: the band $[F_n(x) - \varepsilon_n, F_n(x) + \varepsilon_n]$ with

contains $F$ uniformly in $x$ with probability $\ge 1 - \alpha$.

Corollary 2 (the DKW envelope). With $\varepsilon_n$ as above, $\Pr\bigl(F_n(x) - \varepsilon_n \le F(x) \le F_n(x) + \varepsilon_n \text{ for all } x\bigr) \ge 1 - \alpha$ regardless of $F$. (DVO1956; MAS1990.)

A common ML workflow: you suspect the distribution of a production feature has shifted since training. Collect $n_1$ training samples and $n_2$ production samples. Compute the ECDFs $F_{n_1}^{\text{train}}$ and $F_{n_2}^{\text{prod}}$ and their sup-norm distance $D = \sup_x |F_{n_1}^{\text{train}}(x) - F_{n_2}^{\text{prod}}(x)|$. The one-sample DKW envelopes give simultaneous confidence bands for the true train and production CDFs; if those bands overlap at every $x$, no detectable shift. For the formal two-sample test — including the asymptotic null distribution — see §29.8 (with the two-sample variant a one-line extension to ksTwoSample in nonparametric.ts). Production systems run this test nightly to flag shift on each logged feature.

CLT says $\sqrt n\, (F_n(x) - F(x)) \Rightarrow \mathcal{N}(0, F(x)(1 - F(x)))$ for each fixed $x$ — a pointwise asymptotic statement. DKW says $\sup_x |F_n(x) - F(x)|$ is bounded non-asymptotically with probability $\ge 1 - \alpha$ — a uniform finite-sample statement. Both are useful: CLT gives the pointwise limiting distribution that §29.8 Thm 7 upgrades to uniform (via Donsker); DKW gives the uniform bound that does not require any asymptotic machinery. The Bahadur proof (§29.6) exploits both: CLT at the point $\xi_p$ plus a uniform-over-neighborhood fluctuation bound that DKW is the finite-sample analog of.

For $x$ at the median of $F$, CLT gives asymptotic variance $F(x)(1 - F(x)) = 1/4$. DKW’s exponential rate $2 e^{-2n\varepsilon^2}$ at $\varepsilon = 1/\sqrt{2n}$ gives probability $\ge 2/e \approx 0.74$ that the sup-norm gap is below $1/\sqrt{2n}$, which matches the CLT tail probability to within constants. For $x$ deep in the tails where $F(x) \approx 0$ or $1$, CLT gives much smaller pointwise variance than $1/4$, so the DKW uniform bound is loose at tails — a property exploited by the tail-aware inequalities of Talagrand (Talagrand 1994) developed at formalml concentration-inequalities (forthcoming).

Let $X_1, X_2, \dots, X_n$ be iid with CDF $F$. Fix $p \in (0, 1)$ and assume $F$ is differentiable at $\xi_p = F^{-1}(p)$ with $F'(\xi_p) = f(\xi_p) > 0$ and that $F$ is twice-differentiable in a neighborhood of $\xi_p$ with bounded second derivative. Let $\hat\xi_p$ be the Type-7 sample quantile (Def 3). Then

where the remainder satisfies $R_n = O_p\!\bigl(n^{-3/4}\,\sqrt{\log n}\bigr)$ (Kiefer 1967 rate).

Corollary 1 (asymptotic normality). Under the same hypotheses,

(BAH1966; KIE1967; SER1980 §2.5; vdV2000 §21.)

For $X_i \sim \mathcal{N}(0, 1)$ at $p = 0.5$: $\xi_{0.5} = 0$, $f(0) = 1/\sqrt{2\pi}$. Bahadur Cor 1 gives

The limiting standard error of the sample median is $\sqrt{\pi/2} \approx 1.253$, compared to $1.0$ for the sample mean. The median is less efficient than the mean under exact Normal data by the factor $2/\pi \approx 0.637$ — the asymptotic relative efficiency ARE(median, mean) = $2/\pi$. Under heavier tails this ratio flips; see Ex 8.

For $X_i \sim \mathrm{Cauchy}(0, 1)$ at $p = 0.5$: $\xi_{0.5} = 0$, $f(0) = 1/\pi$. Bahadur Cor 1 gives

The limiting standard error is $\pi/2 \approx 1.571$ — slightly larger than the Normal case, but finite. The sample mean of Cauchy data is undefined asymptotically (the variance of $\bar X$ is infinite; CLT does not apply), so the sample median is infinitely more efficient than the sample mean for Cauchy — ARE is $+\infty$. This is the money-shot example that motivates Bahadur: regularity is about the density at $\xi_p$, not about moments.

The limit variance $p(1-p)/f(\xi_p)^2$ depends on $f(\xi_p)$ — a nuisance parameter the reader typically does not know. Plug-in estimation requires estimating $f$ at $\xi_p$. Topic 30 (Kernel Density Estimation) develops the general KDE machinery that provides consistent estimators $\hat f_n(\xi_p)$ with controlled bandwidth. The density-at-quantile problem is discharged at Topic 30 §30.7 Ex 9, which computes a Silverman-bandwidth KDE at the empirical quantile and plugs the result into Bahadur’s variance formula. Alternatively, the bootstrap (Topic 31: The Bootstrap) sidesteps this entirely by resampling the full empirical process — the plug-in variance is never computed directly.

The Bahadur representation is structurally identical to the Taylor-expansion-of-the-score-function argument that underpins MLE asymptotic normality (Topic 14 §14.6):

where the leading term is a mean-zero empirical average divided by a fixed constant, and the remainder is asymptotically negligible. Bahadur’s remainder rate $n^{-3/4}\sqrt{\log n}$ is worse than MLE’s $o_p(n^{-1/2})$, but the functional form — “linearize into a known empirical average, bound the remainder, apply CLT” — is identical. Every Track 8 result reduces to a Bahadur-style linearization: KDE bandwidth asymptotics (Topic 30), bootstrap-CI validity (Topic 31: The Bootstrap), empirical-process functionals (Topic 32), and quantile-regression (forthcoming on formalml) all use this template.

Bahadur requires $f(\xi_p) > 0$. At extreme $p$ — near 0 or 1 — the density may vanish ($f(\xi_0)$ or $f(\xi_1)$ may be zero or undefined if the support is unbounded). In those regimes the sample quantile’s asymptotic distribution is not Normal at rate $\sqrt n$; it is governed by the Fisher–Tippett–Gnedenko trichotomy (Gumbel, Fréchet, Weibull) at rate determined by the tail behavior of $F$. Extreme-value theory is the subject of formalml extreme-value-theory (forthcoming); Topic 29 restricts to interior $p$ where Bahadur applies.

For iid $X_1, \dots, X_n$ with continuous CDF $F$, and $p \in (0, 1)$, let $r < s$ be integers with $1 \le r \le s \le n$. Then

regardless of $F$. Choosing $r, s$ so that this binomial probability is $\ge 1 - \alpha$ yields a distribution-free exact $(1 - \alpha)$-CI for $\xi_p$. The canonical choice picks the largest $r$ with $\Pr(B < r) \le \alpha / 2$ and the smallest $s$ with $\Pr(B \ge s) \le \alpha / 2$; this gives a two-sided CI with exact coverage $1 - \Pr(B < r) - \Pr(B \ge s) \ge 1 - \alpha$. (SER1980 §2.6.1.)