The D1 Multivariate Wald Test

Why Use a Joint Test?

A model can contain several coefficients whose values matter jointly. For example, a categorical predictor represented by several indicator variables has no single coefficient that captures its entire association with an outcome. Testing each coefficient separately does not answer the joint question “are all these coefficients zero?” A D1 multivariate Wald test uses the pooled estimates and their covariance to test such a joint hypothesis after multiple imputation.

The method requires more than pooling each coefficient’s standard error separately. A joint test depends on the covariances among the estimates, as well as the uncertainty added by variation between imputations. The D1 procedure incorporates both sources of uncertainty and refers its statistic to an approximate F distribution.

Definition: Suppose $m$ completed data sets give estimates $\hat{\boldsymbol{\beta}}_j$ and estimated within-imputation covariance matrices $\mathbf{U}_j$, for $j=1,\ldots,m$. A joint linear null hypothesis is $H_0:\mathbf{R}\boldsymbol{\beta}=\mathbf{r}_0$, where $\mathbf{R}$ has $q$ linearly independent rows and $\mathbf{r}_0$ is a specified vector. The D1 test pools the estimate and covariance of these $q$ restrictions and calculates an approximate F test.

Pooling the Tested Restrictions

First average the coefficient estimates and their estimated within-imputation covariance matrices. The between-imputation covariance describes how much the coefficient estimates vary across the $m$ completed data sets. The total covariance combines the average within-imputation covariance with an adjusted between-imputation covariance.

\bar{\boldsymbol{\beta}}=\frac{1}{m}\sum_{j=1}^{m}\hat{\boldsymbol{\beta}}_j,\qquad \bar{\mathbf{U}}=\frac{1}{m}\sum_{j=1}^{m}\mathbf{U}_j

\mathbf{B}=\frac{1}{m-1}\sum_{j=1}^{m} (\hat{\boldsymbol{\beta}}_j-\bar{\boldsymbol{\beta}}) (\hat{\boldsymbol{\beta}}_j-\bar{\boldsymbol{\beta}})^\mathsf{T}, \qquad \mathbf{T}=\bar{\mathbf{U}}+\left(1+\frac{1}{m}\right)\mathbf{B}

For the null hypothesis, the pooled estimate minus the hypothesized value is $\bar{\mathbf{Q}}=\mathbf{R}\bar{\boldsymbol{\beta}}-\mathbf{r}_0$. The within- and between-imputation covariance matrices in the space of the restrictions are $\bar{\mathbf{U}}_R=\mathbf{R}\bar{\mathbf{U}}\mathbf{R}^\mathsf{T}$ and $\mathbf{B}_R=\mathbf{R}\mathbf{B}\mathbf{R}^\mathsf{T}$. Thus, the total covariance of the tested restrictions is $\mathbf{T}_R=\bar{\mathbf{U}}_R+(1+1/m)\mathbf{B}_R$. Equivalently, $\mathbf{T}_R=\mathbf{R}\mathbf{T}\mathbf{R}^\mathsf{T}$.

Key result: For $q$ linearly independent restrictions, the D1 statistic and its relative increase in variance are $$ F_{\mathrm{D1}}=\frac{1}{q}\bar{\mathbf{Q}}^\mathsf{T}\mathbf{T}_R^{-1}\bar{\mathbf{Q}}, \qquad r_1=\frac{1+1/m}{q}\operatorname{tr}\left(\mathbf{B}_R\bar{\mathbf{U}}_R^{-1}\right). $$ Under the D1 approximation, compare $F_{\mathrm{D1}}$ with an $F$ distribution having $q$ numerator degrees of freedom and $\nu=4+(q+1)/r_1^2$ denominator degrees of freedom. If $r_1=0$, take the denominator degrees of freedom to be infinite.

The trace in $r_1$ summarizes the between-imputation contribution relative to within-imputation uncertainty across the tested directions. The factor $1+1/m$ applies the finite-number correction to the between-imputation component. The p-value is the upper-tail probability $P(F_{q,\nu}\geq F_{\mathrm{D1}})$, evaluated under the approximate reference distribution. It is not the probability that the null hypothesis is true.

Why the Statistic Has This Form

The quadratic form measures the pooled displacement from the null relative to its pooled covariance. To see the role of the covariance, suppose the $q$-vector $\bar{\mathbf{Q}}$ is transformed to coordinates in which $\mathbf{T}_R$ is diagonal. Each squared coordinate is then divided by its corresponding variance; their sum is the squared multivariate distance from the null. The trace-based $r_1$ separately summarizes how much the between-imputation covariance contributes relative to within-imputation covariance. The D1 reference-distribution approximation uses these quantities to account for the uncertainty in a test based on multiply imputed data.

This construction also explains why the full covariance matrix matters. If two tested estimates move together, their joint uncertainty is not the same as the uncertainty obtained by pretending they are independent. Using only the diagonal variances can therefore change the test statistic and the resulting conclusion.

Worked Example: A Two-Coefficient Joint Test

Suppose $m=5$ hypothetical imputations produce pooled coefficient estimates $\bar{\boldsymbol{\beta}}=(0.30,-0.20)^\mathsf{T}$. The hypothesis is that both coefficients are zero, so $\mathbf{R}$ is the $2$-by-$2$ identity matrix, $\mathbf{r}_0=(0,0)^\mathsf{T}$, and $q=2$. Suppose the pooled within-imputation and between-imputation covariance matrices are

\bar{\mathbf{U}}_R= \begin{pmatrix}0.040&0\\0&0.090\end{pmatrix}, \qquad \mathbf{B}_R= \begin{pmatrix}0.010&0\\0&0.010\end{pmatrix}.

These are hypothetical inputs; the within-imputation variances are positive, and the between-imputation matrix is positive semidefinite. The two restrictions are linearly independent, and the total covariance is nonsingular. The finite-number multiplier is $1+1/5=1.2$, giving

\mathbf{T}_R= \begin{pmatrix}0.040&0\\0&0.090\end{pmatrix} +1.2\begin{pmatrix}0.010&0\\0&0.010\end{pmatrix} = \begin{pmatrix}0.052&0\\0&0.102\end{pmatrix}.

The pooled displacement is $\bar{\mathbf{Q}}=(0.30,-0.20)^\mathsf{T}$. Since $\mathbf{T}_R$ is diagonal, the quadratic form is

\bar{\mathbf{Q}}^\mathsf{T}\mathbf{T}_R^{-1}\bar{\mathbf{Q}} =\frac{0.30^2}{0.052}+\frac{(-0.20)^2}{0.102} =1.7308+0.3922=2.1229.

Therefore $F_{\mathrm{D1}}=2.1229/2=1.0615$. The relative increase is

r_1=\frac{1.2}{2}\left(\frac{0.010}{0.040}+\frac{0.010}{0.090}\right) =0.2167, \qquad \nu=4+\frac{3}{0.2167^2}=67.9.

The upper-tail probability for $F_{2,67.9}=1.0615$ is approximately $0.352$. At $\alpha=0.05$, the joint test does not provide evidence against the null that both coefficients are zero. This conclusion concerns the joint hypothesis, not whether each coefficient separately has a particular scalar-test p-value.

Interpreting the Result and Checking the Inputs

The degrees of freedom are approximate and depend on the relative contribution of between-imputation variation. A larger $r_1$ generally yields a smaller denominator degrees of freedom, reflecting greater uncertainty from missing information. Do not substitute a scalar Rubin degrees-of-freedom formula for the D1 multivariate formula: the D1 test is calibrated for a joint set of restrictions.

Before interpreting a result, verify that every completed-data fit estimates the same coefficient vector and that the same restrictions are being tested in every imputation. Check that the within-imputation covariance matrices are valid covariance estimates, that $\bar{\mathbf{U}}_R$ can be inverted for calculating $r_1$, and that $\mathbf{T}_R$ can be inverted for the Wald statistic. The restriction matrix must have rank $q$; redundant rows do not define $q$ independent restrictions. Serious departures from the approximate normal, large-sample setting can also make the reference F distribution less reliable.

Worked Example: A Joint Test With Correlated Estimates

Consider a second hypothetical test of two coefficients, both hypothesized to be zero, with $m=4$, pooled displacement $\bar{\mathbf{Q}}=(0.20,-0.10)^\mathsf{T}$, and

\bar{\mathbf{U}}_R= \begin{pmatrix}0.040&0.010\\0.010&0.090\end{pmatrix}, \qquad \mathbf{B}_R= \begin{pmatrix}0.010&0.004\\0.004&0.016\end{pmatrix}.

The within-imputation matrix has determinant $0.040(0.090)-0.010^2=0.0035>0$, and the between-imputation matrix is positive definite, so both are valid nonsingular covariance inputs. The total covariance is

\mathbf{T}_R=\bar{\mathbf{U}}_R+1.25\mathbf{B}_R =\begin{pmatrix}0.0525&0.015\\0.015&0.110\end{pmatrix}.

Its determinant is $0.0525(0.110)-0.015^2=0.00555$, so its inverse is

\mathbf{T}_R^{-1}=\frac{1}{0.00555} \begin{pmatrix}0.110&-0.015\\-0.015&0.0525\end{pmatrix}.

Substitution into the quadratic form gives

\bar{\mathbf{Q}}^\mathsf{T}\mathbf{T}_R^{-1}\bar{\mathbf{Q}} =\frac{0.110(0.20^2)-2(0.015)(0.20)(-0.10)+0.0525(-0.10)^2}{0.00555} =\frac{0.005525}{0.00555}=0.9955.

Thus $F_{\mathrm{D1}}=0.9955/2=0.4977$. For the degrees of freedom, the inverse of the within-imputation matrix is

\bar{\mathbf{U}}_R^{-1}=\frac{1}{0.0035} \begin{pmatrix}0.090&-0.010\\-0.010&0.040\end{pmatrix}.

The trace is $\operatorname{tr}(\mathbf{B}_R\bar{\mathbf{U}}_R^{-1})=0.4171$, so $r_1=1.25(0.4171)/2=0.2607$ and $\nu=4+3/0.2607^2=48.1$. The p-value is approximately $0.611$; at the 0.05 level, this joint test does not reject the null. The off-diagonal covariance terms were included in both the inverse total covariance and the trace calculation.

Comparison, Diagnostics, and Limitations

Separate scalar tests address separate one-coefficient hypotheses; they are not a substitute for the D1 test of a specified joint null. In particular, separate tests do not use the same joint covariance calculation and do not directly answer whether all the restrictions hold together. An ordinary complete-data Wald test also does not account for between-imputation uncertainty. The D2 method for pooling test statistics and the D3 method for pooling likelihood-ratio tests are other named approaches; they use different pooled testing information and are not interchangeable calculations of the D1 statistic.

A computation may be going wrong if the resulting covariance matrix is asymmetric beyond numerical tolerance, has negative eigenvalues beyond rounding error, or is singular when the formula requires an inverse. These problems can reflect too few imputations to estimate between-imputation covariance stably, redundant restrictions, or unstable model fits in some completed data sets. Review the restriction matrix and the per-imputation fits rather than forcing an inverse or silently dropping a restriction. If $r_1$ is close to zero, the formula gives very large denominator degrees of freedom; this indicates little estimated between-imputation contribution, not proof that missingness has no effect.

D1 is an approximate test, not an exact finite-sample result for every model. Its reliability depends on a meaningful common estimand across imputations, valid covariance estimates, and a sufficiently regular distribution of the pooled estimates. It cannot repair a misspecified analysis model, invalid imputation model, or a poorly identified coefficient. Report the joint null, the D1 statistic, its numerator and denominator degrees of freedom, the p-value, and the tested coefficients or contrasts in context.

Worked Example: A Clear Failure Signal

Suppose a proposed two-restriction test produces $\bar{\mathbf{U}}_R=\begin{pmatrix}0.04&0.02\\0.02&0.01\end{pmatrix}$. Its determinant is $0.04(0.01)-0.02^2=0$, so it is singular. The D1 relative-increase calculation requires $\bar{\mathbf{U}}_R^{-1}$, which does not exist. It would be incorrect to report a D1 p-value by treating this matrix as invertible.

The appropriate response is to investigate why the tested estimates are linearly dependent or why the fitted covariance is degenerate. If the restrictions are redundant, reformulate them as a full-rank set representing the intended null; if the degeneracy comes from unstable fits, address that issue before interpreting a joint test. A numerical generalized inverse does not, by itself, establish that the requested D1 test is valid.

Common Mistake

A common mistake is to calculate each coefficient’s pooled variance, ignore all covariances, and add the squared coefficient-to-standard-error ratios. That is only the right quadratic form when the tested total covariance is diagonal in the coordinates being used. In general the D1 calculation is $\bar{\mathbf{Q}}^\mathsf{T}\mathbf{T}_R^{-1}\bar{\mathbf{Q}}$, with the full matrix inverse. Another mistake is to interpret a non-significant joint p-value as proof that every tested coefficient is zero. The correct interpretation is that the data do not provide sufficient evidence against the stated joint null at the chosen significance level.

Key takeaway: D1 tests $q$ linear restrictions by comparing their pooled displacement with the full pooled covariance matrix, then uses an approximate $F_{q,\nu}$ reference distribution. The restriction rank, covariance matrices, and between-imputation contribution all matter.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

What does $q$ count in a D1 test, and what must be true of the rows of $\mathbf{R}$?
If $m=9$, $\operatorname{tr}(\mathbf{B}_R\bar{\mathbf{U}}_R^{-1})=0.60$, and $q=3$, calculate $r_1$.
Why can two separate scalar tests fail to answer the same question as a D1 joint test?
For $F_{\mathrm{D1}}=2.00$, $q=2$, and $\nu=50$, the upper-tail probability is approximately $0.146$. State the conclusion at $\alpha=0.05$ without overinterpreting it.
What should be checked if $\bar{\mathbf{U}}_R$ is singular?

Solutions. (1) $q$ is the number of independent linear restrictions; $\mathbf{R}$ must have rank $q$. (2) $r_1=(1+1/9)(0.60)/3=0.2222$. (3) Separate tests do not use the joint covariance quadratic form and test individual hypotheses rather than the specified joint null. (4) Since $0.146>0.05$, do not reject the joint null; this is insufficient evidence against it, not proof that it is true. (5) Check for redundant restrictions, degeneracy, and unstable fits; the required inverse does not exist, so do not report the D1 calculation as if it were valid.

What You'll Learn