Rubin's Combining Rules and Variance Decomposition

Why Pooling Requires More Than Averaging

Multiple imputation produces \(M\) completed data sets, not one definitive completed data set. Each gives an estimate of the same parameter, along with an estimate of its sampling covariance if that completed data set had been observed. The estimates can differ because the imputations represent uncertainty about the missing values. A valid summary must account for both the uncertainty within each completed data set and the variability between completed data sets.

The rules below combine these two sources and also account for the fact that only \(M\) imputations are used. They presume proper imputations and a common estimand across completed data sets, as in the earlier material on proper multiple imputation. For the usual variance and inference approximations, the complete-data estimates should also be approximately unbiased and normally distributed, with covariance estimates that estimate their sampling variances adequately.

Scalar Combining Rules

The factor multiplying \(B\) has two parts. The \(B\) term represents uncertainty associated with the missing values. The additional \(B/M\) term allows for the variability in the pooled estimate caused by averaging only \(M\) imputed estimates. Thus, the total variance is not generally just \(W+B\), and it is not \(W+B/M\).

Why the Total Variance Has This Form

The decomposition is an approximation under the usual large-sample multiple-imputation framework, not an exact identity for every finite collection of completed data sets. Its logic is a variance decomposition: uncertainty within completed data sets contributes \(W\); variation in the estimates across imputations contributes \(B\); and estimating the pooled result with a finite number of imputations contributes an additional \(B/M\).

The finite-imputation term can be understood by first imagining an unlimited number of imputations. In that limit, averaging removes Monte Carlo error from the finite set of imputations, but the underlying uncertainty due to the unobserved values remains. With only \(M\) imputations, the average of the imputation-specific estimates has additional simulation variability.

To see the size of that contribution, write the imputation-specific estimates, conditional on the observed data, as draws with a common mean and variance \(V\). The average of \(M\) independent draws then has variance \(V/M\): expanding the variance of the sum gives \(M V\), since cross-covariances are zero, and dividing by \(M^2\) gives \(V/M\). The sample variance \(B\) estimates the across-imputation variance \(V\), so the finite-imputation contribution is estimated by \(B/M\). Adding this to the within-imputation term and the between-imputation uncertainty gives \(W+B+B/M\), the stated formula. In practice, \(W\) and \(B\) are themselves estimated from the completed data sets, so the rule is used as a variance estimator under its stated approximation conditions.

If \(B=0\) and \(W>0\), the estimates are identical across the \(M\) imputations, so \(r=0\); the limiting degrees-of-freedom approximation is infinite. This does not establish that missing-data uncertainty is truly absent: it may simply be too small to detect with the available imputations.

Vector Estimates and Covariance Matrices

For a vector parameter, each completed data set supplies an estimate vector \(\hat{\boldsymbol Q}_m\) and an estimated covariance matrix \(\boldsymbol U_m\). The pooled estimate is the componentwise average. The within-imputation covariance averages the covariance matrices, and the between-imputation covariance records how the estimate vectors vary jointly.

The vector rule follows by applying the scalar decomposition to every linear combination of the parameters. For any fixed vector \(\boldsymbol c\), the scalar estimate \(\boldsymbol c^\mathsf{T}\hat{\boldsymbol Q}_m\) has within-imputation variance \(\boldsymbol c^\mathsf{T}\boldsymbol U_m\boldsymbol c\) and between-imputation variance \(\boldsymbol c^\mathsf{T}\boldsymbol B\boldsymbol c\). Its total variance is therefore \( \boldsymbol c^\mathsf{T}\boldsymbol W\boldsymbol c+(1+1/M)\boldsymbol c^\mathsf{T}\boldsymbol B\boldsymbol c =\boldsymbol c^\mathsf{T}\boldsymbol T\boldsymbol c. \) This verifies that \(\boldsymbol T\) gives the Rubin variance rule for every linear contrast, not just for individual coordinates. Because \(\boldsymbol W\) and \(\boldsymbol B\) are covariance matrices, each is positive semidefinite; their nonnegative weighted sum \(\boldsymbol T\) is also positive semidefinite.

More Imputations and a Common Pooling Error

Increasing \(M\) reduces the Monte Carlo term \(B/M\), but it does not remove the between-imputation uncertainty \(B\). In the limit as \(M\) grows, the total variance approaches \(W+B\). Thus, adding imputations can stabilize estimates of the pooled quantities and improve Monte Carlo precision, but it does not make missing information disappear.

A common mistake is to average the estimates and then report the average complete-data standard error, \(\sqrt{W}\), as if the imputed values had been known. That calculation ignores between-imputation variability. Another mistake is to report \(\sqrt{W+B}\) for a finite \(M\), omitting the extra \(B/M\). Use \(W+(1+1/M)B\) for the total variance, and distinguish the variance \(T\) from its square root, the standard error. For vector parameters, apply the same rule to the full covariance matrices rather than pooling only their diagonal entries when a joint contrast is required.