Variance of Cohen's $d$ for within subjects designs

Question

Based on this paper by Lakens 2013, I understand there are two recommended ways for calculating Cohen's $d$ for within subjects measures. 

1. Cohen's $d_{rm}$ (which assumes the correlation between measures is known) and can be calculated as: $$ \text{Cohen's }d_{rm}=( M_{diff}/\sqrt{(SD_1^2+SD_2^2-2*r*SD_1*SD_2})*\sqrt{2(1-r)} $$ Where $M_{diff}$ is the difference in means, $SD_1$ and $SD_2$ are the standard deviations of these means and $r$ is the correlation between measures.

   the variance of Cohen's $d_{rm}$ can be calculated using the following: $$ V(d_{rm})=(1/n+d_{rm}^2/2n)2(1-r) $$ where $n$ is the sample size. Again this assumes the correlation is known. 

2. Cohen's $d_{av}$ (which ignores any correlation between measures and used the average of the standard deviations). $$ \text{Cohen's }d_{av}=M_{diff}/((SD_1+SD_2)/2) $$ My question is how to calculate the variance of Cohen's $d_{av}$ $(V(d_{av}))$? 

Can it be calculate in a similar manner to the calculation for  independent groups - simply by substituting $d$ for $d_{av}$? $$ V(d_{av}) = (n_1+n_2/n_1*n_2)+(d_{av}^2/2(n_1+n_2)) $$ if so some elements of this equation appear unclear as the sample size is the same for both observations.

Calculating variance of Cohen's d for repeated measures designs? - ResearchGate (accessed Jan 13, 2017).

Answer 1

There is no know derivation of the variance of $d_{av}$. In fact, this is not how one should compute the $d$ value for a within-subjects design (neither with $d_{rm}$ or the $d_{av}$).

There are two approaches for computing a $d$ value for a within-subjects design. The first uses change score standardization and is given by $$d_c = \frac{M_d}{SD_d},$$ where $M_d$ is the mean change and $SD_d$ is the SD of the change scores (which is equal to $SD_d = \sqrt{SD_1^2 + SD_2^2 - 2 \times r \times SD_1SD_2}$). The large-sample variance of $d_c$ is $$Var[d_c] = \frac{1}{n} + \frac{d_c^2}{2 \times n}.$$

The second approach uses raw score standardization. Here, we standardize based on the SD of either the pre- or the post-test scores (typically, the SD of the pre-test scores is used). So, we compute $$d_r = \frac{M_d}{SD_1}.$$ The large-sample variance of $d_r$ is $$Var[d_r] = \frac{2(1-r)}{n} + \frac{d_r^2}{2 \times n}.$$

We do not take the average of $SD_1$ and $SD_2$ in the second approach (and if anything, we would compute $\sqrt{(SD_1^2 + SD_2^2)/2}$ but again we also do not do that). If you would do this, then it becomes nearly impossible to derive the variance, because $SD_1^2 + SD_2^2$ does not follow a (scaled) $\chi^2$ distribution, since the pre- and post-test scores are not independent. If you look at the answer given here (which is about the independent samples case), you will see that the derivation of the variance involves the non-central t-distribution, which is also involved in deriving the two variance equations above. However, to get that non-central t, you need (the square-root of) a (scaled) $\chi^2$ distribution in the denominator.

Answer 2

Structural equation modeling (SEM) is your friend in this case. We need the sample means and the covariance matrix as inputs. We may label the parameters (m1, m2, sd1, and sd2) and define $d_{av}$=(m2-m1)/((sd1+sd2)/2). SEM automatically calculates $Var(d_{av})$ by taken the correlation between the pre- and post-test scores into account.

The followings are the R code using the lavaan package. You may refer to Chapter 3 in Cheung (2015) for the details.

library("lavaan")

## Sample covariance matrix on pre- and post-test scores
lower <- '10
          8 12'
( Cov <- getCov(lower, diag=TRUE, names=c("x_pre","x_post")) )

## Means of the pre- and post-tests
Mean <- c(10, 13)

## Sample size
N <- 50

model1 <- '# Label the sds with sd1 and sd2
           eta_pre =~ sd1*x_pre
           eta_post =~ sd2*x_post
           # Fix the error variances at 0
           x_pre ~~ 0*x_pre
           x_post ~~ 0*x_post
           # Label the means with m1 and m2 
           x_pre ~ m1*1
           x_post ~ m2*1
           # Define d_av
           d_av := (m2-m1)/((sd1+sd2)/2)'

## Fit the model
fit1 <- cfa(model1, sample.cov=Cov, sample.mean=Mean, 
            sample.nobs=N, std.lv=TRUE, 
            sample.cov.rescale=FALSE)

## Display d_av and its SE
parameterEstimates(fit1)[12, -c(1,2,3)]

## Output
   label   est    se   z pvalue ci.lower ci.upper
12  d_av 0.905 0.131 6.9      0    0.648    1.163

Reference

Cheung, M. W.-L. (2015). Meta-analysis: A structural equation modeling approach. Chichester, West Sussex: John Wiley & Sons, Inc.