# How to get proper effect size from LMM?

I now have a problem with effect size of LMM. Someone insisted that I should have effect sizes after p values. I then thought that I can use 'estimate' in the output of pairwise comparison using lsmeans function. He is not convinced by using estimate and suggested me to use cohen's d. Then I look at the calculation of cohen'd; it only uses means and SD of groups; it does not come from variance in LMM.

Another thing is also about effect size, I have used mixed function in afex in R which gives me overall significant values of parameters, but does not provide me effect size. Then I emailed asking the author of that function how I can have effect size, and he replied that there is no agreement on how to calculate effect size in LMM, so he didn't provide in that function, but again my colleague insisted that I should have it. I feel like if I have to us effect size calculated from ANOVA table, then it is not effect size of LMM that I calculate. This is output from LMM in afex.

        > mixed3 <- mixed(peak_Mid ~ (1|item) + (1+vowel3|speaker) + sex*vowel3*Language, data=data1.frame, na.action=na.omit)
Fitting 9 (g)lmer() models:
[.........]
Obtaining 8 p-values:
[Note: method with signature ‘sparseMatrix#ANY’ chosen for function ‘kronecker’,
target signature ‘dgCMatrix#ngCMatrix’.
"ANY#sparseMatrix" would also be valid
........]
> summary(mixed3)
Effect        stat ndf       ddf F.scaling      p.value      stat.U ndf.U     ddf.U F.scaling.U    p.value.U
1         (Intercept) 9500.922104   1  70.40672 1.0000000 7.529698e-77 9500.922104     1  70.40672          NA 7.529698e-77
2                 sex   15.980281   1  71.52842 1.0000000 1.538529e-04   15.980281     1  71.52842          NA 1.538529e-04
3              vowel3    8.596702   2  27.40531 0.9916348 1.264905e-03    8.669222     2  27.40531          NA 1.209863e-03

4            Language    3.996819   2  70.74337 0.9909675 2.267036e-02    4.033250     2  70.74337          NA 2.194066e-02
5          sex:vowel3    1.746398   2  75.92257 0.9870432 1.813334e-01    1.769323     2  75.92257          NA 1.774036e-01
6        sex:Language    4.136050   2 170.78334 0.9964821 1.761500e-02    4.150652     2 170.78334          NA 1.737140e-02
7     vowel3:Language    1.573332   4  66.15951 0.9799283 1.917146e-01    1.605559     4  66.15951          NA 1.832701e-01
8 sex:vowel3:Language    1.239002   4 195.29430 0.9894859 2.956981e-01    1.252168     4 195.29430          NA 2.903144e-01
>


Do you have any suggestions on if LMM has any way to calculate effect size rather than using estimate?

## 2 Answers

As I understand it, Cohen's $d = (\mu_i-\mu_j)/\sigma_y$ where $\sigma_y$ is the SD of the response variable. You have estimates of the numerator from lsmeans, so what you need is an estimate of $\sigma_y$.

1. Idea 1: Your mixed model appears to have random effects for item and speaker, and the coefficient of vowel3. So the random terms in the model say that $$\sigma_y^2 = \sigma^2_I + \sigma^2_S + v_3^2\sigma^2_V + \sigma^2_E$$ Wher the $\sigma^2$s here denote the variances of items, speakers, vowel3 slopes, and error, and $v_3$ is a particular value of vowel3. If you can use the estimates of these variances from the model-fitting output and add them up and take a square root. You need to decide on what $v_3$ to use -- maybe the mean of vowel3?
2. Idea 2: Compute err <- y - predict(model), where the predictions are based only on the fixed effects. Then use sd(err) as an estimate of $\sigma_y$. This implicitly weights the estimate according to the distribution of the predictor combinations, but maybe this would be considered appropriate to do.

Disclaimer: Personally, I don't see what's so important or relevant about effect sizes in most situations, especially where it is difficult to decide what you even mean by them. So my answer is admittedly focused on finding a way to make the reviewer happy.

• Hi Russ, thank you for your suggestions. I have tried the second way (as the first way is too complicated for me) using 'err <- data1.frame$sex - predict(mixed3)', but it shows 'Error in UseMethod("predict") : no applicable method for 'predict' applied to an object of class "mixed"'. Is it because the second formula cannot apply with mixed function in afex? – user3288202 Sep 7 '14 at 8:15 • I haven't used the afex package, but evidently it doesn't implement a predict method for its objects. You may need to fit the same model using lmer in the lme4 package, and then use its predict method. I think older versions of lme4 didn't have the predict method either, and maybe afex is based on one of those. – Russ Lenth Sep 7 '14 at 11:33 • Hi Russ, it does not work in lme4 either: > err <- data1.frame$peak_Mid - predict(lmer50) Warning message: In data1.frame\$peak_Mid - predict(lmer50) : longer object length is not a multiple of shorter object length – user3288202 Sep 8 '14 at 7:32
• I can't make much sense of the above. But I think you need to make sure you're asking for the right predictions. And I suggest it might help to unravel this if you get the predictions in one variable and then see how many values you got, and then you can compare that with the number of y values. Perhaps some cases got thrown out due to incomplete observations, in which case you'll have fewer predictions than observations. The complete.cases function can be used to identify the ones that were actually used. – Russ Lenth Sep 8 '14 at 17:59

You should use Cohen's f square as effect size. Refer to this paper: Selya, Arielle S., et al. "A practical guide to calculating Cohen’s f2, a measure of local effect size, from PROC MIXED." Frontiers in psychology 3 (2012): 111.