Approximate and Asymptotic Confidence Intervals

Let $X_1, X_2, \dots, X_n$ be an independent, identically distributed sample from a distribution $P_\theta$ with unknown parameter $\theta \in \Theta \subseteq \mathbb{R}$. Suppose $T_n = T_n(\mathbf{X})$ is a consistent point estimator of $\theta$ satisfying a Central Limit Theorem:

where $V(\theta) > 0$ is the asymptotic variance function. The quantity $\frac{\sqrt{n}(T_n - \theta)}{\sqrt{V(\theta)}}$ converges in distribution to $\mathcal{N}(0, 1)$, but it cannot be directly used as an empirical pivot if $V(\theta)$ depends on the unknown parameter $\theta$.

Slutsky's Theorem and Asymptotic Pivots

To eliminate the unknown parameter from the denominator, we replace $V(\theta)$ with a consistent estimator $\hat{V}_n = V(T_n)$ or an empirical variance estimator such that $\hat{V}_n \xrightarrow{P} V(\theta)$. By Slutsky's Theorem: If $U_n = \sqrt{n}(T_n - \theta) \xrightarrow{d} \mathcal{N}(0, V(\theta))$ and $W_n = \sqrt{\hat{V}_n} \xrightarrow{P} \sqrt{V(\theta)}$, then the ratio:

The studentized ratio $Z_n$ is an asymptotically pivotal quantity. For any pre-assigned confidence level $1 - \alpha$, let $z_{1 - \alpha/2}$ denote the $(1 - \alpha/2)$-quantile of the standard normal distribution:

Inverting the inequality yields the canonical large-sample confidence interval:

The coverage probability of this interval converges to $1 - \alpha$ as $n \to \infty$. However, in finite samples, the true coverage probability $C_n(\theta) = \mathbb{P}_\theta(L_n \le \theta \le U_n)$ deviates from $1 - \alpha$ by an approximation error of order $\mathcal{O}(1/\sqrt{n})$ or $\mathcal{O}(1/n)$.

Edgeworth Expansions and Cornish–Fisher Asymptotic Inversions

The Central Limit Theorem guarantees that the standardized studentized statistic $T_n = \frac{\sqrt{n}(\hat{\theta} - \theta)}{\hat{\sigma}}$ converges in distribution to $\mathcal{N}(0, 1)$ at rate $\mathcal{O}(1/\sqrt{n})$:

To achieve higher-order coverage accuracy, mathematical statisticians employ Edgeworth Expansions, which expand the true CDF $F_n(x)$ in terms of the Hermite polynomials $H_k(x)$ and population cumulants:

where $\kappa_3 = \frac{\mathbb{E}[(X - \mu)^3]}{\sigma^3}$ is the skewness coefficient, and $H_2(x) = x^2 - 1$, $H_3(x) = x^3 - 3x$, $H_5(x) = x^5 - 10x^3 + 15x$. Notice that the leading error term $-\phi(x) \frac{\kappa_3}{6\sqrt{n}} (x^2 - 1)$ is driven entirely by population skewness! If the population is skewed ($\kappa_3 \ne 0$), the true distribution of $T_n$ is asymmetric, and a naive Gaussian interval $T_n \pm z_{1 - \alpha/2}$ experiences asymmetric tail coverage error of order $\mathcal{O}(1/\sqrt{n})$. Inverting the Edgeworth expansion yields the Cornish–Fisher Expansion for the exact quantile $w_\alpha$:

Using Cornish–Fisher quantiles replaces the standard $z_{\alpha}$ cutoffs with skewness-adjusted critical values, reducing interval coverage error from $\mathcal{O}(1/\sqrt{n})$ down to $\mathcal{O}(1/n)$ (second-order accuracy).

The most ubiquitous asymptotic interval in applied statistics is the Wald Confidence Interval, which evaluates the standard error directly at the Maximum Likelihood Estimate $\hat{\theta}_n$.

Definition of the Wald Interval

Let $\hat{\theta}_n$ be the Maximum Likelihood Estimator of $\theta$. By Theorem 5.3:

Estimating the standard error using the observed Fisher Information $J_n(\hat{\theta}) = -\ell_n''(\hat{\theta})$:

The Wald $(1 - \alpha)$ Confidence Interval is:

The Pathology of Parameterization Sensitivity

While computationally effortless, the Wald confidence interval suffers from a notorious mathematical flaw: it is fundamentally not invariant under non-linear reparameterizations! Suppose an investigator reparameterizes the model using a smooth, strictly monotonic transformation $\psi = g(\theta)$. By the invariance property of MLEs, the MLE of $\psi$ is $\hat{\psi} = g(\hat{\theta})$. By the Delta Method, the standard error of $\hat{\psi}$ is:

The Wald interval for $\psi$ is:

Now, invert this interval back into the original parameter space $\theta = g^{-1}(\psi)$:

In general, $I_\theta^* \ne I_\theta$! Different scientific teams analyzing the exact same dataset can obtain radically different confidence intervals simply by choosing whether to parameterize their likelihood in terms of variance $\sigma^2$, standard deviation $\sigma$, precision $\tau = 1/\sigma^2$, or log-variance $\ln(\sigma^2)$. Furthermore, if the log-likelihood surface is asymmetric or non-quadratic, the Wald interval can extend beyond natural parameter boundaries, producing absurd negative variance estimates or probabilities exceeding $1$.

To overcome the parameterization fragility of the Wald interval, C. R. Rao and Edwin B. Wilson proposed inverting the Score Statistic instead of the point estimator.

The Score Statistic Inversion

Recall the Score function $S_n(\theta) = \frac{\partial \ell_n(\theta)}{\partial \theta}$ and expected Fisher Information $I_n(\theta)$. By Bartlett's first identity and the Central Limit Theorem:

The Score Confidence Region gathers all candidate parameter values $\theta$ whose score statistic is within the normal critical threshold:

Key Mathematical Advantage over Wald

Notice the profound structural difference between the Score and Wald methods:

• In the Wald method, the denominator $\text{SE}(\hat{\theta})$ is evaluated at the single point estimate $\hat{\theta}$ and held frozen as $\theta$ varies.
• In the Score method, the denominator $\sqrt{I_n(\theta)}$ is evaluated dynamically at the hypothesized parameter $\theta$!

Evaluating information at $\theta$ automatically adjusts the standard error to the local curvature of the likelihood across the entire parameter space. For asymmetric distributions, inverting the quadratic score equation naturally produces asymmetric intervals that cannot violate physical boundary constraints.

The third member of the likelihood trinity is the Likelihood Ratio (LR) Confidence Interval, obtained by inverting the Likelihood Ratio Test statistic.

Theorem: Inversion of Wilks' Likelihood Ratio Test

Let $\ell_n(\theta)$ be the log-likelihood function and let $\hat{\theta}_n$ be the unconstrained MLE. By Wilks' Theorem (1938), under true parameter $\theta_0$:

The $100(1 - \alpha)\%$ Likelihood Ratio Confidence Region is defined as the set of all parameter values $\theta$ whose drop in log-likelihood does not exceed the critical Chi-square threshold:

Equivalently, in terms of the raw likelihood function:

For a $95\%$ confidence interval, $z_{0.975} = 1.96$, so $z^2 = 3.8415$ and $\frac{z^2}{2} = 1.9207$. The $95\%$ LR interval consists of all parameter values where the log-likelihood is within $1.92$ units of its maximum peak:

Invariance of Likelihood Ratio Intervals

Unlike the Wald interval, the Likelihood Ratio confidence interval is strictly invariant under one-to-one transformations! If $\psi = g(\theta)$ is any monotonic reparameterization, the likelihood value satisfies $L^*(g(\theta)) = L(\theta)$ and $L^*(g(\hat{\theta})) = L(\hat{\theta})$. Therefore:

The Likelihood Ratio interval for $\psi$ is precisely the direct mapping of the interval for $\theta$:

This functional invariance makes the Likelihood Ratio interval the gold standard in non-linear modeling, dose-response analysis, and survival analysis.

Asymptotic Confidence Regions in Multiparameter Models: Profile Deviance Contours

In multidimensional models $\boldsymbol{\theta} = (\theta_1, \dots, \theta_p)^T$, establishing simultaneous confidence sets for a subset of parameters $\boldsymbol{\psi} = (\theta_1, \dots, \theta_r)^T$ in the presence of nuisance parameters $\boldsymbol{\lambda} = (\theta_{r+1}, \dots, \theta_p)^T$ relies on the Profile Deviance Statistic:

By Wilks' Theorem for partitioned parameters:

The $100(1 - \alpha)\%$ simultaneous confidence region for $\boldsymbol{\psi}$ is the contour set:

Geometric Characteristics:

1. Curved and Asymmetric Topography: Unlike the elliptical contours produced by the multivariate Wald method $(\hat{\boldsymbol{\psi}} - \boldsymbol{\psi})^T [\widehat{\text{Var}}(\hat{\boldsymbol{\psi}})]^{-1} (\hat{\boldsymbol{\psi}} - \boldsymbol{\psi}) \le \chi_r^2$, profile deviance contours follow the exact non-linear ridges and valleys of the likelihood surface.

2. Banana-Shaped Contours: In non-linear models (such as Michaelis–Menten enzyme kinetics or compartment models in pharmacokinetics), the parameters exhibit strong non-linear collinearity. The profile deviance contours naturally bend into banana-shaped or crescent manifolds, capturing parameter trade-offs that quadratic Wald approximations misrepresent.

3. Invariance: The profile deviance region is invariant under any invertible transformation $\boldsymbol{\phi} = \mathbf{g}(\boldsymbol{\psi})$.

The estimation of an unknown success probability $p \in (0, 1)$ based on $X \sim \text{Binomial}(n, p)$ is one of the most widely used procedures in science. However, the classical textbook interval exhibits catastrophic coverage failures.

1. The Classical Wald Interval (The "Standard" Method)

The MLE is the sample proportion $\hat{p} = X / n$. Its estimated standard error is $\sqrt{\frac{\hat{p}(1 - \hat{p})}{n}}$. The classical Wald interval is:

Catastrophic Flaws of the Wald Interval (Brown, Cai, and DasGupta, 2001):

• Zero-Width Failure: If $X = 0$ or $X = n$, then $\hat{p}(1 - \hat{p}) = 0$, producing a degenerate confidence interval of length zero $[0, 0]$ or $[1, 1]$, which asserts absolute certainty!
• Oscillation and Systematic Under-Coverage: Because the binomial distribution is discrete, the true coverage probability $C_n(p) = \sum_{x=0}^n \binom{n}{x} p^x (1 - p)^{n-x} \mathbf{1}_{\{L(x) \le p \le U(x)\}}$ oscillates wildly and drops far below the nominal $95\%$ target—frequently falling to $80\%$ or lower, even for sample sizes as large as $n = 100$!

2. The Wilson Score Interval (1927)

Inverting the Score statistic $\frac{(\hat{p} - p)^2}{p(1 - p) / n} \le z^2$:

This is a quadratic inequality in $p$:

Solving for the roots via the quadratic formula yields the celebrated Wilson Score Interval:

where the center of the interval is the shrinkage proportion:

For $95\%$ confidence ($z \approx 2$), the center $\tilde{p} \approx \frac{X + 2}{n + 4}$! This is the famous Agresti–Coull "Add-Two-Successes-and-Two-Failures" rule: by adding two successes and two failures to the sample, the simple Wald formula produces an interval whose coverage matches the exact Wilson score interval.

3. The Clopper–Pearson Exact Interval (1934)

To guarantee that coverage never falls below $1 - \alpha$ for any $p$, C. J. Clopper and E. S. Pearson inverted the exact binomial cumulative distribution function:

Using the identity between the incomplete Beta integral and the binomial sum, the exact endpoints can be computed via Beta quantiles:

The Clopper–Pearson interval is strictly conservative ($C_n(p) \ge 1 - \alpha$ for all $p$), making it the standard choice in regulatory safety applications where under-coverage is unacceptable.

Asymptotic Confidence Intervals for Odds Ratios and Relative Risk in 2x2 Tables

In clinical epidemiology and medical genetics, the association between an exposure (or treatment) and a binary health outcome is summarized in a $2 \times 2$ contingency table:

Two primary metrics quantify the association:

1. Relative Risk (Risk Ratio, RR):

2. Odds Ratio (OR):

Both $\widehat{\text{RR}}$ and $\widehat{\text{OR}}$ are strictly non-negative ratios whose sampling distributions are severely skewed to the right for small to moderate sample sizes. Direct Wald intervals on the raw ratios produce catastrophic under-coverage and can violate the non-negativity bound. To achieve rapid Gaussian convergence, statisticians apply the logarithmic transformation:

By the Multivariate Delta Method, using the multinomial variance-covariance structure, the asymptotic standard error of the log-odds ratio is given by Woolf's Formula (1955):

The asymptotic $100(1 - \alpha)\%$ confidence interval for $\ln(\text{OR})$ is:

Exponentiating both endpoints yields the exact multiplicative confidence interval for the Odds Ratio:

Similarly, for the Relative Risk, the standard error of $\ln(\widehat{\text{RR}})$ is:

Exponentiating yields an asymmetric, range-preserving confidence interval for the Relative Risk that naturally respects the $(0, \infty)$ domain.

For an i.i.d. sample $X_1, \dots, X_n \sim \text{Poisson}(\lambda)$, the sample mean is $\bar{X}$ and the sum is $T = \sum X_i \sim \text{Poisson}(n\lambda)$.

1. Wald Interval:

Fails when $\bar{X} = 0$, producing $[0, 0]$.

2. Score Interval:

Inverting $\frac{(\bar{X} - \lambda)^2}{\lambda / n} \le z^2 \iff n \lambda^2 - (2n\bar{X} + z^2)\lambda + n\bar{X}^2 \le 0$:

When $\bar{X} = 0$, the score interval smoothly evaluates to $\left[ 0, \, \frac{z^2}{n} \right]$, providing a non-trivial, valid upper bound.

3. Exact Garwood Interval:

Using the relationship between the Poisson sum and the Chi-square distribution:

This guarantees exact coverage $\ge 1 - \alpha$ for all sample sizes.

When the asymptotic variance depends on the unknown parameter $\theta$, standard normal approximations can be accelerated by transforming $\hat{\theta}$ via a variance-stabilizing transformation $g(\theta) = \int \sqrt{I_1(\theta)} \, d\theta$.

General Procedure:

1. Identify the transformation $g$ such that $\text{Var}(g(\hat{\theta})) \approx \frac{c^2}{n}$, a constant independent of $\theta$.
2. Construct the symmetric Gaussian interval for the transformed parameter $\psi = g(\theta)$:

1. Invert the endpoints back using $g^{-1}$:

Classical Examples:

• Binomial Proportion ($g(p) = \arcsin(\sqrt{p})$):

Inverting yields endpoints bounded strictly inside $[0, 1]$.

• Poisson Rate ($g(\lambda) = \sqrt{\lambda}$):

• Correlation Coefficient ($g(\rho) = \frac{1}{2}\ln\frac{1+\rho}{1-\rho} = \text{arctanh}(\rho)$):

Inverting via $\tanh$ guarantees endpoints strictly inside $[-1, 1]$.

Woodruff's Method for Nonparametric Quantile Confidence Intervals

In econometric surveys of household income and wealth, distributions are heavily skewed and Pareto-tailed, making mean estimation unstable. Policy analysts rely instead on population quantiles $\xi_p = F^{-1}(p)$ (such as the median $\xi_{0.5}$ or 90th percentile $\xi_{0.9}$). Ralph S. Woodruff (1952) established an ingenious, distribution-free method for constructing confidence intervals for any quantile $\xi_p$ without requiring density estimation. Let $\hat{F}_n(t) = \frac{1}{n} \sum_{i=1}^n \mathbf{1}_{\{X_i \le t\}}$ be the empirical cumulative distribution function. At the true quantile $\xi_p$, the indicator $\mathbf{1}_{\{X_i \le \xi_p\}}$ is a Bernoulli variable with success probability $p$. By the Central Limit Theorem:

A $100(1 - \alpha)\%$ confidence interval for the probability value $p$ is:

Woodruff inverted this vertical probability interval horizontally through the sample quantile function $\hat{F}_n^{-1}$:

For an ordered sample $X_{(1)} \le X_{(2)} \le \dots \le X_{(n)}$, the endpoints are simply the order statistics:

where the integer ranks are:

Woodruff's interval requires zero kernel density estimation, avoids bandwidth selection dilemmas, and guarantees exact non-parametric coverage across all continuous distribution shapes.

When analytical distributions cannot be derived, resampling methods provide distribution-free confidence intervals.

1. The Standard Normal Bootstrap Interval

Assumes the bootstrap distribution is symmetric and approximately normal.

2. The Percentile Bootstrap Interval (Efron, 1981)

Sort the $B$ bootstrap estimates in ascending order: $\hat{\theta}^{*(1)} \le \hat{\theta}^{*(2)} \le \dots \le \hat{\theta}^{*(B)}$. The $100(1 - \alpha)\%$ percentile interval is defined by the empirical quantiles:

The percentile interval is range-preserving and transformation-invariant, but exhibits coverage error of order $\mathcal{O}(1/\sqrt{n})$ if the sampling distribution is skewed.

3. The Bias-Corrected and Accelerated (BCa) Interval

To achieve second-order accuracy (coverage error $\mathcal{O}(1/n)$), Bradley Efron (1987) introduced the BCa interval, which adjusts the quantile endpoints for both median bias ($z_0$) and skewness/acceleration ($a$):

where:

• The bias-correction parameter $z_0 = \Phi^{-1}\left( \frac{1}{B} \sum_{b=1}^B \mathbf{1}_{\{\hat{\theta}^{*b} < \hat{\theta}\}} \right)$ measures discrepancy between the bootstrap median and the point estimate.
• The acceleration parameter $a = \frac{\sum_{i=1}^n (\bar{\theta}_{(\cdot)} - \hat{\theta}_{(-i)})^3}{6 [\sum_{i=1}^n (\bar{\theta}_{(\cdot)} - \hat{\theta}_{(-i)})^2]^{3/2}}$ is computed from jackknife influence values.

The BCa interval automatically corrects for skewness and parameter-dependent variance, achieving second-order asymptotic correctness without requiring analytical Edgeworth corrections.

When evaluating the quality of an asymptotic confidence interval, mathematicians distinguish between first-order accurate and second-order accurate procedures:

• A confidence interval $[L_n, U_n]$ is said to be first-order accurate if its one-sided coverage errors decay at rate $\mathcal{O}(1/\sqrt{n})$:

The standard Wald interval $\hat{\theta} \pm z_{1-\alpha/2} \text{SE}(\hat{\theta})$ and standard normal bootstrap intervals are only first-order accurate.

• An interval is second-order accurate if the $\mathcal{O}(1/\sqrt{n})$ term vanishes identically across all $\theta \in \Theta$:

Second-order accuracy means that the interval achieves tenfold higher coverage precision in practical sample sizes. The Wilson Score interval, the Likelihood Ratio interval, and the BCa bootstrap interval are all second-order accurate.

Bartlett Correction for Likelihood Ratio Intervals

In 1937, M. S. Bartlett discovered that the expectation of the likelihood ratio statistic $w = 2[\ell_n(\hat{\theta}) - \ell_n(\theta_0)]$ satisfies:

where $b$ is a constant depending on the third and fourth log-likelihood cumulants. Bartlett showed that scaling the likelihood ratio statistic by its expected bias factor:

accelerates the Chi-square convergence rate from $\mathcal{O}(1/n)$ down to $\mathcal{O}(1/n^2)$! Inverting the Bartlett-corrected likelihood ratio statistic:

yields confidence intervals whose coverage errors decay at the rapid rate $\mathcal{O}(1/n^2)$, providing near-exact coverage even in tiny samples ($n = 10$).

Can a Bayesian credible interval achieve exact frequentist coverage? In 1963, B. L. Welch and H. W. Peers investigated what prior distributions $\pi(\theta)$ generate Bayesian posterior credible intervals whose frequentist coverage probability matches the nominal level $1 - \alpha$ to higher orders of asymptotic approximation:

Welch and Peers proved that in a single-parameter regular model, the unique prior that eliminates the $\mathcal{O}(1/\sqrt{n})$ frequentist coverage bias is Jeffreys' invariant prior:

Under Jeffreys' prior, the Bayesian $(1 - \alpha)$ one-sided credible limit is second-order accurate: its frequentist coverage error decays at rate $\mathcal{O}(1/n)$ rather than the usual $\mathcal{O}(1/\sqrt{n})$! This celebrated result, known as the Welch–Peers Theorem, establishes a deep and satisfying mathematical bridge between the two competing philosophies: objective Bayesian credible intervals constructed using Jeffreys' prior are simultaneously the most accurate frequentist confidence intervals conceivable.

In 1988, Art B. Owen invented Empirical Likelihood (EL), a revolutionary inference method that combines the non-parametric distribution-free flexibility of the bootstrap with the likelihood-ratio elegance and Bartlett correctability of parametric models.

Definition of Empirical Likelihood

Let $X_1, X_2, \dots, X_n$ be an i.i.d. continuous sample. Without specifying any parametric density family, we place probability weights $p_i = \mathbb{P}(X = X_i)$ on each observed data point, subject to $p_i \ge 0$ and $\sum_{i=1}^n p_i = 1$. The non-parametric likelihood function is:

Unconstrained, $L(F)$ is maximized by the empirical distribution $p_i = 1/n$, yielding maximum value $(1/n)^n$. Now, consider estimating a general parameter $\boldsymbol{\theta}$ defined through an estimating equation $\mathbb{E}[m(X; \boldsymbol{\theta})] = \mathbf{0}$ (such as the mean $m(X; \mu) = X - \mu$). The Empirical Profile Likelihood Ratio is:

Using Lagrange multipliers, the optimal weights are $p_i = \frac{1}{n} \frac{1}{1 + \boldsymbol{\lambda}^T m(X_i; \boldsymbol{\theta})}$, where $\boldsymbol{\lambda}$ is the Lagrange multiplier vector satisfying $\sum_{i=1}^n \frac{m(X_i; \boldsymbol{\theta})}{1 + \boldsymbol{\lambda}^T m(X_i; \boldsymbol{\theta})} = \mathbf{0}$.

Theorem: Nonparametric Wilks' Theorem (Owen, 1988)

Under the true parameter value $\boldsymbol{\theta}_0$:

The $100(1 - \alpha)\%$ Empirical Likelihood Confidence Region is:

Superb Properties of Empirical Likelihood:

1. Completely Non-Parametric: Requires zero parametric assumptions about the underlying distribution family (no normality required).

2. Data-Determined Shape: EL confidence regions naturally adapt to the skewness and curvature of the data. Unlike symmetric Wald intervals, EL contours naturally deform into asymmetric, data-driven shapes.

3. Range Preserving: EL intervals for bounded parameters (e.g. positive variances) automatically respect physical boundaries without ad-hoc truncation.

4. Bartlett Correctable: Peter Hall (1990) proved that Owen's empirical likelihood is Bartlett-correctable, reducing coverage error from $\mathcal{O}(1/n)$ to $\mathcal{O}(1/n^2)$.

Higher-Order Tail Calibrations and Higher-Order Asymptotics of Confidence Intervals

While standard first-order asymptotic confidence intervals rely upon the normal approximation $Z = \frac{\hat{\theta} - \theta}{\text{SE}(\hat{\theta})} \approx \mathcal{N}(0, 1)$, empirical coverage probabilities in moderate samples ($n \in [20, 100]$) frequently exhibit substantial coverage errors of order $\mathcal{O}(1/\sqrt{n})$ due to skewness in the underlying distribution.

To eliminate these asymmetric tail coverage errors, higher-order asymptotic theory develops Cornish–Fisher inversions and Edgeworth expansions of the standardized pivotal distribution:

where $\kappa_3$ and $\kappa_4$ denote the standardized skewness and kurtosis of the score function increments. Inverting this expansion analytically gives rise to skewness-corrected confidence bounds, where the upper and lower critical values are modified by asymmetric perturbation terms proportional to $\hat{\kappa}_3 / (6\sqrt{n})$.

Furthermore, Bartlett corrections apply directly to the likelihood ratio statistic $w(\theta) = 2(\ell(\hat{\theta}) - \ell(\theta))$. By dividing $w(\theta)$ by its expected value scale factor $1 + b/n$, where $b$ is computed from joint moments of the log-likelihood derivatives:

The resulting Bartlett-corrected likelihood ratio confidence interval achieves an exceptional coverage error rate of order $\mathcal{O}(1/n^2)$, outperforming ordinary Wald intervals by two entire orders of magnitude. This remarkable error reduction is achieved without requiring any computer-intensive resampling.