Cohen's d is a way to describe the effect size relative to the standard deviation of the data.

For instance in the case of the difference between the means of two populations

$$\begin{array}{} \text{absolute effect size} &=& \bar{x_1} - \bar{x_2} \\ \text{relative effect size} &=& \frac{\bar{x_1} - \bar{x_2}}{\hat\sigma} &=& \text{Cohen's d} \\ \end{array}$$

If we want to apply Rubin's Rules to pool the results of multiple types of imputation of the same data, then should we apply the rules to the absolute effect size or to cohen's d?

Say, we have the following two approaches. We have some data and different imputations of it. We could do one of the following two:

  • Compute the different $d$ values and their standard error for the different imputations and apply the Rubin's Rules to it to get a pooled $d$ and it's variance.
  • Compute the different absolute effects $\bar{x_1} - \bar{x_2}$ and the different population variance estimates $\hat{\sigma}$ for the different imputations and apply the Rubin's Rules to each seperately. From those two results compute a pooled $d$ and it's variance.

Can the second approach be done (or maybe some other way of applying the rules to the absolute effects) and could it be more accurate?


1 Answer 1


Below are results for a comparison made by 1000 simulations (R-code below)

The steps taken are:

  • generate two vectors We use two samples of sizes 20 and normal distributions $X \sim N(0,1)$ and $Y \sim N(0.2,1)$ For this case Cohen's d would be equal to 0.2

  • randomly delete data We selected randomly 10 values out of the 40 values.

  • impute the data For this I used the mice r-package. The standard method that it uses is predictive mean matching.

  • compute Cohen's d This is done in multiple ways.

    1. In one case I apply Rubin's rule to the $d$ computed for each imputation.
    2. In another case, I apply Rubin's rule to the difference of the means and to the estimate of the standard deviation. Finally, to compute $d$ I take the ratio of those two.
    3. A third case. In this simplified case imputation is not at all needed. We can compute Cohen's d for two samples of unequal size and there is no need to fill the missing data.

    Imputing the data is usefull when the method can not deal with missing data (for instance a linear regression where some of the regressor values $X$ are missing). Also possible might be that one wishes to create more balanced data and use imputation to generate equal representation of classes in both variables (e.g. demographic data), but in this case there might be potentially alternatives by applying different weights based on the classes or by using a linear mixed effects model.

The results are

Method           Mean of 1000 simulations  Mean of squared error
1 relative       0.2124757                 0.1502107
2 absolute       0.2118989                 0.1493271
3 no imputation  0.2033886                 0.1302398



nd <- 20
nimp <- 10
effect <- 0.2
m <- 5

### function to create data with imputations from mice
set <- function(imp, n) {
  data_imp <- imp$data
  data_imp[imp$where[,1],1] <- imp$imp[[1]][,n]
  data_imp[imp$where[,2],2] <- imp$imp[[2]][,n]

### function to do the simulation
sim <- function() {
  ### generate data
  data <- data.frame(x = rnorm(nd,0,1), y = rnorm(nd,effect,1))

  ### sample nimp points to remove (make equal to NA)
  imp <- sample(0:(nd*2-1),nimp)
  for (i in imp) {
    k <- floor(i/30)+1
    l <- (i%%30)+1
    data[l,k] <- NA
  ### perform 5 imputations with mice
  imp2 <- mice(data[], m = m, printFlag = F)
  ### vectors to contain results from for-loop
  ve_d <- rep(0,5)
  sd_d <- rep(0,5)
  ve_m <- rep(0,5)
  sd_m <- rep(0,5)
  ve_v <- rep(0,5)
  sd_v <- rep(0,5)
  for (i in 1:m) {  ### repeatedly compute statistics for different imputations
    m1 <- set(imp2,i)
    par1 <- mean(m1[,2])-mean(m1[,1])            ### difference of means
    par2 <- sqrt(0.5*(var(m1[,1])+var(m1[,2])))  ### pooled variance estimate
    par3 <- par1/par2                            ### Cohen's d
    nu <- 2*nd-2
    ### standard deviations of par1, par2, par3
    sd1 <- par2*sqrt(2/nd)
    sd2 <- par2 * sqrt(gamma(nu/2)/gamma((nu+1)/2) * nu/2- 1)
    sd3 <- sqrt((nd+nd)/(nd*nd) * (nu)/(nu-2)+ par3^2 * (nu/(nu-2)-(1-3/(4*nu-1))^-2))
    ve_m[i] <- par1
    sd_m[i] <- sd1
    ve_v[i] <- par2
    sd_v[i] <- sd2
    ve_d[i] <- par3
    sd_d[i] <- sd3    
  ### Rubin's Rule applied to Cohen's d
  est_d <- mean(ve_d)
  var_d <- mean(sd_d^2) + var(ve_d) * (1+1/m)

  ### Rubin's Rule applied to difference
  est_m <- mean(ve_m)
  var_m <- mean(sd_m^2) + var(ve_m) * (1+1/m)
  ### Rubin's Rule applied to variance
  est_v <- mean(ve_v)
  var_v <- mean(sd_v^2) + var(ve_v) * (1+1/m)
  ### Straightforward Cohen's d
  x <- data$x[!is.na(data$x)]
  y <- data$y[!is.na(data$y)]
  nx <- length(x)
  ny <- length(y)
  cohend <- (mean(y)-mean(x))/sqrt( ((nx-1)*var(x)+ (ny-1)*var(y))/(nx+ny-2))

  ### return the results
  r <- list(est_d = est_d, var_d = var_d,
            est_m = est_m, var_m = var_m,
            est_v = est_v, var_v = var_v,
            cohend = cohend)

results <- replicate(10^3,sim(), simplify = TRUE)

# mean
# variance

# error

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.