Deriving Box’s M Test For Equal Covariance Matrices

Why Compare Covariance Matrices Across Groups?

In a multivariate study, groups can have similar mean vectors but different patterns of variability and association. For example, treatment groups might have the same average values on several outcomes while showing different outcome correlations. Box’s M tests whether the groups share a common covariance matrix. Its construction uses the within-group scatter matrices and compares their determinants with the determinant of a pooled covariance estimate.

The test is derived under independent multivariate normal sampling. The means are allowed to differ freely across groups, so the hypothesis concerns covariance matrices alone. The derivation below uses the likelihood of the groups’ residual scatter matrices, which have Wishart distributions after the group means are estimated.

Definition: Suppose there are $g\geq 2$ independent groups, each with $n_i$ observations in $p$ dimensions. In group $i$, assume $X_{ij}\sim N_p(\mu_i,\Sigma_i)$, independently for $j=1,\ldots,n_i$. Let $\bar X_i$ be the group sample mean and define the residual scatter and unbiased sample covariance by $$ W_i=\sum_{j=1}^{n_i}(X_{ij}-\bar X_i)(X_{ij}-\bar X_i)^\mathsf{T}, \qquad S_i=\frac{W_i}{n_i-1}. $$ Box’s M tests $H_0:\Sigma_1=\cdots=\Sigma_g$ against the alternative that the group covariance matrices are not all equal. Assume each $\Sigma_i$ is positive definite and $n_i-1\geq p$, so each $S_i$ is positive definite with probability one.

Write $\nu_i=n_i-1$ for the residual degrees of freedom and $\nu=\sum_{i=1}^g\nu_i$. The pooled covariance estimate under the equality constraint is a degrees-of-freedom-weighted average, not an unweighted average of the group covariances.

Definition: The pooled covariance estimate is $$ S_p=\frac{\sum_{i=1}^g W_i}{\nu} =\frac{\sum_{i=1}^g \nu_i S_i}{\sum_{i=1}^g \nu_i}. $$ The conventional Box’s M statistic is $$ M=\nu\log|S_p|-\sum_{i=1}^g \nu_i\log|S_i|. $$ Here $|A|$ denotes the determinant of $A$, and $\log$ is the natural logarithm.

Deriving the Pooled Estimate and the Likelihood Ratio

For multivariate normal observations, each residual scatter matrix has a Wishart distribution: $W_i\sim W_p(\Sigma_i,\nu_i)$. These matrices are independent across groups. As a function of $\Sigma_i$, the Wishart likelihood is proportional to $|\Sigma_i|^{-\nu_i/2}\exp\{-\tfrac12\operatorname{tr}(\Sigma_i^{-1}W_i)\}$. Terms that do not involve the covariance parameters can be omitted when maximizing.

Result: Under the unrestricted model, the maximizer of the residual-scatter likelihood in group $i$ is $\hat\Sigma_i=S_i$. Under $H_0$, the common-covariance maximizer is $\hat\Sigma=S_p$.

To verify these estimates, consider first one group with scatter $W$ and degrees of freedom $\nu_0$. Its covariance-dependent log-likelihood is $$ \ell(\Sigma)=-\frac{\nu_0}{2}\log|\Sigma|-\frac12\operatorname{tr}(\Sigma^{-1}W)+C, $$ where $C$ does not depend on $\Sigma$. Differentiating with respect to the positive-definite matrix $\Sigma$, the likelihood equation is $\nu_0\Sigma^{-1}=\Sigma^{-1}W\Sigma^{-1}$. Multiplying on the left and right by $\Sigma$ gives $W=\nu_0\Sigma$. Thus the maximizer is $W/\nu_0$, which gives $S_i=W_i/\nu_i$ in each unrestricted group. Under $H_0$, the same calculation uses the summed scatter $\sum_i W_i$ and total degrees of freedom $\nu$, giving $S_p$. The log-likelihood is maximized there: equivalently, in the parameterization $A=\Sigma^{-1}$, it is a strictly concave function of $A$ on the positive-definite cone. $\square$

The likelihood-ratio statistic compares the constrained maximum likelihood with the unrestricted maximum. At the unrestricted estimates, the trace contribution from group $i$ is $\operatorname{tr}(S_i^{-1}W_i)=\nu_i p$. At the constrained estimate, it is $\operatorname{tr}(S_p^{-1}\sum_i W_i)=\nu p$. These contributions are equal because $\nu=\sum_i\nu_i$, so they cancel in the likelihood ratio.

Derivation: Let $L_0$ and $L_1$ denote the maximized residual-scatter likelihoods under $H_0$ and the unrestricted alternative, respectively. Then $$ -2\log\left(\frac{L_0}{L_1}\right) =\nu\log|S_p|-\sum_{i=1}^g\nu_i\log|S_i| =M. $$

This identity explains the determinant comparison in Box’s M. A covariance matrix with a large determinant has a large generalized variance: loosely, it describes greater overall multivariate spread. The statistic is zero when all sample covariance matrices are identical, and it becomes larger as their determinants and orientations depart from what a common covariance matrix would imply. It is nonnegative because the constrained model is nested within the unrestricted one, so its maximized likelihood cannot be larger.

Degrees of Freedom and Asymptotic Calibration

A symmetric $p\times p$ covariance matrix has $p$ diagonal entries and $p(p-1)/2$ distinct off-diagonal entries. Its number of free parameters is therefore $p(p+1)/2$. The unrestricted model has $g$ covariance matrices, while the null model has one. The difference in parameter counts gives the degrees of freedom for the likelihood-ratio test.

Theorem: For fixed $p$ and $g$, under $H_0$, if the independent multivariate normal groups satisfy the stated positive-definiteness conditions and every $n_i$ tends to infinity, the residual-scatter likelihood-ratio statistic has the first-order asymptotic distribution $$ M\xrightarrow{d}\chi^2_q, \qquad q=(g-1)\frac{p(p+1)}{2}. $$

The parameter-count difference is $g\,p(p+1)/2-p(p+1)/2=(g-1)p(p+1)/2$. To obtain the asymptotic calibration, apply Wilks’ theorem to the Wishart likelihood: under the null the true common covariance is an interior point of the positive-definite parameter space, and the null imposes $q$ independent smooth restrictions on the unrestricted covariance parameters. With fixed $p,g$ and group sizes increasing, the regular likelihood-ratio limit is chi-square with $q$ degrees of freedom. $\square$

This chi-square statement is a large-sample approximation, not an assertion that $M$ has an exact chi-square distribution for finite group sizes. Box’s M is particularly sensitive to departures from multivariate normality, and its finite-sample calibration uses a correction to the raw statistic. The correction is a separate refinement; the derivation here identifies the uncorrected likelihood-ratio form and its first-order reference distribution.

Worked Calculations

Worked Example: Two Groups With Opposing Variability

There are two groups with $n_1=n_2=5$ observations in $p=2$ dimensions. Suppose their sample covariance matrices are $S_1=\begin{pmatrix}1&0\\0&4\end{pmatrix}$ and $S_2=\begin{pmatrix}4&0\\0&1\end{pmatrix}$. Compute the pooled estimate, Box’s M, and its degrees of freedom.

Each group has $\nu_i=5-1=4$, so $\nu=8$. The pooled estimate is $S_p=(4S_1+4S_2)/8=(S_1+S_2)/2$, giving $S_p=\begin{pmatrix}2.5&0\\0&2.5\end{pmatrix}$. Its determinant is $2.5(2.5)=6.25$, while $|S_1|=|S_2|=4$.

Consequently, $M=8\log(6.25)-4\log(4)-4\log(4)=8\log(6.25/4)=8\log(1.5625)\approx3.570$. The degrees of freedom are $q=(2-1)\,2(2+1)/2=3$. The two groups have the same determinant, but their covariance matrices differ in which coordinate has greater variance; the pooled determinant comparison captures that discrepancy.

Worked Example: Unequal Group Sizes and Weighted Pooling

Consider two groups in $p=2$ dimensions with $n_1=3$, $n_2=5$, and $S_1=\begin{pmatrix}2&0\\0&1\end{pmatrix}$, $S_2=\begin{pmatrix}1&0\\0&3\end{pmatrix}$. Calculate the pooled estimate and $M$.

Here $\nu_1=2$, $\nu_2=4$, and $\nu=6$. Weighting by residual degrees of freedom gives $S_p=(2S_1+4S_2)/6$. Its first diagonal entry is $(2\cdot2+4\cdot1)/6=8/6=4/3$; its second is $(2\cdot1+4\cdot3)/6=14/6=7/3$. Thus $|S_p|=(4/3)(7/3)=28/9$.

The group determinants are $|S_1|=2$ and $|S_2|=3$, so $M=6\log(28/9)-2\log2-4\log3\approx6(1.13498)-2(0.69315)-4(1.09861)\approx1.029$. There are $q=3$ degrees of freedom. The larger group receives twice the weight of the smaller group because it contributes twice as many residual degrees of freedom.

Worked Example: The Univariate Case

For three independent univariate normal groups, suppose each group has $n_i=3$ and sample variances $S_1=1$, $S_2=2$, and $S_3=4$. Compute the pooled variance, $M$, and degrees of freedom.

Each $\nu_i=2$, so $\nu=6$. The pooled variance is $S_p=(2\cdot1+2\cdot2+2\cdot4)/6=14/6=7/3$. For $p=1$, determinants are just the variances themselves. Therefore $M=6\log(7/3)-2\log1-2\log2-2\log4$. Since $\log1=0$ and $\log2+\log4=\log8$, $M=6(0.84730)-2(2.07944)\approx5.0838-4.1589=0.925$.

The degrees of freedom are $q=(3-1)\,1(1+1)/2=2$. This is the covariance-equality likelihood-ratio construction in one dimension: the group covariance matrices reduce to ordinary variances.

Interpretation and Common Pitfalls

Box’s M asks whether the observed group covariance matrices are sufficiently similar to support a common covariance model. This question matters whenever an analysis relies on common within-group covariance, such as some multivariate comparisons of group means. A rejection indicates evidence against covariance equality; it does not identify which group or which covariance feature is responsible. The test also does not establish that the groups are multivariate normal, an assumption used in its derivation.

Common pitfall: It is incorrect to pool the $S_i$ by taking their unweighted arithmetic mean when group sizes differ. For the second worked example, that would give diagonal entries $1.5$ and $2$, rather than the correct $4/3$ and $7/3$. The correct weights are $\nu_i=n_i-1$, because each $S_i$ summarizes that many residual degrees of freedom. Another potential confusion is to replace these weights by $n_i$: the conventional Box’s M formula uses the residual-scatter Wishart likelihood after fitting the group means. It is not the same formula as the raw-data Gaussian profile likelihood statistic formed with maximum-likelihood covariance estimates $W_i/n_i$.

Check Your Understanding

Use the definitions and derivation above to answer the following.

Why does the pooled covariance estimate weight $S_i$ by $n_i-1$, rather than by equal weights?
For $g=4$ groups and $p=3$ outcomes, what are the degrees of freedom for the uncorrected Box’s M calibration?
If every group covariance estimate is identical, what value does $M$ take? Verify this from its determinant formula.
Why do the trace terms cancel when forming the residual-scatter likelihood ratio?
What assumption about the groups’ distributions underlies the Wishart likelihood derivation?