Worked D1 Test For A Joint Regression Hypothesis

Why a Joint Test Needs Pooling

A regression model may include several coefficients whose joint relevance is of interest. For example, a researcher might ask whether two treatment indicators are both zero, or whether a group of covariates contributes to a model. With multiple imputation, each completed data set produces its own coefficient estimates and covariance matrix. A joint test must account for both the uncertainty within each completed-data analysis and the variation among the estimates across imputations.

The D1 method is a multivariate Wald test for a set of linear restrictions. It pools the relevant coefficient estimates and their covariance information, then uses an approximate F reference distribution. It is not the average of the separate tests from each imputed data set.

Definition: Let $m$ be the number of imputations and let $\boldsymbol q_i$ be the $k$-dimensional estimate of the quantities being tested in imputation $i$, with estimated covariance matrix $\boldsymbol U_i$. Define the pooled estimate $\bar{\boldsymbol q}$, within-imputation covariance $\bar{\boldsymbol U}$, and between-imputation covariance $\boldsymbol B$ by $$ \bar{\boldsymbol q}=\frac{1}{m}\sum_{i=1}^m\boldsymbol q_i,\qquad \bar{\boldsymbol U}=\frac{1}{m}\sum_{i=1}^m\boldsymbol U_i,\qquad \boldsymbol B=\frac{1}{m-1}\sum_{i=1}^m(\boldsymbol q_i-\bar{\boldsymbol q})(\boldsymbol q_i-\bar{\boldsymbol q})^\mathsf T. $$ For the null hypothesis $H_0:\boldsymbol q=\boldsymbol 0$, the D1 relative increase in variance is $r_1=(1+1/m)\operatorname{tr}(\boldsymbol B\bar{\boldsymbol U}^{-1})/k$. The test statistic is $F_{\mathrm{D1}}=\bar{\boldsymbol q}^\mathsf T\bar{\boldsymbol U}^{-1}\bar{\boldsymbol q}/\{k(1+r_1)\}$, with numerator degrees of freedom $df_1=k$ and approximate denominator degrees of freedom $df_2=4+\{(k+1)/(k+3)\}k/r_1^2$, when $r_1>0$.

Setting Up the Joint Hypothesis

Suppose the regression coefficient vector in imputation $i$ is $\hat{\boldsymbol\beta}_i$, with estimated covariance matrix $\boldsymbol V_i$. A joint linear hypothesis can be written $H_0:\boldsymbol R\boldsymbol\beta=\boldsymbol c$, where $\boldsymbol R$ has $k$ rows. For each imputation, form the restricted estimate and covariance $ \boldsymbol q_i=\boldsymbol R\hat{\boldsymbol\beta}_i-\boldsymbol c $ and $ \boldsymbol U_i=\boldsymbol R\boldsymbol V_i\boldsymbol R^\mathsf T $. The test is then $H_0:\boldsymbol q=\boldsymbol 0$. This setup ensures that each imputation tests the same restrictions on the same scale.

The usual D1 calculation requires that the estimates and covariance matrices refer to the same estimand and parameterization in every imputation. The within-imputation covariance matrices should be valid and nonsingular for the tested quantities, and the imputation-specific estimates should be approximately multivariate normal for the Wald approximation to be reasonable. The multiple-imputation procedure must also be suitable for the missingness and analysis models; a pooling formula cannot repair biased imputations or a badly specified regression model.

The D1 approximation can be understood by comparing $\boldsymbol B$ with $\bar{\boldsymbol U}$. The matrix $(1+1/m)\boldsymbol B$ represents the added variance associated with uncertainty across imputations. The eigenvalues of $\bar{\boldsymbol U}^{-1/2}(1+1/m)\boldsymbol B\bar{\boldsymbol U}^{-1/2}$ describe the relative increase in different directions of the joint hypothesis. Their average is $$ \frac{1}{k}\operatorname{tr}\!\left\{\bar{\boldsymbol U}^{-1/2}(1+1/m)\boldsymbol B\bar{\boldsymbol U}^{-1/2}\right\} =\frac{1+1/m}{k}\operatorname{tr}(\boldsymbol B\bar{\boldsymbol U}^{-1})=r_1. $$

D1 uses this average increase as a scalar adjustment. The Wald quadratic form based on $\bar{\boldsymbol U}$ is divided by $1+r_1$, so the added variance reduces the test statistic. The denominator degrees of freedom approximate the additional uncertainty in estimating that inflation. This is an approximation: a single average cannot retain every direction-specific feature of a multivariate variance increase.

Worked Example: A One-Parameter D1 Test

Suppose three hypothetical imputations give estimates of one regression coefficient of $0.8$, $1.0$, and $1.2$, with estimated variance $0.25$ in each imputation. Test $H_0:\beta=0$. The same coefficient and variance scale are used in all three analyses, and the estimates are treated as approximately normal.

Here $m=3$ and $k=1$. The pooled estimate is $\bar q=(0.8+1.0+1.2)/3=1.0$, and $\bar U=(0.25+0.25+0.25)/3=0.25$. The between-imputation variance is $B=\{(0.8-1.0)^2+(1.0-1.0)^2+(1.2-1.0)^2\}/(3-1)=0.04$.

Therefore $r_1=(1+1/3)(0.04/0.25)=0.2133$. The statistic and degrees of freedom are $$ F_{\mathrm{D1}}=\frac{1.0^2/0.25}{1+0.2133}=3.297,\qquad df_1=1,\qquad df_2=4+\frac{2}{4}\frac{1}{0.2133^2}\approx14.99. $$

Against an $F$ distribution with 1 and approximately 14.99 degrees of freedom, the two-sided one-parameter test has $p\approx0.090$. At significance level $0.05$, this is not sufficient evidence against $H_0$. The pooled coefficient is $1.0$ coefficient units above zero, but this test does not establish that the population coefficient is zero.

A Joint Test With Corrected Between-Imputation Covariance

The following hypothetical example uses two coefficients and two imputations. It also illustrates why every element of the between-imputation covariance matrix must be calculated from the imputation-specific vectors, rather than inferred from a verbal description.

Worked Example: Testing Two Regression Coefficients Together

Suppose the two coefficient estimates are $\boldsymbol q_1=(1.5+\sqrt{10}/2,\ 1.5-1/\sqrt{10})^\mathsf T$ and $\boldsymbol q_2=(1.5-\sqrt{10}/2,\ 1.5+1/\sqrt{10})^\mathsf T$. In both imputations, suppose the covariance matrix for these estimates is the $2$ by $2$ identity matrix. Test $H_0:\boldsymbol q=\boldsymbol 0$. These are hypothetical values, and the coefficients are assumed to be on the same scale in both imputations.

There are $m=2$ imputations and $k=2$ restrictions. Their mean is $\bar{\boldsymbol q}=(1.5,1.5)^\mathsf T$, and $\bar{\boldsymbol U}=\boldsymbol I$. The deviations from the mean are $(\sqrt{10}/2,-1/\sqrt{10})^\mathsf T$ and its negative. Thus $$ \boldsymbol B= \begin{pmatrix} 5&-1\\ -1&0.2 \end{pmatrix}. $$ For example, $B_{11}=2(\sqrt{10}/2)^2=5$, while $B_{12}=2(\sqrt{10}/2)(-1/\sqrt{10})=-1$. The second diagonal element is $B_{22}=2(1/\sqrt{10})^2=0.2$. The matrix is positive semidefinite: its determinant is $5(0.2)-(-1)^2=0$, and its diagonal entries are nonnegative.

Since $\bar{\boldsymbol U}^{-1}=\boldsymbol I$, the trace in $r_1$ is $5+0.2=5.2$. Hence $r_1=(1+1/2)(5.2/2)=3.9$. The pooled Wald quadratic form is $\bar{\boldsymbol q}^\mathsf T\bar{\boldsymbol U}^{-1}\bar{\boldsymbol q}=1.5^2+1.5^2=4.5$. Therefore $$ F_{\mathrm{D1}}=\frac{4.5}{2(1+3.9)}=0.4592,\qquad df_1=2,\qquad df_2=4+\frac{3}{5}\frac{2}{3.9^2}\approx4.079. $$

The upper-tail $p$-value for $F(2,4.079)$ at $0.4592$ is approximately $0.661$. The data therefore do not provide evidence against the joint null at the $0.05$ level. The large $r_1$ reflects substantial between-imputation variation relative to within-imputation covariance; it increases the variance scale and consequently reduces the D1 statistic.

Accounting for Correlated Coefficients

The off-diagonal elements of the covariance matrices matter. A joint Wald statistic uses the full inverse covariance matrix, not just separate standard errors. The next example shows the computation when both within-imputation and between-imputation covariance have nonzero off-diagonal terms.

Worked Example: A Joint Test With Correlated Estimates

For three hypothetical imputations, suppose the estimate vectors are $(0.4,0.2)^\mathsf T$, $(0.5,0.1)^\mathsf T$, and $(0.3,0.3)^\mathsf T$. In each imputation, suppose the covariance matrix is $\begin{pmatrix}0.16&0.06\\0.06&0.09\end{pmatrix}$. Test the joint null that both coefficients are zero. The covariance matrix is positive definite because its determinant is $0.16(0.09)-0.06^2=0.0108>0$.

The pooled estimate is $\bar{\boldsymbol q}=(0.4,0.2)^\mathsf T$, and $\bar{\boldsymbol U}$ is the stated covariance matrix. The estimate deviations are $(0,0)^\mathsf T$, $(0.1,-0.1)^\mathsf T$, and $(-0.1,0.1)^\mathsf T$. Dividing the sum of their outer products by $m-1=2$ gives $\boldsymbol B=\begin{pmatrix}0.01&-0.01\\-0.01&0.01\end{pmatrix}$.

The inverse within-imputation covariance is approximately $\bar{\boldsymbol U}^{-1}=\begin{pmatrix}8.333&-5.556\\-5.556&14.815\end{pmatrix}$. The trace product is $\operatorname{tr}(\boldsymbol B\bar{\boldsymbol U}^{-1})=0.3426$, so $r_1=(1+1/3)(0.3426/2)=0.2284$. The Wald quadratic form is approximately $(0.4,0.2)\bar{\boldsymbol U}^{-1}(0.4,0.2)^\mathsf T=1.037$. Thus $$ F_{\mathrm{D1}}=\frac{1.037}{2(1+0.2284)}=0.4221,\qquad df_1=2,\qquad df_2=4+\frac{3}{5}\frac{2}{0.2284^2}\approx27.00. $$

The upper-tail $p$-value for $F(2,27.00)$ is approximately $0.660$. At the $0.05$ level, the joint null is not rejected. Ignoring the off-diagonal entries would change the quadratic form and could change the conclusion; the covariance structure is part of the joint test.

Diagnostics, Alternatives, and Limitations

Before interpreting a D1 result, check the inputs and the approximation. Each imputation must provide the same $k$ estimates in the same order, together with their corresponding covariance matrices. The pooled within-imputation covariance must be invertible. The between-imputation covariance should be symmetric and positive semidefinite, apart from tiny numerical discrepancies. Investigate coding or convergence problems if estimates are extreme, if covariance matrices are singular, or if between-imputation variation is unexpectedly large.

D1 is most useful when the analysis supplies coefficient estimates and covariance matrices and the target is a joint Wald hypothesis. Its approximate F reference distribution is not an exact finite-sample result. In particular, a small number of imputations or a large relative increase in variance can make the degrees-of-freedom approximation consequential. The formula for $df_2$ requires $r_1>0$; if the between-imputation covariance is exactly zero, the formula is undefined. In that limiting situation there is no estimated between-imputation variance to adjust for, and inference should use the appropriate complete-data test reference distribution rather than inserting $r_1=0$ into the displayed expression.

D1 differs from approaches that combine test statistics rather than pooled coefficient and covariance information. D2 combines imputation-specific Wald-type statistics using a separate combining rule, while D3 is based on likelihood-ratio information. These methods are not interchangeable by definition: their inputs, approximations, and suitability depend on the model and available output. Averaging imputation-specific $p$-values is not a substitute for a pooled test, because $p$-values discard the estimates and covariance information needed to represent pooled uncertainty.

Common pitfall: It is tempting to read $1+r_1$ as an inflation of the Wald statistic. The formula places this factor in the denominator: it inflates the variance scale for the extra uncertainty and therefore reduces the test statistic. Also, do not substitute separately pooled standard errors for the full covariance matrices in a joint test; the off-diagonal covariances contribute to the result.

Check Your Understanding

Use the D1 definitions above. The solutions are included so that each calculation can be checked.

If a test concerns three linear restrictions, what is $df_1$? Solution: $df_1=k=3$, the number of restrictions.
For $m=4$, $k=1$, $\bar U=0.5$, and $B=0.1$, calculate $r_1$. Solution: $r_1=(1+1/4)(0.1/0.5)=1.25(0.2)=0.25$.
In that same one-parameter case, if $\bar q=1$, calculate the D1 statistic. Solution: The unadjusted quadratic form is $1^2/0.5=2$. Therefore $F_{\mathrm{D1}}=2/(1+0.25)=1.60$.
What happens to $F_{\mathrm{D1}}$ when $r_1$ increases while the pooled estimate and within-imputation covariance stay fixed? Solution: The denominator $k(1+r_1)$ increases, so the statistic decreases. This reflects greater uncertainty from between-imputation variation.
Why is averaging the imputation-specific $p$-values not the D1 procedure? Solution: A $p$-value alone does not retain the pooled coefficient vector or its within- and between-imputation covariance. D1 calculates its statistic from that pooled estimate and covariance information, then uses an approximate F reference distribution.