Nordheim Statistics 411 Spring 2015 (April 9 or 14) Example of Multiple Imputation Using “MICE” > library(mice) > ttm This is the data set (with missing data). x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 NA NA 6.9 2.5 5 33 0.17 12.0 1.8 6 28 NA NA 2.4 7 50 0.18 8.5 3.7 8 NA 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 NA 0.09 9.1 3.1 13 NA 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 NA 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 NA 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 NA 2.4 > md.pattern(ttm) This describes the pattern of missing data. y x3 x2 x1 12 1 1 1 1 0 4 1 1 1 0 1 1 1 1 0 1 1 1 1 0 1 1 1 1 1 1 0 0 2 1 1 0 0 1 2 0 2 3 5 10 > imp=mice(ttm) This performs the imputation (using defaults for all variables). > imp Number of multiple imputations: 5 Missing cells per column: x1 x2 x3 y 5 3 2 0 Imputation methods: x1 x2 x3 y "pmm" "pmm" "pmm" "" The default imputation method is “predictive mean matching”. > complete(imp,1) This is the first set of imputed values. x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 32 0.11 6.9 2.5 5 33 0.17 12.0 1.8 6 28 0.20 8.1 2.4 7 50 0.18 8.5 3.7 8 28 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 50 0.09 9.1 3.1 13 31 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 0.20 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 31 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 11.2 2.4> complete(imp,2) 2nd set x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 28 0.20 6.9 2.5 5 33 0.17 12.0 1.8 6 28 0.09 5.8 2.4 7 50 0.18 8.5 3.7 8 28 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 45 0.09 9.1 3.1 13 32 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 0.20 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 32 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 9.1 2.4 > complete(imp,3) 3rd set x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 33 0.21 6.9 2.5 5 33 0.17 12.0 1.8 6 28 0.25 10.7 2.4 7 50 0.18 8.5 3.7 8 22 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 44 0.09 9.1 3.1 13 28 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 0.17 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 33 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 9.1 2.4 > complete(imp,4) 4th set x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 28 0.20 6.9 2.5 5 33 0.17 12.0 1.8 6 28 0.21 8.1 2.4 7 50 0.18 8.5 3.7 8 22 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 44 0.09 9.1 3.1 13 28 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 0.25 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 33 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 8.5 2.4> complete(imp,5) 5th set x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 32 0.20 6.9 2.5 5 33 0.17 12.0 1.8 6 28 0.23 9.2 2.4 7 50 0.18 8.5 3.7 8 32 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 40 0.09 9.1 3.1 13 33 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 0.14 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 33 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 9.9 2.4 > fit=with(imp, lm(y~x1+x2+x3)) The coefficient estimates for each data set > fit [[1]] Coefficients: (Intercept) x1 x2 x3 0.14905 0.09531 5.26476 -0.19746 Coefficients: (Intercept) x1 x2 x3 0.2325 0.1000 3.8854 -0.2022 [[3]] Coefficients: (Intercept) x1 x2 x3 0.1738 0.0942 3.9233 -0.1681 [[4]] Coefficients: (Intercept) x1 x2 x3 0.57364 0.09394 3.06969 -0.19706 [[5]] Coefficients: (Intercept) x1 x2 x3 -0.2995 0.1045 4.5086 -0.1723 > pool(fit) Averaging across the 5 data sets. Call: pool(object = fit) Pooled coefficients: (Intercept) x1 x2 x3 0.16591619 0.09758009 4.13033501 -0.18742304 Fraction of information about the coefficients missing due to nonresponse: (Intercept) x1 x2 x3 0.3628758 0.3323977 0.3412963 0.2541291 > summary(pool(fit)) est se t df Pr(>|t|) lo 95 (Intercept) 0.16591619 0.69557295 0.2385317 9.399039 8.165884e-01 -1.39745459 x1 0.09758009 0.01084748 8.9956480 10.026907 4.078506e-06 0.07341919 x2 4.13033501 1.90469969 2.1684967 9.841864 5.573497e-02 -0.12286165 x3 -0.18742304 0.04721189 -3.9698268 11.687462 1.955724e-03 -0.29059484 hi 95 nmis fmi lambda (Intercept) 1.72928698 NA 0.3628758 0.2403405 x1 0.12174098 5 0.3323977 0.2113117 x2 8.38353168 3 0.3412963 0.2197852 x3 -0.08425124 2 0.2541291 0.1365530> pool(fit2) Two additional fits Pooled coefficients: (Intercept) x1 x2 x3 0.15240851 0.09764443 3.76957000 -0.18312963 > pool(fit3) Pooled coefficients: (Intercept) x1 x2 x3 0.18978318 0.09632006 4.82807417 -0.19774438 ------------------- > tt=read.csv("imptest.csv",header=T) The original data set (with no missing data). > tt x1 x2 x3 y 1 22 0.14 7.3 1.9 2 31 0.23 11.2 2.3 3 45 0.09 5.8 3.9 4 29 0.24 6.9 2.5 5 33 0.17 12.0 1.8 6 28 0.12 6.8 2.4 7 50 0.18 8.5 3.7 8 27 0.08 9.3 1.3 9 38 0.20 13.2 2.3 10 52 0.11 10.7 4.3 11 43 0.17 7.8 3.2 12 38 0.09 9.1 3.1 13 31 0.25 8.1 2.7 14 44 0.13 9.9 2.9 15 33 0.19 10.4 2.0 16 47 0.20 6.2 4.8 17 40 0.14 7.5 3.5 18 36 0.11 8.5 2.0 19 25 0.16 9.2 1.0 20 32 0.21 11.4 2.4 > out=lm(y~x1+x2+x3,tt) The parameter estimates with the original data. > summary(out) Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 0.13481 0.62994 0.214 0.83324 x1 0.09956 0.01062 9.372 6.72e-08 *** x2 3.34201 1.78102 1.876 0.07895 . x3 -0.17521 0.04494 -3.898 0.00128 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 0.3899 on