Generalized Likelihood Ratio Tests (GLRT), Wilks' Theorem & the Asymptotic Trinity

While Uniformly Most Powerful (UMP), Uniformly Most Powerful Unbiased (UMPU), and Uniformly Most Powerful Invariant (UMPI) tests provide exact finite-sample optimality within exponential families and group-symmetric models, many complex multiparameter models in econometrics, biostatistics, and physical sciences do not admit complete boundary sufficient statistics or simple transformation groups. In 1928, Jerzy Neyman and Egon Pearson introduced a universal constructive methodology for testing arbitrary composite hypotheses: the Generalized Likelihood Ratio Test (GLRT) (often called simply the Likelihood Ratio Test, LRT).

Definition of the Generalized Likelihood Ratio Statistic $\Lambda(\mathbf{x})$

Let $\mathbf{X} = (X_1, \dots, X_n)^T$ have joint probability density or mass function $f(\mathbf{x}; \boldsymbol{\theta})$ indexed by a $k$-dimensional parameter vector $\boldsymbol{\theta} = (\theta_1, \dots, \theta_k)^T \in \Theta \subseteq \mathbb{R}^k$, and denote the likelihood function by $L(\boldsymbol{\theta}; \mathbf{x}) = f(\mathbf{x}; \boldsymbol{\theta})$. Consider testing the general composite null hypothesis against the composite alternative:

Instead of comparing two fixed parameter points as in the simple Neyman-Pearson lemma, the Generalized Likelihood Ratio statistic $\Lambda(\mathbf{x})$ compares the supremum of the likelihood function restricted to the null parameter space $\Theta_0$ against the unrestricted supremum of the likelihood function across the entire parameter space $\Theta = \Theta_0 \cup \Theta_1$:

where:

• $\hat{\boldsymbol{\theta}}_0(\mathbf{x}) = \arg\max_{\boldsymbol{\theta} \in \Theta_0} L(\boldsymbol{\theta}; \mathbf{x})$ is the restricted Maximum Likelihood Estimator (MLE) computed under the null hypothesis constraints $\boldsymbol{\theta} \in \Theta_0$.
• $\hat{\boldsymbol{\theta}}(\mathbf{x}) = \arg\max_{\boldsymbol{\theta} \in \Theta} L(\boldsymbol{\theta}; \mathbf{x})$ is the unrestricted Maximum Likelihood Estimator (MLE) computed over the full parameter space $\Theta$.

Fundamental Analytical Properties of $\Lambda(\mathbf{x})$

1. Bounded Range $[0, 1]$ and Rejection Direction: Because the null subset $\Theta_0$ is contained within the full parameter space $\Theta$ ($\Theta_0 \subset \Theta$), the supremum of a non-negative function over a smaller domain can never exceed its supremum over a larger domain. Therefore, for every sample realization $\mathbf{x} \in \mathcal{X}$:

When the data $\mathbf{x}$ are highly consistent with the null hypothesis $\Theta_0$, the unrestricted MLE $\hat{\boldsymbol{\theta}}(\mathbf{x})$ falls inside or very close to $\Theta_0$, making $L(\hat{\boldsymbol{\theta}}_0; \mathbf{x}) \approx L(\hat{\boldsymbol{\theta}}; \mathbf{x})$ and $\Lambda(\mathbf{x}) \approx 1$. Conversely, when the data strongly favor the alternative space $\Theta_1$, the maximum likelihood achievable inside $\Theta_0$ is orders of magnitude smaller than the unrestricted maximum likelihood $L(\hat{\boldsymbol{\theta}}; \mathbf{x})$, driving $\Lambda(\mathbf{x})$ toward $0$. Consequently, a Generalized Likelihood Ratio Test of level $\alpha$ always rejects $H_0$ for small values of $\Lambda(\mathbf{x})$:

where the critical constant $\lambda_\alpha \in (0, 1)$ is calibrated to satisfy the size constraint $\sup_{\boldsymbol{\theta} \in \Theta_0} \mathbb{P}_{\boldsymbol{\theta}}\left(\Lambda(\mathbf{X}) \le \lambda_\alpha\right) = \alpha$.

2. Sufficiency Factorization: If $\mathbf{T}(\mathbf{X})$ is a sufficient statistic for $\boldsymbol{\theta}$, then by the Neyman-Fisher Factorization Theorem $L(\boldsymbol{\theta}; \mathbf{x}) = g(\mathbf{T}(\mathbf{x}); \boldsymbol{\theta}) h(\mathbf{x})$ with $h(\mathbf{x}) > 0$. Substituting this factorization into $\Lambda(\mathbf{x})$, the parameter-free factor $h(\mathbf{x})$ cancels identically between the numerator and denominator, proving that $\Lambda(\mathbf{x})$ depends on the sample exclusively through the sufficient statistic $\mathbf{T}(\mathbf{x})$.

3. Reparameterization Invariance: By the functional invariance of Maximum Likelihood Estimators, if $\boldsymbol{\psi} = \mathbf{q}(\boldsymbol{\theta})$ is any bijective reparameterization of $\Theta$, the supremum of the likelihood over $\Psi_0 = \mathbf{q}(\Theta_0)$ and $\Psi = \mathbf{q}(\Theta)$ is completely unchanged. Thus, $\Lambda(\mathbf{x})$ is strictly invariant to coordinate reparameterizations.

Profile Likelihood Formulation in the Presence of Nuisance Parameters

In multiparameter statistical models, the parameter vector typically partitions into a vector of interest $\boldsymbol{\vartheta} \in \mathbb{R}^r$ and a vector of nuisance parameters $\boldsymbol{\eta} \in \mathbb{R}^{k-r}$, so that $\boldsymbol{\theta} = (\boldsymbol{\vartheta}, \boldsymbol{\eta})$. When testing a point null hypothesis on the parameters of interest, $H_0 : \boldsymbol{\vartheta} = \boldsymbol{\vartheta}_0$ versus $H_1 : \boldsymbol{\vartheta} \neq \boldsymbol{\vartheta}_0$ while leaving the nuisance vector $\boldsymbol{\eta}$ completely unconstrained, maximizing the likelihood over $\Theta_0 = \{(\boldsymbol{\vartheta}_0, \boldsymbol{\eta}) : \boldsymbol{\eta} \in \Omega\}$ is equivalent to evaluating the profile likelihood function $L_P(\boldsymbol{\vartheta}; \mathbf{x})$:

where $\hat{\boldsymbol{\eta}}(\boldsymbol{\vartheta})$ denotes the conditional Maximum Likelihood Estimator of the nuisance vector $\boldsymbol{\eta}$ for a fixed value of $\boldsymbol{\vartheta}$. Because the unrestricted supremum of $L(\boldsymbol{\vartheta}, \boldsymbol{\eta}; \mathbf{x})$ over $(\boldsymbol{\vartheta}, \boldsymbol{\eta}) \in \Theta$ coincides with the supremum of the profile likelihood $L_P(\boldsymbol{\vartheta}; \mathbf{x})$ over $\boldsymbol{\vartheta}$ (attained at $\hat{\boldsymbol{\vartheta}}$), the Generalized Likelihood Ratio statistic reduces directly to the profile likelihood ratio:

By the Test-Confidence Set Inversion Duality established in Unit 1, inverting the family of profile likelihood ratio tests for all candidate values $\boldsymbol{\vartheta}_0$ generates the profile likelihood confidence region $\{\boldsymbol{\vartheta} : -2 \ln [L_P(\boldsymbol{\vartheta}; \mathbf{x}) / L_P(\hat{\boldsymbol{\vartheta}}; \mathbf{x})] \le \chi^2_{r, 1-\alpha}\}$. Unlike Wald confidence intervals that force symmetric elliptical boundaries around $\hat{\boldsymbol{\vartheta}}$, profile likelihood confidence sets automatically adapt to skewness, nonlinear parameter curvature, and natural domain boundaries.

In Gaussian linear models, the Generalized Likelihood Ratio statistic $\Lambda(\mathbf{X})$ can be evaluated in closed form and shown to be a strictly monotone function of the classical UMPU and UMPI $t$- and $F$-statistics.

1. Exact GLRT for the One-Sample Normal Mean (Unknown Variance)

Let $X_1, \dots, X_n \overset{\text{iid}}{\sim} \mathcal{N}(\mu, \sigma^2)$ with $\boldsymbol{\theta} = (\mu, \sigma^2) \in \Theta = \mathbb{R} \times (0, \infty)$, and consider testing the composite null hypothesis $H_0 : \mu = \mu_0$ ($\Theta_0 = \{(\mu_0, \sigma^2) : \sigma^2 > 0\}$) against the two-sided alternative $H_1 : \mu \neq \mu_0$. The joint likelihood function is:

• Maximization over the Null Space $\Theta_0$: Fixing $\mu = \mu_0$ and differentiating $\ln L(\mu_0, \sigma^2)$ with respect to $\sigma^2$ yields the restricted MLE $\hat{\sigma}_0^2 = \frac{1}{n}\sum_{i=1}^n (x_i - \mu_0)^2$. Substituting $(\mu_0, \hat{\sigma}_0^2)$ back into $L$ gives the maximized null likelihood:

• Maximization over the Full Parameter Space $\Theta$: Maximizing over both $\mu \in \mathbb{R}$ and $\sigma^2 > 0$ yields the unrestricted MLEs $\hat{\mu} = \bar{x}$ and $\hat{\sigma}^2 = \frac{1}{n}\sum_{i=1}^n (x_i - \bar{x})^2$, with maximized full likelihood:

Taking the ratio and applying the orthogonal sum-of-squares decomposition $\sum_{i=1}^n (x_i - \mu_0)^2 = \sum_{i=1}^n (x_i - \bar{x})^2 + n(\bar{x} - \mu_0)^2$:

where $t = \frac{\sqrt{n}(\bar{x} - \mu_0)}{s}$ is Student's $t$-statistic! Because $u \mapsto \left(1 + \frac{u}{n-1}\right)^{-n/2}$ is a strictly decreasing function of $u = t^2 \ge 0$, the GLRT rejection region $\{\Lambda(\mathbf{X}) \le \lambda_\alpha\}$ is mathematically identical to $\{|t(\mathbf{X})| \ge t_{n-1, 1-\alpha/2}\}$, reproducing the exact UMPU and UMPI two-sided Student's $t$-test!

2. Two-Sample Pooled $t$-Test and the Behrens-Fisher Welch-Satterthwaite Approximation

Now suppose $X_1, \dots, X_{n_1} \overset{\text{iid}}{\sim} \mathcal{N}(\mu_1, \sigma^2)$ and $Y_1, \dots, Y_{n_2} \overset{\text{iid}}{\sim} \mathcal{N}(\mu_2, \sigma^2)$ are independent samples sharing a common unknown variance $\sigma^2 > 0$. To test $H_0 : \mu_1 = \mu_2$ versus $H_1 : \mu_1 \neq \mu_2$, let $\nu = n_1 + n_2 - 2$ and let $S_p^2 = \frac{(n_1 - 1)S_X^2 + (n_2 - 1)S_Y^2}{n_1 + n_2 - 2}$ denote the pooled sample variance. Under $\Theta$, the unrestricted MLE of $\sigma^2$ is $\hat{\sigma}^2 = \frac{\nu S_p^2}{n_1 + n_2}$. Under $H_0 : \mu_1 = \mu_2 = \mu$, the restricted MLE of the common mean is the grand mean $\hat{\mu}_0 = \frac{n_1 \bar{X} + n_2 \bar{Y}}{n_1 + n_2}$, and the restricted residual sum of squares decomposes via algebra as:

Consequently, the two-sample GLRT statistic is:

When the population variances $\sigma_1^2 \neq \sigma_2^2$ cannot be assumed equal (the Behrens-Fisher Problem), Yuri Linnik (1968) proved that no exact non-randomized similar test with smooth critical boundary exists. In practice, one uses Welch's unpooled $t$-statistic $t_W = \frac{\bar{X} - \bar{Y}}{\sqrt{S_X^2/n_1 + S_Y^2/n_2}}$, where the linear combination of independent Chi-Square variances $V = \frac{S_X^2}{n_1} + \frac{S_Y^2}{n_2}$ is approximated by a scaled Chi-Square distribution via matching first and second moments (Welch-Satterthwaite effective degrees of freedom):

3. General Linear Hypothesis $F$-Test in the Gaussian Linear Model

More generally, let $\mathbf{Y} \sim \mathcal{N}_n(\mathbf{X}\boldsymbol{\beta}, \sigma^2 \mathbf{I}_n)$ be an $n$-dimensional response vector with full-rank $n \times p$ design matrix $\mathbf{X}$ ($p < n$), and consider testing $r$ linearly independent constraints $H_0 : \mathbf{C}\boldsymbol{\beta} = \mathbf{d}$, where $\mathbf{C} \in \mathbb{R}^{r \times p}$ has rank $r \le p$. Let $\text{SSE}_{\Omega} = \|\mathbf{Y} - \mathbf{X}\hat{\boldsymbol{\beta}}\|^2$ denote the unrestricted residual sum of squares and let $\text{SSE}_{\omega} = \|\mathbf{Y} - \mathbf{X}\hat{\boldsymbol{\beta}}_0\|^2$ denote the restricted residual sum of squares under $H_0$. Because the maximized Gaussian likelihood is inversely proportional to $(\text{SSE})^{n/2}$, the exact GLRT statistic is:

where $F = \frac{(\text{SSE}_{\omega} - \text{SSE}_{\Omega}) / r}{\text{SSE}_{\Omega} / (n - p)} \sim F_{r, n-p}$ under $H_0$. Thus every classical ANOVA and regression $F$-test is an exact finite-sample Generalized Likelihood Ratio Test!

Testing equality of variances across two or more normal populations is essential both for comparing dispersion directly and for validating the homoscedasticity assumption underlying pooled $t$-tests and Analysis of Variance (ANOVA).

1. Two-Sample GLRT for Equality of Variances ($H_0 : \sigma_1^2 = \sigma_2^2$)

Let $X_1, \dots, X_{n_1} \overset{\text{iid}}{\sim} \mathcal{N}(\mu_1, \sigma_1^2)$ and $Y_1, \dots, Y_{n_2} \overset{\text{iid}}{\sim} \mathcal{N}(\mu_2, \sigma_2^2)$ be independent normal samples. To test $H_0 : \sigma_1^2 = \sigma_2^2$ versus $H_1 : \sigma_1^2 \neq \sigma_2^2$ with $\mu_1, \mu_2$ unknown, the unrestricted MLEs are $\hat{\sigma}_1^2 = \frac{n_1 - 1}{n_1} S_X^2$ and $\hat{\sigma}_2^2 = \frac{n_2 - 1}{n_2} S_Y^2$, whereas under $H_0 : \sigma_1^2 = \sigma_2^2 = \sigma^2$, the restricted MLE is $\hat{\sigma}_0^2 = \frac{(n_1 - 1)S_X^2 + (n_2 - 1)S_Y^2}{n_1 + n_2}$. Forming the Generalized Likelihood Ratio statistic yields:

Notice that $\Lambda$ is a unimodal function of $F = S_X^2 / S_Y^2$ that vanishes as $F \to 0$ and as $F \to \infty$, achieving its global maximum at $F = \frac{n_1(n_2-1)}{n_2(n_1-1)}$. Thus the GLRT rejects $H_0$ in two tails: $F < c_1$ or $F > c_2$. When $n_1 \neq n_2$, the raw GLRT using $n_i$ divisors is slightly biased; replacing $n_1, n_2$ by the degrees of freedom $\nu_1 = n_1 - 1$ and $\nu_2 = n_2 - 1$ yields the modified likelihood ratio test, which coincides with the exact UMPU test!

2. Bartlett's Test for Homogeneity of $k$ Normal Variances (1937)

Suppose we have $k \ge 2$ independent normal groups $X_{ij} \sim \mathcal{N}(\mu_i, \sigma_i^2)$ ($i = 1, \dots, k$; $j = 1, \dots, n_i$), with group sample variances $S_i^2 = \frac{1}{\nu_i}\sum_{j=1}^{n_i}(X_{ij} - \bar{X}_{i\cdot})^2$ and degrees of freedom $\nu_i = n_i - 1$. Let $\nu = \sum_{i=1}^k \nu_i = N - k$ denote the total error degrees of freedom, and let $S_p^2$ denote the pooled sample variance:

To test the multi-sample homoscedasticity hypothesis:

Maurice Bartlett replaced the raw sample sizes $n_i$ in the $k$-sample GLRT by their unbiased degrees of freedom $\nu_i = n_i - 1$, defining the modified likelihood ratio statistic $\Lambda^* = \prod_{i=1}^k \left(\frac{S_i^2}{S_p^2}\right)^{\nu_i / 2}$. Taking $-2 \ln \Lambda^*$ gives the numerator:

Notice that $\frac{M}{\nu} = \ln(S_p^2) - \sum_{i=1}^k \frac{\nu_i}{\nu} \ln(S_i^2)$ is the logarithm of the weighted arithmetic mean of the sample variances minus the logarithm of their weighted geometric mean! By the Arithmetic-Geometric Mean Inequality, $M \ge 0$, with $M = 0$ if and only if $S_1^2 = S_2^2 = \dots = S_k^2$.

To improve the finite-sample convergence of $M$ to the asymptotic $\chi^2_{k-1}$ distribution, Bartlett derived the exact expectation $\mathbb{E}_0[M] \approx (k - 1) C$ via Stirling's expansion of the digamma function, where $C$ is Bartlett's correction factor:

At level $\alpha$, Bartlett's test rejects $H_0$ if $B > \chi^2_{k-1, 1-\alpha}$. In multivariate analysis of $p$-dimensional normal populations $\mathcal{N}_p(\boldsymbol{\mu}_i, \boldsymbol{\Sigma}_i)$, George Box (1949) generalized Bartlett's test to Box's $M$-test for equality of $k$ covariance matrices $H_0 : \boldsymbol{\Sigma}_1 = \dots = \boldsymbol{\Sigma}_k$, replacing $\ln(S_i^2)$ by the log-determinant $\ln |\mathbf{S}_i|$ with degrees of freedom $\frac{1}{2} p(p+1)(k-1)$.

3. Sensitivity to Kurtosis and the Levene / Brown-Forsythe Robust Alternatives

George Box (1953) demonstrated a celebrated cautionary result regarding normal-theory variance tests: while the ANOVA $F$-test for equality of means is remarkably robust to non-normality, the $F$-test and Bartlett's test for equality of variances are extrordinarily sensitive to population kurtosis $\gamma_2 = \frac{\mathbb{E}[(X-\mu)^4]}{\sigma^4} - 3$. Specifically, because the asymptotic variance of a sample variance $S_i^2$ is $\text{Var}(S_i^2) \approx \frac{\sigma^4}{n_i}(2 + \gamma_2)$ rather than the Gaussian value $\frac{2\sigma^4}{n_i}$, under any heavy-tailed distribution with $\gamma_2 > 0$, Bartlett's statistic $B$ converges under $H_0$ to $\left(1 + \frac{1}{2}\gamma_2\right) \chi^2_{k-1}$, causing severe inflation of Type I error rates! To remark on Box's famous aphorism: "To make the preliminary test on variances is rather like putting to sea in a rowing boat to find out whether conditions are sufficiently calm for an ocean liner to leave port!"

To achieve robustness against non-normality, Howard Levene (1960) and Morton Brown and Alan Forsythe (1974) replaced squared residuals by absolute deviations from the group median $\tilde{X}_i = \text{Median}(X_{i1}, \dots, X_{in_i})$, defining transformed observations $Z_{ij} = |X_{ij} - \tilde{X}_i|$ and applying a standard one-way ANOVA $F$-test to $Z_{ij}$, which maintains accurate Type I error control across heavy-tailed and skewed distributions.

When an exact finite-sample distribution for the Generalized Likelihood Ratio statistic $\Lambda(\mathbf{X})$ is unavailable, Samuel S. Wilks (1938) discovered a universal asymptotic theorem that elevates the GLRT to the centerpiece of large-sample parametric inference: under mild regularity conditions, $-2 \ln \Lambda(\mathbf{X})$ converges in distribution under $H_0$ to a central Chi-Square distribution whose degrees of freedom equal the number of independent parameter constraints imposed by $H_0$.

Regularity Conditions for Wilks' Theorem

For Wilks' asymptotic approximation to hold, the statistical family $\{f(\mathbf{x}; \boldsymbol{\theta}) : \boldsymbol{\theta} \in \Theta\}$ must satisfy the classical Cramér-Wald regularity conditions:

1. Identifiability and Common Support: Distinct parameters $\boldsymbol{\theta}_1 \neq \boldsymbol{\theta}_2$ define distinct probability distributions, and the support $\{\mathbf{x} : f(\mathbf{x}; \boldsymbol{\theta}) > 0\}$ does not depend on $\boldsymbol{\theta}$.

2. Interior Parameter Point: The true null parameter $\boldsymbol{\theta}_0 \in \Theta_0$ is an interior point of the $k$-dimensional open parameter space $\Theta \subseteq \mathbb{R}^k$, and the null subspace $\Theta_0$ is locally a smooth $(k - r)$-dimensional submanifold near $\boldsymbol{\theta}_0$.

3. Smoothness and Dominated Third Derivatives: The log-density $\ln f(\mathbf{x}; \boldsymbol{\theta})$ is three times continuously differentiable with respect to $\boldsymbol{\theta}$ in a neighborhood of $\boldsymbol{\theta}_0$, differentiation up to second order can be interchanged with integration, and the third partial derivatives are bounded in absolute value by an integrable function $M(\mathbf{x})$ with $\mathbb{E}_{\boldsymbol{\theta}_0}[M(\mathbf{X})] < \infty$.

4. Positive-Definite Fisher Information: The $k \times k$ single-observation Fisher Information matrix $\mathbf{I}_1(\boldsymbol{\theta}_0) = \mathbb{E}_{\boldsymbol{\theta}_0}[\nabla \ln f \nabla \ln f^T] = -\mathbb{E}_{\boldsymbol{\theta}_0}[\nabla^2 \ln f]$ is finite and strictly positive-definite.

Theorem 4.1: Wilks' Likelihood Ratio Theorem

Let $\mathbf{X}_1, \dots, \mathbf{X}_n$ be an i.i.d. sample from a regular parametric family $\{f(\mathbf{x}; \boldsymbol{\theta}) : \boldsymbol{\theta} \in \Theta\}$ satisfying the regularity conditions above. Suppose the null hypothesis $H_0 : \boldsymbol{\theta} \in \Theta_0$ restricts $\boldsymbol{\theta}$ to a smooth $(k - r)$-dimensional submanifold $\Theta_0 \subset \Theta$ defined by $r$ independent smooth constraints ($1 \le r \le k$):

Then under any true parameter $\boldsymbol{\theta}_0 \in \text{int}(\Theta_0)$, as sample size $n \to \infty$, the log-likelihood ratio test statistic satisfies:

Complete Second-Order Taylor Expansion Proof (Simple Null Case $r = k$):

To see why the factor of $-2$ and the $\chi^2_r$ distribution arise naturally, consider first the case where $H_0 : \boldsymbol{\theta} = \boldsymbol{\theta}_0$ specifies all $r = k$ parameters ($\Theta_0 = \{\boldsymbol{\theta}_0\}$). Let $\ell_n(\boldsymbol{\theta}) = \sum_{i=1}^n \ln f(\mathbf{X}_i; \boldsymbol{\theta})$ denote the log-likelihood, and let $\hat{\boldsymbol{\theta}}_n$ be the unrestricted MLE, which satisfies the first-order score equation $\nabla \ell_n(\hat{\boldsymbol{\theta}}_n) = \mathbf{0}$.

Expand $\ell_n(\boldsymbol{\theta}_0)$ in a second-order multivariable Taylor series around the unrestricted MLE $\hat{\boldsymbol{\theta}}_n$:

where $\tilde{\boldsymbol{\theta}}_n$ lies on the line segment connecting $\hat{\boldsymbol{\theta}}_n$ and $\boldsymbol{\theta}_0$. Because $\nabla \ell_n(\hat{\boldsymbol{\theta}}_n) = \mathbf{0}$, the linear term vanishes identically! Rearranging and multiplying by $-2$:

Now invoke two fundamental asymptotic theorems of Maximum Likelihood Estimation under $H_0 : \boldsymbol{\theta} = \boldsymbol{\theta}_0$:

1. Consistency and Law of Large Numbers for Observed Information: Because $\hat{\boldsymbol{\theta}}_n \xrightarrow{P} \boldsymbol{\theta}_0$, we also have $\tilde{\boldsymbol{\theta}}_n \xrightarrow{P} \boldsymbol{\theta}_0$. By the Weak Law of Large Numbers and continuity of the Hessian, the normalized negative Hessian converges in probability to the Fisher Information matrix:

2. Asymptotic Normality of the MLE: Under $H_0$, the scaled estimation error converges in distribution to a multivariate normal vector with covariance $\mathbf{I}_1(\boldsymbol{\theta}_0)^{-1}$:

By Slutsky's Theorem and the Continuous Mapping Theorem, the quadratic form converges in distribution to:

Since $\mathbf{Z} \sim \mathcal{N}_k(\mathbf{0}, \mathbf{I}_1(\boldsymbol{\theta}_0)^{-1})$, the standardized vector $\mathbf{Y} = \mathbf{I}_1(\boldsymbol{\theta}_0)^{1/2} \mathbf{Z} \sim \mathcal{N}_k(\mathbf{0}, \mathbf{I}_k)$ consists of $k$ independent standard normal random variables, so $\mathbf{Z}^T \mathbf{I}_1(\boldsymbol{\theta}_0) \mathbf{Z} = \mathbf{Y}^T \mathbf{Y} = \sum_{j=1}^k Y_j^2 \sim \chi^2_k$!

In the general composite case where $r < k$ parameters are fixed under $H_0$ and $k - r$ nuisance parameters are estimated via $\hat{\boldsymbol{\theta}}_0$, expanding both $\ell_n(\hat{\boldsymbol{\theta}})$ and $\ell_n(\hat{\boldsymbol{\theta}}_0)$ around the true $\boldsymbol{\theta}_0$ projects the $k$-dimensional score vector onto the orthogonal complement of the $(k - r)$-dimensional nuisance tangent space, leaving an idempotent quadratic form of rank $k - (k - r) = r$, yielding $\chi^2_r$.

Beside Wilks' Likelihood Ratio Test ($W_{\text{LR}}$), two companion large-sample testing procedures were introduced by Abraham Wald (1943) and C. Radhakrishna Rao (1948): the **Wald Test ($W_W$)** and the **Rao Score Test ($W_S$)** (known in econometrics as the **Lagrange Multiplier, LM, Test**). Together, $(W_{\text{LR}}, W_W, W_S)$ constitute the celebrated **Asymptotic Trinity** of parametric inference. ### Mathematical Formulation of the Three Trinity Statistics Consider testing the $r$-dimensional nonlinear or linear parameter constraint $H_0 : \mathbf{h}(\boldsymbol{\theta}) = \mathbf{0}$ (where $\mathbf{h} : \mathbb{R}^k \to \mathbb{R}^r$ has full-rank $r \times k$ Jacobian matrix $\mathbf{H}(\boldsymbol{\theta}) = \frac{\partial \mathbf{h}(\boldsymbol{\theta})}{\partial \boldsymbol{\theta}^T}$ of rank $r$). In the scalar parameter case ($k = r = 1$) testing $H_0 : \theta = \theta_0$ versus $H_1 : \theta \neq \theta_0$ with log-likelihood $\ell_n(\theta)$, score $S_n(\theta) = \ell_n'(\theta)$, and total Fisher information $I_n(\theta) = n I_1(\theta)$: 1. **The Likelihood Ratio Statistic ($W_{\text{LR}}$)**: Measures the **vertical drop** in the log-likelihood surface from its peak $\ell_n(\hat{\theta})$ to its null height $\ell_n(\theta_0)$: $$W_{\text{LR}} = 2 \Big[ \ell_n(\hat{\theta}) - \ell_n(\theta_0) \Big]$$ *Computational requirement*: Requires computing **both** the unrestricted MLE $\hat{\boldsymbol{\theta}}$ and the restricted MLE $\hat{\boldsymbol{\theta}}_0$. 2. **The Wald Statistic ($W_W$)**: Measures the **horizontal squared distance** between the unrestricted MLE $\hat{\theta}$ and the null value $\theta_0$, weighted by the curvature (Fisher information) at $\hat{\theta}$: $$W_W = n I_1(\hat{\theta}) \left(\hat{\theta} - \theta_0\right)^2 = \frac{\left(\hat{\theta} - \theta_0\right)^2}{\widehat{\text{Var}}(\hat{\theta})}$$ In the general $k$-parameter case testing $H_0 : \mathbf{h}(\boldsymbol{\theta}) = \mathbf{0}$: $$W_W = \mathbf{h}(\hat{\boldsymbol{\theta}})^T \Big[ \mathbf{H}(\hat{\boldsymbol{\theta}}) \, \mathbf{I}_n(\hat{\boldsymbol{\theta}})^{-1} \, \mathbf{H}(\hat{\boldsymbol{\theta}})^T \Big]^{-1} \mathbf{h}(\hat{\boldsymbol{\theta}})$$ *Computational requirement*: Requires computing **only the unrestricted MLE $\hat{\boldsymbol{\theta}}$**! 3. **The Rao Score / Lagrange Multiplier Statistic ($W_S$)**: Measures the **squared slope (gradient)** of the log-likelihood surface evaluated directly at the null hypothesis $\theta_0$ (where the slope should be near zero if $H_0$ is true), standardized by the null Fisher information $I_n(\theta_0)$: $$W_S = \frac{\left[S_n(\theta_0)\right]^2}{I_n(\theta_0)} = \frac{\left[\ell_n'(\theta_0)\right]^2}{n I_1(\theta_0)}$$ In the general $k$-parameter case with restricted MLE $\hat{\boldsymbol{\theta}}_0$: $$W_S = \mathbf{S}_n(\hat{\boldsymbol{\theta}}_0)^T \, \mathbf{I}_n(\hat{\boldsymbol{\theta}}_0)^{-1} \, \mathbf{S}_n(\hat{\boldsymbol{\theta}}_0)$$ *Computational requirement*: Requires computing **only the restricted null MLE $\hat{\boldsymbol{\theta}}_0$**! Notice why $W_S$ is called the **Lagrange Multiplier Test**: when maximizing $\ell_n(\boldsymbol{\theta}) - \boldsymbol{\lambda}^T \mathbf{h}(\boldsymbol{\theta})$ under $H_0$, the first-order optimality condition at $\hat{\boldsymbol{\theta}}_0$ is $\mathbf{S}_n(\hat{\boldsymbol{\theta}}_0) = \mathbf{H}(\hat{\boldsymbol{\theta}}_0)^T \hat{\boldsymbol{\lambda}}$, so $W_S$ is literally a quadratic form in the shadow price (Lagrange multiplier) vector $\hat{\boldsymbol{\lambda}}$ measuring how strongly the data pull against the null constraint! ### Asymptotic Equivalence and Local Non-Central Chi-Square Power **Theorem 4.2 (Asymptotic Equivalence of the Trinity)**: Under standard regularity conditions, under both the null hypothesis $H_0 : \boldsymbol{\theta} = \boldsymbol{\theta}_0$ and under any sequence of contiguous local Pitman alternatives $H_{1,n} : \boldsymbol{\theta}_n = \boldsymbol{\theta}_0 + \mathbf{h}/\sqrt{n}$: $$W_{\text{LR}} - W_W \xrightarrow{P} 0, \qquad W_{\text{LR}} - W_S \xrightarrow{P} 0, \qquad W_W - W_S \xrightarrow{P} 0$$ Consequently: 1. Under $H_0$, all three statistics converge in distribution to the same central Chi-Square distribution: $W_{\text{LR}}, W_W, W_S \xrightarrow{d} \chi^2_r$. 2. Under the local alternative $\boldsymbol{\theta}_n = \boldsymbol{\theta}_0 + \mathbf{h}/\sqrt{n}$, all three statistics converge in distribution to the **non-central Chi-Square distribution** $\chi^2_r(\delta)$ with $r$ degrees of freedom and non-centrality parameter $\delta = \mathbf{h}^T \mathbf{I}_1(\boldsymbol{\theta}_0) \mathbf{h}$, achieving identical first-order Pitman asymptotic power! 3. When the parameter of interest is partitioned from nuisance parameters as $\boldsymbol{\theta} = (\boldsymbol{\vartheta}, \boldsymbol{\eta})$, the block-inverse Schur complement of the Fisher Information matrix replaces $\mathbf{I}_1(\boldsymbol{\theta}_0)$ by the **efficient Fisher Information matrix** $\mathbf{I}_{\vartheta\vartheta \cdot \eta} = \mathbf{I}_{\vartheta\vartheta} - \mathbf{I}_{\vartheta\eta} \mathbf{I}_{\eta\eta}^{-1} \mathbf{I}_{\eta\vartheta}$, quantifying the exact asymptotic power reduction induced by estimating nuisance parameters. ### Confidence Set Inversion of the Asymptotic Trinity Because every point null test $H_0 : \boldsymbol{\theta} = \boldsymbol{\theta}_0$ can be inverted into a $100(1-\alpha)\%$ confidence region by collecting all parameter values $\boldsymbol{\theta}_0$ that are not rejected at level $\alpha$, the Asymptotic Trinity generates three distinct large-sample confidence regions: - **Wald Confidence Ellipsoid**: Collecting all $\boldsymbol{\theta}_0$ such that $(\hat{\boldsymbol{\theta}} - \boldsymbol{\theta}_0)^T \mathbf{I}_n(\hat{\boldsymbol{\theta}})(\hat{\boldsymbol{\theta}} - \boldsymbol{\theta}_0) \le \chi^2_{p, 1-\alpha}$ produces an exact **symmetric geometric ellipsoid** centered at $\hat{\boldsymbol{\theta}}$. In one dimension, this yields the familiar symmetric interval $\hat{\theta} \pm z_{1-\alpha/2} / \sqrt{I_n(\hat{\theta})}$. However, symmetric Wald intervals can spill outside the legal parameter space (for instance, producing negative lower bounds for a Binomial proportion $p$ near zero or a variance $\sigma^2$) and violate transformation invariance. - **Profile Likelihood Ratio Confidence Region**: Collecting all $\boldsymbol{\theta}_0$ such that $-2\ln\Lambda(\boldsymbol{\theta}_0) = 2[\ell_n(\hat{\boldsymbol{\theta}}) - \ell_n(\boldsymbol{\theta}_0)] \le \chi^2_{p, 1-\alpha}$ produces a **likelihood contour region** whose shape automatically adapts to the skewness and natural boundaries of the log-likelihood surface! Moreover, because the likelihood ratio is strictly invariant under smooth reparameterizations $\boldsymbol{\psi} = \mathbf{g}(\boldsymbol{\theta})$, the profile likelihood confidence region transforms equivariantly without recalibration. - **Score (Wilson-Type) Confidence Region**: Collecting all $\boldsymbol{\theta}_0$ such that $\mathbf{S}_n(\boldsymbol{\theta}_0)^T \mathbf{I}_n(\boldsymbol{\theta}_0)^{-1} \mathbf{S}_n(\boldsymbol{\theta}_0) \le \chi^2_{p, 1-\alpha}$ accounts for how the theoretical variance changes as a function of the hypothesized parameter $\boldsymbol{\theta}_0$. In the Binomial case, solving $W_S(p_0) \le z_{1-\alpha/2}^2$ as a quadratic equation in $p_0$ yields **Edwin Bidwell Wilson's (1927) score interval**, which remains strictly inside $(0, 1)$ even when zero successes ($x = 0$) are observed and maintains dramatically superior finite-sample coverage compared to the Wald interval.

A critical step in the proof of Wilks' Theorem (Section 4.4) is the assumption that the true null parameter $\boldsymbol{\theta}_0$ is an interior point of the full parameter space $\Theta$, so that the unrestricted MLE $\hat{\boldsymbol{\theta}}_n$ can fluctuate freely in all directions around $\boldsymbol{\theta}_0$ and satisfy the unconstrained zero-gradient equation $\nabla \ell_n(\hat{\boldsymbol{\theta}}_n) = \mathbf{0}$. When $\boldsymbol{\theta}_0$ lies on the topological boundary of the parameter space $\Theta$, or when parameters become unidentifiable under $H_0$, Wilks' $\chi^2_r$ distribution breaks down completely. Understanding these non-regular regimes is essential for modern applied statistics, econometrics, biostatistics, and machine learning.

Theorem 4.3: Chernoff's Boundary Theorem and the $\frac{1}{2}\chi^2_0 + \frac{1}{2}\chi^2_1$ Mixture Distribution

Suppose we test a scalar parameter that is physically or mathematically constrained to be non-negative ($\Theta = [0, \infty)$), such as testing whether a random-effects variance component $\sigma_u^2$ in a linear or generalized linear mixed model is zero, testing for overdispersion in a Negative Binomial vs. Poisson regression model, or testing a non-negative signal amplitude above background noise in particle physics:

Let $\hat{\theta}_{\text{uncon}}$ denote the unconstrained maximizer of the smooth extension of $\ell_n(\theta)$ over $\mathbb{R}$, which under $H_0 : \theta = 0$ satisfies $\sqrt{n I_1(0)} \, \hat{\theta}_{\text{uncon}} \xrightarrow{d} Z \sim \mathcal{N}(0, 1)$. Because the true parameter space is restricted to the non-negative half-line $\Theta = [0, \infty)$, the constrained maximum likelihood estimator $\hat{\theta}_n = \arg\max_{\theta \ge 0} \ell_n(\theta)$ projects negative unconstrained estimates directly onto the boundary point $\theta = 0$:

Consequently, under $H_0 : \theta = 0$, two distinct sample-space regimes each occur with asymptotic probability $1/2$:

1. When $Z \le 0$ (asymptotic probability $1/2$): The unconstrained maximizer is negative ($\hat{\theta}_{\text{uncon}} \le 0$), so the constrained MLE on $[0, \infty)$ is truncated at the boundary $\hat{\theta}_n = 0 = \theta_0$! Because the unrestricted MLE on $\Theta$ coincides with the restricted null parameter ($\hat{\theta}_n = \theta_0 = 0$), the likelihood ratio is identically $\Lambda(\mathbf{X}) = L(0) / L(0) = 1$, and the log-likelihood ratio test statistic is identically zero:

2. When $Z > 0$ (asymptotic probability $1/2$): The unconstrained maximizer falls strictly in the interior ($\hat{\theta}_n = \hat{\theta}_{\text{uncon}} > 0$), satisfying the zero-score equation $\ell_n'(\hat{\theta}_n) = 0$, and the second-order Taylor expansion yields $W_{\text{LR}} \approx Z^2 \sim \chi^2_1$ (since the square of a standard normal variable conditioned on being positive, $Z^2 \mid (Z > 0)$, has the exact same $\chi^2_1$ distribution as $Z^2$ by symmetry of the Gaussian density!).

Combining these two equiprobable cases establishes Herman Chernoff's Boundary Theorem (1954) (later extended to general convex cones by Self and Liang, 1987): under the boundary null hypothesis $H_0 : \theta = 0$ with $\Theta = [0, \infty)$, the asymptotic distribution of $-2 \ln \Lambda(\mathbf{X})$ is a 50:50 mixture of a degenerate point mass at zero ($\chi^2_0 \equiv \delta_0$) and a standard $\chi^2_1$ distribution (often denoted $\bar{\chi}^2_{01}$ or $\frac{1}{2}\chi^2_0 + \frac{1}{2}\chi^2_1$):

For any positive threshold $c > 0$, the asymptotic Type I error probability is exactly one-half of the standard $\chi^2_1$ upper-tail probability:

Therefore, to achieve significance level $\alpha = 0.05$, one must set $\frac{1}{2}\mathbb{P}(\chi^2_1 > c_\alpha) = 0.05 \implies \mathbb{P}(\chi^2_1 > c_\alpha) = 2\alpha = 0.10$, yielding the critical value $c_{0.05} = \chi^2_{1, 0.90} = z_{0.95}^2 = (1.6449)^2 \approx 2.7055$ rather than the standard interior cutoff $\chi^2_{1, 0.95} \approx 3.8415$! Naively comparing $-2\ln\Lambda$ against $3.8415$ when testing a variance component $\sigma_u^2 = 0$ cuts the actual Type I error rate in half to $0.025$ and severely depresses statistical power.

Singular Models, Unidentifiable Parameters Under $H_0$, and Finite Mixture Testing

An even more dramatic failure of Wilks' Theorem occurs in finite mixture models and latent variable models, such as testing whether a population is homogeneous ($H_0 : \text{single normal component}$) or a two-component Gaussian mixture:

Notice that the null hypothesis of a single normal distribution can be represented in two completely different ways within the parameter space:

1. By setting the mixing weight on the boundary: $\pi = 0$ (in which case the second component mean $\mu_2$ completely disappears from the likelihood and becomes unidentifiable under $H_0$!).
2. By setting the component means equal: $\mu_1 = \mu_2$ (in which case the mixing proportion $\pi \in [0, 1]$ disappears from the likelihood and becomes unidentifiable under $H_0$, while the derivative Score vector with respect to $\mu_1$ and $\mu_2$ becomes linearly dependent, rendering the Fisher Information matrix singular!).

Because the parameter $\mu_2$ is unidentifiable when $\pi = 0$, the unconstrained likelihood maximization searches over all possible locations $\mu_2 \in \mathbb{R}$ for the best local noise cluster, causing $-2\ln\Lambda(\mathbf{X})$ to behave asymptotically like the supremum of a squared Gaussian process $\sup_{\mu_2} Z(\mu_2)^2$ rather than a finite-degree-of-freedom Chi-Square variable! In such singular models (analyzed systematically by Dacunha-Castelle and Gassiat, 1999, and Sumio Watanabe's Singular Learning Theory), practitioners must calibrate the critical value of $-2\ln\Lambda(\mathbf{X})$ using parametric bootstrap simulation under $H_0$ rather than consulting standard Chi-Square tables.

Higher-Order Bartlett Correction, Signed Root LRT ($r^*$), and the Hauck-Donner Phenomenon

Three additional finite-sample phenomena distinguish the Asymptotic Trinity in regular interior models:

• General Bartlett Adjustment (Lawley, 1956): Under interior regularity conditions, Donald Lawley proved that all cumulants of $W_{\text{LR}} = -2 \ln \Lambda(\mathbf{X})$ are scaled by the exact same $O(n^{-1})$ factor! Specifically, the finite-sample expectation of $W_{\text{LR}}$ admits the expansion $\mathbb{E}_0[W_{\text{LR}}] = r \left(1 + \frac{b}{n} + O(n^{-2})\right)$. Dividing $W_{\text{LR}}$ by its rescaled expectation factor $1 + b/n$ produces the Bartlett-corrected statistic $W_{\text{LR}}^* = \frac{W_{\text{LR}}}{1 + b/n}$, which simultaneously matches all cumulants of the $\chi^2_r$ distribution to order $O(n^{-2})$ and reduces the distributional tail error from $O(n^{-1})$ to $O(n^{-2})$ across all quantiles!
• Signed Root Likelihood Ratio Statistic and Barndorff-Nielsen's $r^*$: For a scalar parameter $\theta \in \mathbb{R}$ when testing a one-sided hypothesis $H_0 : \theta \le \theta_0$ versus $H_1 : \theta > \theta_0$, the unsigned statistic $W_{\text{LR}} \ge 0$ loses directional sign information. One therefore defines the signed root likelihood ratio statistic $r(\theta_0) = \text{sign}(\hat{\theta} - \theta_0)\sqrt{W_{\text{LR}}} \xrightarrow{d} \mathcal{N}(0, 1)$, and Ole Barndorff-Nielsen's higher-order magic formula $r^* = r + \frac{1}{r}\ln(q/r)$ achieves third-order normal accuracy $O(n^{-3/2})$ even for very small sample sizes!
• The Hauck-Donner Effect in Wald Tests (1977): Unlike $W_{\text{LR}}$ and $W_S$, the Wald statistic $W_W = (\hat{\theta} - \theta_0)^2 I_n(\hat{\theta})$ evaluates the Fisher information at the unrestricted MLE $\hat{\theta}$. Walter Hauck and Allan Donner (1977) discovered that in binomial logistic regression, contingency tables with near-complete separation, and Poisson log-linear models, when the true effect is very large and $|\hat{\theta} - \theta_0|$ becomes extreme, the estimated information $I_n(\hat{\theta})$ decays toward zero exponentially faster than $(\hat{\theta} - \theta_0)^2$ grows. Consequently, the Wald statistic $W_W$ decreases back toward zero as the evidence against $H_0$ becomes stronger, failing to reject a blatantly false null hypothesis! Because the Likelihood Ratio Test $W_{\text{LR}}$ and the Score Test $W_S$ are completely immune to the Hauck-Donner pathology and invariant to nonlinear parameter transformations, modern statistical software uses profile likelihood ratio tests and score tests whenever boundary proximity or large effect sizes are possible.