Studies

This section analyzes the estimators properties using mathematical proofs. Most proofs are adapted from various textbooks and research papers, but only essential references are provided.

Unlike the main part of the manual, the studies require knowledge of classic statistical methods. Well-known facts and commonly accepted notation are used without special introduction. The studies provide detailed analyses of estimator properties for practitioners interested in rigorous proofs and numerical simulation results.

Summary Estimator Properties

This section compares the toolkits robust estimators against traditional statistical methods to demonstrate their advantages across diverse conditions. While traditional estimators often work well under ideal conditions, the toolkits estimators maintain reliable performance across diverse real-world scenarios.

Average Estimators:

Mean (arithmetic average):

Median:

Center (Hodges-Lehmann estimator):

Dispersion Estimators:

Standard Deviation:

Median Absolute Deviation (around the median):

Spread (Shamos scale estimator):

Breakdown

Heavy-tailed distributions naturally produce extreme outliers that completely distort traditional estimators. A single extreme measurement from the $\underline{\operatorname{Power}}$ distribution can make the sample mean arbitrarily large. Real-world data can also contain corrupted measurements from instrument failures, recording errors, or transmission problems. Both natural extremes and data corruption create the same challenge: extracting reliable information when some measurements are too influential.

The breakdown point (Huber 2009) is the fraction of a sample that can be replaced by arbitrarily large values without making an estimator arbitrarily large. The theoretical maximum is $50\%$; no estimator can guarantee reliable results when more than half the measurements are extreme or corrupted. In such cases, summary estimators are not applicable, and a more sophisticated approach is needed.

A $50\%$ breakdown point is rarely needed in practice, as more conservative values also cover practical needs. Additionally, a high breakdown point corresponds to low precision (information is lost by neglecting part of the data). The optimal practical breakdown point should be between $0\%$ (no robustness) and $50\%$ (low precision).

The $\operatorname{Center}$ and $\operatorname{Spread}$ estimators achieve $29\%$ breakdown points, providing substantial protection against realistic contamination levels while maintaining good precision. Below is a comparison with traditional estimators.

Asymptotic breakdown points for average estimators:

Asymptotic breakdown points for dispersion estimators:

Drift

Drift measures estimator precision by quantifying how much estimates scatter across repeated samples. It is based on the $\operatorname{Spread}$ of estimates and therefore has a breakdown point of approximately $29\%$.

Drift is useful for comparing the precision of several estimators. To simplify the comparison, one of the estimators can be chosen as a baseline. A table of squared drift values, normalized by the baseline, shows the required sample size adjustment factor for switching from the baseline to another estimator. For example, if $\operatorname{Center}$ is the baseline and the rescaled drift square of $\operatorname{Median}$ is $1.5$, this means that $\operatorname{Median}$ requires $1.5$ times more data than $\operatorname{Center}$ to achieve the same precision. See the From Statistical Efficiency to Drift section for details.

Squared Asymptotic Drift of Average Estimators (values are approximated):

Rescaled to $\operatorname{Center}$ (sample size adjustment factors):

Squared Asymptotic Drift of Dispersion Estimators (values are approximated):

Rescaled to $\operatorname{Spread}$ (sample size adjustment factors):

Invariance

Invariance properties determine how estimators respond to data transformations. These properties are crucial for analysis design and interpretation:

• Location-invariant estimators are invariant to additive shifts: $T(\mathbf{x}+k)=T(\mathbf{x})$
• Scale-invariant estimators are invariant to positive rescaling: $T(k \cdot \mathbf{x})=T(\mathbf{x})$ for $k>0$
• Equivariant estimators change predictably with transformations, maintaining relative relationships

Choosing estimators with appropriate invariance properties ensures that results remain meaningful across different measurement scales, units, and data transformations. For example, when comparing datasets collected with different instruments or protocols, location-invariant estimators eliminate the need for data centering, while scale-invariant estimators eliminate the need for normalization.

Location-invariance: An estimator $T$ is location-invariant if adding a constant to the measurements leaves the result unchanged:

Location-equivariance: An estimator $T$ is location-equivariant if it shifts with the data:

Scale-invariance: An estimator $T$ is scale-invariant if multiplying by a positive constant leaves the result unchanged:

Scale-equivariance: An estimator $T$ is scale-equivariant if it scales proportionally with the data:

Reframings

From Statistical Efficiency to Drift

Statistical efficiency measures estimator precision (Serfling 2009). When multiple estimators target the same quantity, efficiency determines which provides more reliable results.

Efficiency measures how tightly estimates cluster around the true value across repeated samples. For an estimator $T$ applied to samples from distribution $X$, absolute efficiency is defined relative to the optimal estimator $T^*$:

Relative efficiency compares two estimators by taking the ratio of their variances:

Under $\underline{\operatorname{Additive}}$ (Normal) distributions, this approach works well. The sample mean achieves optimal efficiency, while the median operates at roughly 64% efficiency.

However, this variance-based definition creates four critical limitations:

• Absolute efficiency requires knowing the optimal estimator, which is difficult to determine. For many distributions, deriving the minimum-variance unbiased estimator requires complex mathematical analysis. Without this reference point, absolute efficiency cannot be computed.
• Relative efficiency only compares estimator pairs, preventing systematic evaluation. This limits understanding of how multiple estimators perform relative to each other. Practitioners cannot rank estimators comprehensively or evaluate individual performance in isolation.
• The approach depends on variance calculations that break down when variance becomes infinite or when distributions have heavy tails. Many real-world distributions, such as those with power-law tails, exhibit infinite variance. When the variance is undefined, efficiency comparisons become impossible.
• Variance is not robust to outliers, which can corrupt efficiency calculations. A single extreme observation can greatly inflate variance estimates. This sensitivity can make efficient estimators look inefficient and vice versa.

The $\operatorname{Drift}$ concept provides a robust alternative. Drift measures estimator precision using $\operatorname{Spread}$ instead of variance, providing reliable comparisons across a wide range of distributions.

For an average estimator $T$, random variable $X$, and sample size $n$:

This formula measures estimator variability compared to data variability. $\operatorname{Spread}[T(X_1, \ldots, X_n)]$ captures the median absolute difference between estimates across repeated samples. Multiplying by $\sqrt{n}$ removes sample size dependency, making drift values comparable across different sample sizes. Dividing by $\operatorname{Spread}[X]$ creates a scale-free measure that provides consistent drift values across different distribution parameters and measurement units.

Dispersion estimators use a parallel formulation:

Here $\operatorname{RelSpread}$ normalizes by the estimators typical value for fair comparison.

Drift offers four key advantages:

• For estimators with $\sqrt{n}$ convergence rates, drift remains finite and comparable across distributions; for heavier tails, drift may diverge, flagging estimator instability.
• It provides absolute precision measures rather than only pairwise comparisons.
• The robust $\operatorname{Spread}$ foundation resists outlier distortion that corrupts variance-based calculations.
• The $\sqrt{n}$ normalization makes drift values comparable across different sample sizes, enabling direct comparison of estimator performance regardless of sample size.

Under $\underline{\operatorname{Additive}}$ (Normal) conditions, drift matches traditional efficiency. The sample mean achieves drift near 1.0; the median achieves drift around 1.25. This consistency validates drift as a proper generalization of efficiency that extends to realistic data conditions where traditional efficiency fails.

When switching from one estimator to another while maintaining the same precision, the required sample size adjustment follows:

This applies when estimator $T_1$ has lower drift than $T_2$.

The ratio of squared drifts determines the data requirement change. If $T_2$ has drift 1.5 times higher than $T_1$, then $T_2$ requires $(1.5)^2 = 2.25$ times more data to match $T_1$s precision. Conversely, switching to a more precise estimator allows smaller sample sizes.

For asymptotic analysis, $\operatorname{Drift}(T, X)$ denotes the limiting value as $n \to \infty$. With a baseline estimator, rescaled drift values enable direct comparisons:

The standard drift definition assumes $\sqrt{n}$ convergence rates typical under $\underline{\operatorname{Additive}}$ (Normal) conditions. For broader applicability, drift generalizes to:

The instability parameter adapts to estimator convergence rates. The toolkit uses $\text{instability} = \frac{1}{2}$ throughout because this choice provides natural intuition and mental representation for the $\underline{\operatorname{Additive}}$ (Normal) distribution. Rather than introduce additional complexity through variable instability parameters, the fixed $\sqrt{n}$ scaling offers practical convenience while maintaining theoretical rigor for the distribution classes most common in applications.

From Confidence Level to Misrate

Traditional statistics expresses uncertainty through confidence levels: 95% confidence interval, 99% confidence, 99.9% confidence. This convention emerged from early statistical practice when tables printed confidence intervals for common levels like 90%, 95%, and 99%.

The confidence level approach creates practical problems:

• Cognitive difficulty with high confidence. Distinguishing between 99.999% and 99.9999% confidence requires mental effort. The difference matters — one represents a 1-in-100,000 error rate, the other 1-in-1,000,000 — but the representation obscures this distinction.
• Asymmetric scale. The confidence level scale compresses near 100%, where most practical values cluster. Moving from 90% to 95% represents a 2× change in error rate, while moving from 99% to 99.9% represents a 10× change, despite similar visual spacing.
• Indirect interpretation. Practitioners care about error rates, not success rates. Whats the chance Im wrong? matters more than Whats the chance Im right? Confidence level forces mental subtraction to answer the natural question.
• Unclear defaults. Traditional practice offers no clear default confidence level. Different fields use different conventions (95%, 99%, 99.9%), creating inconsistency and requiring arbitrary choices.

The $\mathrm{misrate}$ parameter provides a more natural representation. Misrate expresses the probability that computed bounds fail to contain the true value:

This simple inversion provides several advantages:

• Direct interpretation. $\mathrm{misrate} = 0.01$ means 1% chance of error or wrong 1 time in 100. $\mathrm{misrate} = 10^{-6}$ means wrong 1 time in a million. No mental arithmetic required.
• Linear scale for practical values. $\mathrm{misrate} = 0.1$ (10%), $\mathrm{misrate} = 0.01$ (1%), $\mathrm{misrate} = 0.001$ (0.1%) form a natural sequence. Scientific notation handles extreme values cleanly: $10^{-3}$, $10^{-6}$, $10^{-9}$.
• Clear comparisons. $10^{-5}$ versus $10^{-6}$ immediately shows a 10× difference in error tolerance. 99.999% versus 99.9999% confidence obscures this same relationship.
• Pragmatic default. The toolkit recommends $\mathrm{misrate} = 10^{-3}$ (one-in-a-thousand error rate) as a reasonable default for everyday analysis. For critical decisions where errors are costly, use $\mathrm{misrate} = 10^{-6}$ (one-in-a-million).

The terminology shift from confidence level to misrate parallels other clarifying renames in this toolkit. Just as $\underline{\operatorname{Additive}}$ better describes the distributions formation than Normal, and $\operatorname{Center}$ better describes the estimators purpose than Hodges-Lehmann, $\mathrm{misrate}$ better describes the quantity practitioners actually reason about: the probability of error.

Traditional confidence intervals become bounds in this framework, eliminating statistical jargon in favor of descriptive terminology. $\operatorname{ShiftBounds}(\mathbf{x}, \mathbf{y}, \mathrm{misrate})$ clearly indicates: it provides bounds on the shift, with a specified error rate. No background in classical statistics required to understand the concept.

From Mann-Whitney U-test to Pairwise Margin

The Mann-Whitney $U$ test (also known as the Wilcoxon rank-sum test) ranks among the most widely used non-parametric statistical tests, testing whether two independent samples come from the same distribution. Under $\underline{\operatorname{Additive}}$ (Normal) conditions, it achieves nearly the same precision as the Students $t$-test, while maintaining reliability under diverse distributional conditions where the $t$-test fails.

The test operates by comparing all pairs of measurements between the two samples. Given samples $\mathbf{x} = (x_1, \ldots, x_n)$ and $\mathbf{y} = (y_1, \ldots, y_m)$, the Mann-Whitney $U$ statistic counts how many pairs satisfy $x_i > y_j$:

If the samples come from the same distribution, $U$ should be near $n \frac{m}{2}$ (roughly half the pairs favor $\mathbf{x}$, half favor $\mathbf{y}$). Large deviations from $n \frac{m}{2}$ suggest the distributions differ.

The test answers: Could this $U$ value arise by chance if the samples were truly equivalent? The $p$-value quantifies this probability. If $p < 0.05$, traditional practice declares the difference statistically significant.

This approach creates several problems for practitioners:

• Binary thinking. The test produces a yes/no answer: reject or fail to reject the null hypothesis. Practitioners typically want to know the magnitude of difference, not just whether one exists.
• Arbitrary thresholds. The 0.05 threshold has no universal justification, yet it dominates practice and creates a false dichotomy between $p = 0.049$ and $p = 0.051$.
• Hypothesis-centric framework. The test assumes a null hypothesis of no difference and evaluates evidence against it. Real questions rarely concern exact equality; practitioners want to know how different? rather than different or not?
• Inverted logic. The natural question is what shifts are consistent with my data? The test answers is this specific shift (zero) consistent with my data?

The toolkit inverts this framework. Instead of testing whether a hypothesized shift is plausible, we compute which shifts are plausible given the data. This inversion transforms hypothesis testing into bounds estimation.

The mathematical foundation remains the same. The distribution of pairwise comparisons under random sampling determines which order statistics of pairwise differences form reliable bounds. The Mann-Whitney $U$ statistic measures pairwise comparisons ($x_i > y_j$). The $\operatorname{Shift}$ estimator uses pairwise differences ($x_i - y_j$). These quantities are mathematically related: a pairwise difference $x_i - y_j$ is positive exactly when $x_i > y_j$. The toolkit renames this comparison count from $U$ to $\text{Dominance}(\mathbf{x}, \mathbf{y})$, clarifying its purpose: measuring how often one sample dominates the other in pairwise comparisons.

The distribution of $\text{Dominance}$ determines which order statistics form reliable bounds. Define the margin function:

This function computes how many extreme pairwise differences could occur by chance with probability $\mathrm{misrate}$, based on the distribution of pairwise comparisons.

The $\operatorname{PairwiseMargin}$ function requires knowing the distribution of pairwise comparisons under sampling. Two computational approaches exist:

• Exact computation (Löfflers algorithm, 1982). Uses a recurrence relation to compute the exact distribution of pairwise comparisons for small samples. Practical for combined sample sizes up to several hundred.
• Approximation (Edgeworth expansion, 1955). Refines the normal approximation with correction terms based on higher moments of the distribution. Provides accurate results for large samples where exact computation becomes impractical.

The toolkit automatically selects the appropriate method based on sample sizes, ensuring both accuracy and computational efficiency.

This approach naturally complements $\operatorname{Center}$ and $\operatorname{Spread}$:

• $\operatorname{Center}(\mathbf{x})$ uses the median of pairwise averages $\frac{x_i + x_j}{2}$
• $\operatorname{Spread}(\mathbf{x})$ uses the median of pairwise differences $\lvert x_i - x_j \rvert$
• $\operatorname{Shift}(\mathbf{x}, \mathbf{y})$ uses the median of pairwise differences $x_i - y_j$
• $\operatorname{ShiftBounds}(\mathbf{x}, \mathbf{y}, \mathrm{misrate})$ uses order statistics of the same pairwise differences

All procedures build on pairwise operations. This structural consistency reflects the mathematical unity underlying robust statistics: pairwise operations provide natural robustness while maintaining computational feasibility and statistical efficiency.

The inversion from hypothesis testing to bounds estimation represents a philosophical shift in statistical practice. Traditional methods ask should I believe this specific hypothesis? Pragmatic methods ask what should I believe, given this data? Bounds provide actionable answers: they tell practitioners which values are plausible, enabling informed decisions without arbitrary significance thresholds.

Traditional Mann-Whitney implementations apply tie correction when samples contain repeated values. This correction modifies variance calculations to account for tied observations, changing $p$-values and confidence intervals in ways that depend on measurement precision. The toolkit deliberately omits tie correction. Continuous distributions produce theoretically distinct values; observed ties result from finite measurement precision and digital representation. When measurements appear identical, this reflects rounding of underlying continuous variation, not true equality in the measured quantity. Treating ties as artifacts of discretization rather than distributional features simplifies computation while maintaining accuracy. The exact and approximate methods compute comparison distributions without requiring adjustments for tied values, eliminating a source of complexity and potential inconsistency in statistical practice.

Historical Development

The mathematical foundations emerged through decades of refinement. Mann and Whitney (1947) established the distribution of pairwise comparisons under random sampling, creating the theoretical basis for comparing samples through rank-based methods. Their work demonstrated that comparison counts follow predictable patterns regardless of the underlying population distributions.

The original computational approaches suffered from severe limitations. Mann and Whitney proposed a slow exact method requiring exponential resources and a normal approximation that proved grossly inaccurate for practical use. The approximation works reasonably in distribution centers but fails catastrophically in the tails where practitioners most need accuracy. For moderate sample sizes, approximation errors can exceed factors of $10^11$.

Fix and Hodges (1955) addressed the approximation problem through higher-order corrections. Their expansion adds terms based on the distributions actual moments rather than assuming perfect normality. This refinement reduces tail probability errors from orders of magnitude to roughly 1%, making approximation practical for large samples where exact computation becomes infeasible.

Löffler (1982) solved the exact computation problem through algorithmic innovation. The naive recurrence requires quadratic memory— infeasible for samples beyond a few dozen measurements. Löffler discovered a reformulation that reduces memory to linear scale, making exact computation practical for combined sample sizes up to several hundred.

Despite these advances, most statistical software continues using the 1947 approximation. The computational literature contains the solutions, but software implementations lag decades behind theoretical developments. This toolkit implements both the exact method for small samples and the refined approximation for large samples, automatically selecting the appropriate approach based on sample sizes.

The shift from hypothesis testing to bounds estimation requires no new mathematics. The same comparison distributions that enable hypothesis tests also determine which order statistics form reliable bounds. Traditional applications ask is zero plausible? and answer yes or no. This toolkit asks which values are plausible? and answers with an interval. The perspective inverts while the mathematical foundation remains identical.

Notes

On Bootstrap for Center Bounds

A natural question arises: can bootstrap resampling improve $\operatorname{CenterBounds}$ coverage for asymmetric distributions where the weak symmetry assumption fails?

The idea is appealing. The signed-rank approach computes bounds from order statistics of Walsh averages using a margin derived from the Wilcoxon distribution, which assumes symmetric deviations from the center. Bootstrap makes no symmetry assumption: resample the data with replacement, compute $\operatorname{Center}$ on each resample, and take quantiles of the bootstrap distribution as bounds. This should yield valid bounds regardless of distributional shape.

This manual deliberately does not provide a bootstrap-based alternative to $\operatorname{CenterBounds}$. The reasons are both computational and statistical.

Computational cost

$\operatorname{CenterBounds}$ computes bounds in $O(n \log n)$ time: a single pass through the Walsh averages guided by the signed-rank margin. No resampling, no iteration.

A bootstrap version requires $B$ resamples (typically $B = 10000$ for stable tail quantiles), each computing $\operatorname{Center}$ on the resample. $\operatorname{Center}$ itself costs $O(n \log n)$ via the fast selection algorithm on the implicit pairwise matrix. The total cost becomes $O(B \cdot n \log n)$ per call — roughly $10000 \times$ slower than the signed-rank approach.

For $n = 5$, each call to $\operatorname{Center}$ operates on $15$ Walsh averages. The bootstrap recomputes this $10000$ times. The computation is not deep — it is merely wasteful. For $n = 100$, there are $5050$ Walsh averages per resample, and $10000$ resamples produce $5 \times 10^7$ selection operations per bounds call. In a simulation study that evaluates bounds across many samples, this cost becomes prohibitive.

Statistical quality

Bootstrap bounds are nominal, not exact. The percentile method has well-documented undercoverage for small samples: requesting 95% confidence ($\mathrm{misrate} = 0.05$) typically yields 85–92% actual coverage for $n < 30$. This is inherent to the bootstrap percentile method — the quantile estimates from $B$ resamples are biased toward the sample and underrepresent tail behavior. Refined methods (BCa, bootstrap-$t$) partially address this but add complexity and still provide only asymptotic guarantees.

Meanwhile, $\operatorname{CenterBounds}$ provides exact distribution-free coverage under symmetry. For $n = 5$ requesting $\mathrm{misrate} = 0.1$, the signed-rank method delivers exactly $10\%$ misrate. A bootstrap method, requesting the same $10\%$, typically delivers $12$–$15\%$ misrate. The exact method is simultaneously faster and more accurate.

Behavior under asymmetry

Under asymmetric distributions, the signed-rank margin is no longer calibrated: the Wilcoxon distribution assumes symmetric deviations, and asymmetry shifts the actual distribution of comparison counts.

However, the coverage degradation is gradual, not catastrophic. Mild asymmetry produces mild coverage drift. The bounds remain meaningful — they still bracket the pseudomedian using order statistics of the Walsh averages — but the actual misrate differs from the requested value.

This is the same situation as bootstrap, which also provides only approximate coverage. The practical difference is that the signed-rank approach achieves this approximate coverage in $O(n \log n)$ time, while bootstrap achieves comparable approximate coverage in $O(B \cdot n \log n)$ time.

Why not both?

One might argue for providing both methods: the signed-rank approach as default, and a bootstrap variant for cases where symmetry is severely violated.

This creates a misleading choice. If the bootstrap method offered substantially better coverage under asymmetry, the complexity would be justified. But for the distributions practitioners encounter ($\underline{\operatorname{Multiplic}}$, $\underline{\operatorname{Exp}}$, and other moderate asymmetries), the coverage difference between the two approaches is small relative to the $10000 \times$ cost difference. For extreme asymmetries where the signed-rank coverage genuinely breaks down, the sign test provides an alternative foundation for median bounds (see the study on misrate efficiency of MedianBounds), but its $O(n^{-\frac{1}{2}})$ efficiency convergence makes it impractical for moderate sample sizes.

The toolkit therefore provides $\operatorname{CenterBounds}$ as the single bounds estimator. The weak symmetry assumption means the method performs well under approximate symmetry and degrades gracefully under moderate asymmetry. There is no useful middle ground that justifies a $10000 \times$ computational penalty for marginally different approximate coverage.

On Misrate Efficiency of MedianBounds

This study analyzes $MedianBounds$, a bounds estimator for the population median based on the sign test, and explains why pragmastat omits it in favor of $\operatorname{CenterBounds}$.

Definition

where $k$ is the largest integer satisfying $2 \cdot \Pr(B \leq k-1) \leq \mathrm{misrate}$ and $B \sim \text{Binomial}(n, 0.5)$. The interval brackets the population median using order statistics, with $\mathrm{misrate}$ controlling the probability that the true median falls outside the bounds.

$MedianBounds$ requires no symmetry assumption — only weak continuity — making it applicable to arbitrarily skewed distributions. This is its principal advantage over $\operatorname{CenterBounds}$, which assumes weak symmetry.

Sign test foundation

The method is equivalent to inverting the sign test. Under weak continuity, each observation independently falls above or below the true median with probability $\frac{1}{2}$. The number of observations below the median follows $\text{Binomial}(n, 0.5)$, and the order statistics $x_{(k)}$ and $x_{(n-k+1)}$ form a confidence interval whose coverage is determined exactly by the binomial CDF.

Because the binomial CDF is a step function, the achievable misrate values form a discrete set. The algorithm rounds down to the nearest achievable level, inevitably wasting part of the requested misrate budget. This study derives the resulting efficiency loss and its convergence rate.

Achievable misrate levels

The achievable misrate values for sample size $n$ are:

The algorithm selects the largest $k$ satisfying $m_k \leq \mathrm{misrate}$. The efficiency $\eta = \frac{m_k}{\mathrm{misrate}}$ measures how much of the requested budget is used. Efficiency $\eta = 1$ means the bounds are as tight as the misrate allows; $\eta = 0.5$ means half the budget is wasted, producing unnecessarily wide bounds.

Spacing between consecutive levels

The gap between consecutive achievable misrates is:

For a target misrate $\alpha$, the relevant index $k$ satisfies $m_k \approx \alpha$. By the normal approximation to the binomial, $B \approx cal(N)(\frac{n}{2}, \frac{n}{4})$, the binomial CDF near this index changes by approximately:

where $z_\alpha = \Phi^{-1}(\frac{\alpha}{2})$ is the corresponding standard normal quantile and $\phi$ is the standard normal density. This spacing governs how coarsely the achievable misrates are distributed near the target.

Expected efficiency

The requested misrate $\alpha$ falls at a uniformly random position within a gap of width $\Delta m$. On average, the algorithm wastes $\Delta \frac{m}{2}$, giving expected efficiency:

Define the misrate-dependent constant:

Then the expected efficiency has the form:

The convergence rate is $O(n^{-\frac{1}{2}})$: efficiency improves as the square root of sample size.

Values of $C(\alpha)$

The constant $C(\alpha)$ increases for smaller misrates, meaning tighter error tolerances require proportionally larger samples for efficient bounds:

For $\alpha = 0.1$ and $n = 50$: $\mathbb{E}[\eta] \approx 1 - 2.\frac{06}{\sqrt{50}} \approx 0.71$. Achieving $90\%$ efficiency on average requires $n > (\frac{C(\alpha)}{0}.1)^2$. For $\alpha = 0.1$ this gives $n > 424$; for $\alpha = 0.001$ this gives $n > 1190$.

Comparison with CenterBounds

$\operatorname{CenterBounds}$ uses the signed-rank statistic $W$ with range $[0, \frac{n(n+1)}{2}]$. Under the null hypothesis, $W$ has variance $\sigma^2 = n(n+1)\frac{2n+1}{24} \approx n^\frac{3}{12}$. The CDF spacing at the relevant quantile is:

The expected efficiency for $\operatorname{CenterBounds}$ is therefore:

This converges at rate $O(n^{-\frac{3}{2}})$ — three polynomial orders faster than $MedianBounds$. The difference arises because the signed-rank distribution has $\frac{n(n+1)}{2}$ discrete levels compared to the binomials $n$ levels, providing fundamentally finer resolution.

Why pragmastat omits MedianBounds

The efficiency loss of $MedianBounds$ is not an implementation artifact. It reflects a structural limitation of the sign test: using only the signs of $(X_i - \theta)$ discards magnitude information, leaving only $n$ binary observations to determine coverage. The signed-rank test used by $\operatorname{CenterBounds}$ exploits both signs and ranks, producing $\frac{n(n+1)}{2}$ comparison outcomes and correspondingly finer misrate resolution.

For applications requiring tight misrate control on the median, large samples ($n > 500$) are needed to ensure efficient use of the misrate budget. For smaller samples, the bounds remain valid but conservative: the actual misrate is guaranteed to not exceed the requested value, even though it may be substantially below it.

$\operatorname{CenterBounds}$ with its $O(n^{-\frac{3}{2}})$ convergence achieves near-continuous misrate control even for moderate $n$, at the cost of requiring weak symmetry. For the distributions practitioners typically encounter, this tradeoff favors $\operatorname{CenterBounds}$ as the single bounds estimator in the toolkit. When symmetry is severely violated, the coverage drift of $\operatorname{CenterBounds}$ is gradual — mild asymmetry produces mild drift — making it a robust default without the efficiency penalty inherent to the sign test approach.

Additive (Normal) Distribution

The $\underline{\operatorname{Additive}}$ (Normal) distribution has two parameters: the mean and the standard deviation, written as $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$.

Asymptotic Spread Value

Consider two independent draws $X$ and $Y$ from the $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$ distribution. The goal is to find the median of their absolute difference $\lvert X-Y \rvert$. Define the difference $D = X - Y$. By linearity of expectation, $\mathbb{E}[D] = 0$. By independence, $\operatorname{Var}[D] = 2 \cdot \mathrm{stdDev}^2$. Thus $D$ has distribution $\underline{\operatorname{Additive}}(0, \sqrt{2} \cdot \mathrm{stdDev})$, and the problem reduces to finding the median of $\lvert D \rvert$. The location parameter $\mathrm{mean}$ disappears, as expected, because absolute differences are invariant to shifts.

Let $\tau=\sqrt{2} \cdot \mathrm{stdDev}$, so that $D \sim \underline{\operatorname{Additive}}(0, \tau)$. The random variable $\lvert D \rvert$ then follows the Half-$\underline{\operatorname{Additive}}$ (Folded Normal) distribution with scale $\tau$. Its cumulative distribution function for $z \geq 0$ becomes

where $\Phi$ denotes the standard $\underline{\operatorname{Additive}}$ (Normal) CDF.

The median $m$ is the point at which this cdf equals $\frac{1}{2}$. Setting $F_{\lvert D \rvert}(m)=\frac{1}{2}$ gives

Applying the inverse cdf yields $\frac{m}{\tau} = z_{0.75}$. Substituting back $\tau = \sqrt{2} \cdot \mathrm{stdDev}$ produces

Define $z_{0.75} := \Phi^{-1}(0.75) \approx 0.6744897502$. Numerically, the median absolute difference is approximately $\sqrt{2} \cdot z_{0.75} \cdot \mathrm{stdDev} \approx 0.9538725524 \cdot \mathrm{stdDev}$. This expression depends only on the scale parameter $\mathrm{stdDev}$, not on the mean, reflecting the translation invariance of the problem.

Lemma: Average Estimator Drift Formula

For average estimators $T_n$ with asymptotic standard deviation $a \cdot \mathrm{stdDev} / \sqrt{n}$ around the mean $\mu$, define $\operatorname{RelSpread}[T_n] := \frac{\operatorname{Spread}[T_n]}{\operatorname{Spread}[X]}$. In the $\underline{\operatorname{Additive}}$ (Normal) case, $\operatorname{Spread}[X] = \sqrt{2} \cdot z_{0.75} \cdot \mathrm{stdDev}$.

For any average estimator $T_n$ with asymptotic standard deviation $a \cdot \mathrm{stdDev} / \sqrt{n}$ around the mean $\mu$, the drift calculation follows:

• The spread of two independent estimates: $\operatorname{Spread}[T_n] = \sqrt{2} \cdot z_{0.75} \cdot a \cdot \mathrm{stdDev} / \sqrt{n}$
• The relative spread: $\operatorname{RelSpread}[T_n] = \frac{a}{\sqrt{n}}$
• The asymptotic drift: $\operatorname{Drift}(T, X) = a$

Asymptotic Mean Drift

For the sample mean $\operatorname{Mean}(\mathbf{x}) = \frac{1}{n} \sum_{i=1}^n x_i$ applied to samples from $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$, the sampling distribution of $\operatorname{Mean}$ is also additive with mean $\mathrm{mean}$ and standard deviation $\mathrm{stdDev}/\sqrt{n}$.

Using the lemma with $a = 1$ (since the standard deviation is $\mathrm{stdDev}/\sqrt{n}$):

$\operatorname{Mean}$ achieves unit drift under the $\underline{\operatorname{Additive}}$ (Normal) distribution, serving as the natural baseline for comparison. $\operatorname{Mean}$ is the optimal estimator under the $\underline{\operatorname{Additive}}$ (Normal) distribution: no other estimator achieves lower $\operatorname{Drift}$.

Asymptotic Median Drift

For the sample median $\operatorname{Median}(\mathbf{x})$ applied to samples from $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$, the asymptotic sampling distribution of $\operatorname{Median}$ is approximately $\underline{\operatorname{Additive}}$ (Normal) with mean $\mathrm{mean}$ and standard deviation $\sqrt{\frac{\pi}{2}} \cdot \mathrm{stdDev}/\sqrt{n}$.

This result follows from the asymptotic theory of order statistics. For the median of a sample from a continuous distribution with density $f$ and cumulative distribution $F$, the asymptotic variance is $\frac{1}{4n[f(F^{-1}(0.5))]^2}$. For the $\underline{\operatorname{Additive}}$ (Normal) distribution with standard deviation $\mathrm{stdDev}$, the density at the median (which equals the mean) is $\frac{1}{\mathrm{stdDev} \sqrt{2 \pi}}$. Thus the asymptotic variance becomes $\pi \cdot \mathrm{stdDev}^\frac{2}{2n}$.

Using the lemma with $a = \sqrt{\frac{\pi}{2}}$:

Numerically, $\sqrt{\frac{\pi}{2}} \approx 1.2533$, so the median has approximately 25% higher drift than the mean under the $\underline{\operatorname{Additive}}$ (Normal) distribution.

Asymptotic Center Drift

For the sample center $\operatorname{Center}(\mathbf{x}) = \underset{1 \leq i \leq j \leq n}{\operatorname{Median}} \frac{x_i + x_j}{2}$ applied to samples from $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$, its asymptotic sampling distribution must be determined.

The center estimator computes all pairwise averages (including $i=j$) and takes their median. For the $\underline{\operatorname{Additive}}$ (Normal) distribution, asymptotic theory shows that the center estimator is asymptotically $\underline{\operatorname{Additive}}$ (Normal) with mean $\mathrm{mean}$.

The exact asymptotic variance of the center estimator for the $\underline{\operatorname{Additive}}$ (Normal) distribution is:

This gives an asymptotic standard deviation of:

Using the lemma with $a = \sqrt{\frac{\pi}{3}}$:

Numerically, $\sqrt{\frac{\pi}{3}} \approx 1.0233$, so the center estimator achieves a drift very close to 1 under the $\underline{\operatorname{Additive}}$ (Normal) distribution, performing nearly as well as the mean while offering greater robustness to outliers.

Lemma: Dispersion Estimator Drift Formula

For dispersion estimators $T_n$ with asymptotic center $b \cdot \mathrm{stdDev}$ and standard deviation $a \cdot \mathrm{stdDev} / \sqrt{n}$, define $\operatorname{RelSpread}[T_n] := \frac{\operatorname{Spread}[T_n]}{b \cdot \mathrm{stdDev}}$.

For any dispersion estimator $T_n$ with asymptotic distribution $T_n \sim\text{approx} \underline{\operatorname{Additive}}(b \cdot \mathrm{stdDev}, (a \cdot \mathrm{stdDev})^\frac{2}{n})$, the drift calculation follows:

• The spread of two independent estimates: $\operatorname{Spread}[T_n] = \sqrt{2} \cdot z_{0.75} \cdot a \cdot \mathrm{stdDev} / \sqrt{n}$
• The relative spread: $\operatorname{RelSpread}[T_n] = \sqrt{2} \cdot z_{0.75} \cdot \frac{a}{b \sqrt{n}}$
• The asymptotic drift: $\operatorname{Drift}(T, X) = \sqrt{2} \cdot z_{0.75} \cdot \frac{a}{b}$

Note: The $\sqrt{2}$ factor comes from the standard deviation of the difference $D = T_1 - T_2$ of two independent estimates, and the $z_{0.75}$ factor converts this standard deviation to the median absolute difference.

Asymptotic StdDev Drift

For the sample standard deviation $\operatorname{StdDev}(\mathbf{x}) = \sqrt{\frac{1}{n-1} \sum_{i=1}^n (x_i - \operatorname{Mean}(\mathbf{x}))^2}$ applied to samples from $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$, the sampling distribution of $\operatorname{StdDev}$ is approximately $\underline{\operatorname{Additive}}$ (Normal) for large $n$ with mean $\mathrm{stdDev}$ and standard deviation $\mathrm{stdDev}/\sqrt{2n}$.

Applying the lemma with $a = \frac{1}{\sqrt{2}}$ and $b = 1$:

For the dispersion drift, we use the relative spread formula:

Since $\operatorname{Center}[\operatorname{StdDev}(X_{1:n})] \approx \mathrm{stdDev}$ asymptotically:

Therefore:

Numerically, $z_{0.75} \approx 0.67449$.

Asymptotic MAD Drift

For the median absolute deviation $\operatorname{MAD}(\mathbf{x}) = \operatorname{Median}(\lvert x_i - \operatorname{Median}(\mathbf{x}) \rvert)$ applied to samples from $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$, the asymptotic distribution is approximately $\underline{\operatorname{Additive}}$ (Normal).

For the $\underline{\operatorname{Additive}}$ (Normal) distribution, the population MAD equals $z_{0.75} \cdot \mathrm{stdDev}$. The asymptotic standard deviation of the sample MAD is:

where $c_{\mathrm{mad}} \approx 0.78$.

Applying the lemma with $a = c_{\mathrm{mad}}$ and $b = z_{0.75}$:

Since $\operatorname{Center}[\operatorname{MAD}(X_{1:n})] \approx z_{0.75} \cdot \mathrm{stdDev}$ asymptotically:

Therefore:

Numerically, $\sqrt{2} \cdot c_{\mathrm{mad}} \approx \sqrt{2} \cdot 0.78 \approx 1.10$.

Asymptotic Spread Drift

For the sample spread $\operatorname{Spread}(\mathbf{x}) = \underset{1 \leq i < j \leq n}{\operatorname{Median}} \lvert x_i - x_j \rvert$ applied to samples from $\underline{\operatorname{Additive}}(\mathrm{mean}, \mathrm{stdDev})$, the asymptotic distribution is approximately $\underline{\operatorname{Additive}}$ (Normal).

The spread estimator computes all pairwise absolute differences and takes their median. For the $\underline{\operatorname{Additive}}$ (Normal) distribution, the population spread equals $\sqrt{2} \cdot z_{0.75} \cdot \mathrm{stdDev}$ as derived in the Asymptotic Spread Value section.

The asymptotic standard deviation of the sample spread for the $\underline{\operatorname{Additive}}$ (Normal) distribution is:

where $c_{\mathrm{spr}} \approx 0.72$.

Applying the lemma with $a = c_{\mathrm{spr}}$ and $b = \sqrt{2} \cdot z_{0.75}$:

Since $\operatorname{Center}[\operatorname{Spread}(X_{1:n})] \approx \sqrt{2} \cdot z_{0.75} \cdot \mathrm{stdDev}$ asymptotically:

Therefore:

Numerically, $c_{\mathrm{spr}} \approx 0.72$.

Summary

Summary for average estimators:

The squared drift values indicate the sample size adjustment needed when switching estimators. For instance, switching from $\operatorname{Mean}$ to $\operatorname{Median}$ while maintaining the same precision requires increasing the sample size by a factor of $\frac{\pi}{2} \approx 1.571$ (about 57% more observations). Similarly, switching from $\operatorname{Mean}$ to $\operatorname{Center}$ requires only about 5% more observations.

The inverse squared drift (rightmost column) equals the classical statistical efficiency relative to the $\operatorname{Mean}$. The $\operatorname{Mean}$ achieves optimal performance (unit efficiency) for the $\underline{\operatorname{Additive}}$ (Normal) distribution, as expected from classical theory. The $\operatorname{Center}$ maintains 95.5% efficiency while offering greater robustness to outliers, making it an attractive alternative when some contamination is possible. The $\operatorname{Median}$, while most robust, operates at only 63.7% efficiency under purely $\underline{\operatorname{Additive}}$ (Normal) conditions.

Summary for dispersion estimators:

For the $\underline{\operatorname{Additive}}$ (Normal) distribution, the asymptotic drift values reveal the relative precision of different dispersion estimators:

The squared drift values indicate the sample size adjustment needed when switching estimators. For instance, switching from $\operatorname{StdDev}$ to $\operatorname{MAD}$ while maintaining the same precision requires increasing the sample size by a factor of $1.\frac{22}{0}.45 \approx 2.71$ (more than doubling the observations). Similarly, switching from $\operatorname{StdDev}$ to $\operatorname{Spread}$ requires a factor of $0.\frac{52}{0}.45 \approx 1.16$.

The $\operatorname{StdDev}$ achieves optimal performance for the $\underline{\operatorname{Additive}}$ (Normal) distribution. The $\operatorname{MAD}$ requires about 2.7 times more data to match the $\operatorname{StdDev}$ precision while offering greater robustness to outliers. The $\operatorname{Spread}$ requires about 1.16 times more data to match the $\operatorname{StdDev}$ precision under purely $\underline{\operatorname{Additive}}$ (Normal) conditions while maintaining robustness.