I am using the mice package in R to impute missing data in small study. The study investigates the effect of a behavioral intervention on the frequency of a particular behavior, i.e., count data that can range from anywhere between 0-5 times for most participants and several hundreds for some participants. I have used a package for imputation of count data (https://github.com/kkleinke/countimp) and have successfully imputed the missing data with a reasonable distribution.
Now I want to report the estimated frequency at pre-treatment and post-treatment as well as the pre- to post-change and a significance test of the change, all based on the pooled data. The convention in the field is to report medians and IQR of pre- and post-frequencies as well as the median change in frequency.
Can this be properly done with multiply imputed data? My current idea is to calculate the medians and IQR with standard errors using the quantile regression library for R (quantreg), by employing intercept only regression models for quantiles .25, .50, .75 for each imputation (pre, post, and change score), and then pooling these estimates with standard errors using mice. However, this essentially means that I am taking the mean of medians, which does not feel right.
I'm pasting code that calculates the median (50th percentile) of the pre- and post-frequencies as well as the change in frequency using quantreg and then pools the estimates using mice. After that are calculations of the corresponding means.
# demo.median contains the imputed data # vcov.rq is called by pool to get the covariance estimate of the # quantile regression vcov.rq <- function(fit)summary(fit, se="nid", cov=T)$cov l <- list() l$pre <- summary(pool(with(demo.median, rq(x_PRE ~ 1, tau = .5)))) l$post <- summary(pool(with(demo.median, rq(x_POST ~ 1, tau = .5)))) l$pre_post <- summary(pool(with(demo.median, rq(I(x_PRE - x_POST) ~ 1, tau = .5)))) # Print pooled medians of pre, post and change score do.call(rbind.data.frame, l) l <- list() l$pre <- summary(pool(with(demo.median, lm(x_PRE ~ 1)))) l$post <- summary(pool(with(demo.median, lm(x_POST ~ 1)))) l$pre_post <- summary(pool(with(demo.median, lm(I(x_PRE - x_POST) ~ 1)))) # Print pooled means of pre, post and change score do.call(rbind.data.frame, l)
# Medians est se t df Pr(>|t|) lo 95 hi 95 nmis fmi lambda pre 8.000 1.3250825 6.037360 91.75287 3.286786e-08 5.3681774 10.631823 NA 0.02110754 0.0000000 post 2.875 1.0997478 2.614236 76.09011 1.077545e-02 0.6847043 5.065296 NA 0.16031744 0.1385330 pre_post 2.400 0.8801746 2.726732 52.96378 8.654540e-03 0.6345650 4.165435 NA 0.34977967 0.3256813 # Means est se t df Pr(>|t|) lo 95 hi 95 nmis fmi lambda pre 18.252632 3.417037 5.341655 92.05265 6.619670e-07 11.466153 25.03911 NA 0.02104097 0.000000000 post 12.019474 2.413928 4.979217 91.32515 2.999330e-06 7.224731 16.81422 NA 0.02889543 0.007858808 pre_post 6.233158 2.316257 2.691047 91.26050 8.471778e-03 1.632374 10.83394 NA 0.02957225 0.008535556
My question is: Is this proper? If not, can I achieve my goal of estimating absolute and change score medians from multiply imputed data in any other way?