Hello Everybody,
I do a lot of imputing for applied researchers. I use the excellent
norm library in R. Now I wish to perform ANOVAs (regressions on dummies)
on imputed files and have to pool the multiple outcomes. Unfortunately,
unlike the Windows version of Norm, the R library does not contain a
procedure for the multiparameter method. I am working on an R-version of
the procedure myself. To that end I have been studying the Li,
Raghunathan, & Rubin (1991) JASA article, the 1997 book by Schafer (par
4.3.3), and the comprehensive help function in the Windows version of
Norm. To test my own procedure I performed a multiparameter MI-inference
within the Windows version of Norm as well.
First I had to overcome several difficulties:
A Schafer's book, page 114 and the help function in Norm differ in
the reported estimate of the total covariance matrix. The book seems to
contain an error, it should be (in tex code) \tilde{T}=(1+r_1)\bar{U}
instead of \tilde{T}=(1-r_1)\bar{U} (so, the minus should be a plus).
B In the help function, the r=(1+m^-1)tr(B\bar{U} should also be
divided by k.
Then, the outcomes of my procedure (code is added in the appendices)
and Windows-Norm' differed. Moreover, when using Meng & Rubins (1992) D3
statistic, which I also programmed, the outcomes were very similar to
those of Norm in windows. After staring at my R-code for hours I still
think it is okay; I followed the Li et al formula's, but the code must
be wrong somehow.
My question is the following:
Does anybody see a flaw in the R-code for the mimd.inferende function
(mi-inference for multidimensional estimands) below? In addition, does
anybody have an example for which the multiparameter is valid so I can
test my R-procedure?
I imputed the data 5 times. On each of the files a regression model with
two dummies (three groups) was estimated.
##############
Appendix: R-code
##############
mimd.inference<-function (est, cov.mat)
{
qstar <- est[[1]]
for (i in 2:length(est)) {
qstar <- cbind(qstar, est[[i]])
}
qbar <- apply(qstar, 1, mean)
u <- as.vector(cov.mat[[1]])
for (i in 2:length(cov.mat)) {
u <- cbind(u, as.vector(cov.mat[[i]]))
}
ubar <- matrix(apply(u, 1, mean),nrow=length(est[[1]]))
m <- dim(qstar)[2]
ss <- (est[[1]] - qbar) %*% (est[[1]] - qbar)
for (i in 2:m) {
ss <- cbind(ss,(est[[i]] - qbar) %*% (est[[i]] - qbar))
}
bm <- sum(ss)/(m-1)
rem <- (1 + 1/m) * sum(diag(bm/ubar))/k
tm <- (1 + rem) * ubar
k <- nrow(ubar)
t <- k * (m-1)
if (t > 4)
nu1 <- 4 + (t-4) * (1 + (1-2/t)/rem)^2 else
nu1 <- 0.5 * t * (1 + 1/k) * (1 + 1/rem)^2
D1 <- (t(t(1/tm)%*%qbar)%*%qbar) * 1/k
pval <- 1 - pf(D1, k, nu1)
fminf <- rem/(rem + 1)
result <- list(est = qbar, cov.mat = tm, df1 = k, df2 = nu1, signif
= pval,
r = rem, fminf = fminf)
result
}
##############
Appendix: Outcomes in R
##############
> est<-list(c(-4.553825, -3.815103), c(-3.206370,
-3.818416),c(-3.960845,-4.248486
),c(-2.855024,-3.652531),c(-0.8009951,-2.5596872))
> cov.mat<-list(
+ matrix(c(1.6526685, 0.8221818, 0.8221818, 1.6282424),nrow=2),
+ matrix(c(1.678708, 0.835136, 0.835136, 1.653897),nrow=2),
+ matrix(c(1.6380543, 0.8149114, 0.8149114, 1.6138442),nrow=2),
+ matrix(c(1.6175106, 0.8046912, 0.8046912, 1.5936041),nrow=2),
+ matrix(c(1.5615181, 0.7768356, 0.7768356, 1.5384392),nrow=2))
>
> mimd.inference(est,cov.mat)
$est
[1] -3.075412 -3.618845
$cov.mat
[,1] [,2]
[1,] 4.593598 2.285257
[2,] 2.285257 4.525705
$df1
[1] 2
$df2
[1] 11.97932
$signif
[,1]
[1,] 0.00827624
$r
[1] 1.818691
$fminf
[1] 0.6452254
##############
Appendix: multiple parameters for Norm in Windows
##############
Here follow the estimates of the 5 files. For each of the files, the
first row contain the parameter estimates of the two dummies; the next
two rows contain their covariance matrix.
-4.5538251 -3.8151025
1.6526685 0.8221818
0.8221818 1.6282424
-3.2063700 -3.8184156
1.6787079 0.8351361
0.8351361 1.6538969
-3.9608448 -4.2484857
1.6380543 0.8149114
0.8149114 1.6138442
-2.8550241 -3.6525309
1.6175106 0.8046912
0.8046912 1.5936041
-0.8009951 -2.5596872
1.5615181 0.7768356
0.7768356 1.5384392
#########
Appendix: Outcomes of Norm-Windows
#########
**************************************************
Est. S.E.
QTY_1 -3.07541 1.71355
QTY_2 -3.61884 1.70084
COVARIANCE MATRIX:
QTY_1 QTY_2
QTY_1 2.93626 1.46075
QTY_2 1.46075 2.89286
AVG.REL.INCR. IN VARIANCE DUE TO NONRESPONSE = 0.801728
AVG. RATE OF MISSING INFORMATION = 44.50%
TESTING H0: ALL ESTIMANDS EQUAL TO ZERO
F = 2.61774
DF_1 = 2
DF_2 = 18
PVAL = 0.100482
--
Niels Smits
Research Methodology,
Statistics and Data-analysis
Faculty of Psychology and Education
Free University Amsterdam
Van der Boechorststraat 1
1081 BT Amsterdam
The Netherlands
Tel: +31 (0)20 5988713
Secr: +31 (0)20 5988757
Fax: +31 (0)20 5988758
