mtorchiano / effsize

Effsize - a package for efficient effect size computation
GNU General Public License v2.0
104 stars 11 forks source link

Confidence interval computation correctness #28

Closed mtorchiano closed 6 years ago

mtorchiano commented 6 years ago

I try to calculate the CI for Cohen’s d. I used two packages for this but the results don’t match.

Please consider this example:

data(sleep)
(tTestResult <- t.test(extra ~ group, sleep, paired = TRUE))

if (!require(psych)) install.packages("psych")
psych::cohen.d.ci(apa::cohens_d(tTestResult),
                  n1 = tTestResult$parameter + 1,
                  alpha = .01)

with psych I get the result

         lower    effect     upper
[1,] -2.655439 -1.284558 0.1132465

If I use your function I get

if (!require(effsize)) install.packages("effsize")
effsize::cohen.d(extra ~ group,
                 data = sleep,
                 paired = T,
                 conf.level = 0.99
)
d estimate: -1.284558 (large)
      inf        sup
-2.6983736  0.1292585

However psych calculates with noncentral distribution. Therefore I used it in your function as well

effsize::cohen.d(extra ~ group,
                data = sleep,
                paired = T,
                conf.level = 0.99,
                noncentral = T)

However, now the results are totally messed up. The CI are outside from the effect.

d estimate: -1.284558 (large)
     inf       sup
0.1729418 2.5309425

Can you tell me: Do I use your function wrongly?

mtorchiano commented 6 years ago

Concerning the estimate being outside the CI limits: a bug was found where the sign of the CI limits was flipped.

Concerning the correctness of the values, let us refer to an example available in the literature.

I consider the paper by David C. Howell available online.

The paper reports a case -- taken from Adams, Wright, and Lohr (1996) -- with two samples

Homophobic Nonhomophobic
Mean 24.00 16.5
Variance 148.87 139.16
N 35 29

With the t statistics: t=2.48, df =62 , d=0.62

The Cohen's d CI is then [0.117, 1.124].

We can generate the data as

  set.seed(537)
  hf = generate_data(35,24,sqrt(148.87))
  nhf = generate_data(29,16.5,sqrt(139.16))

The fixed version of effsize package can be called as

cohen.d(hf,nhf,noncentral=TRUE)

and it returns

Cohen's d

d estimate: 0.6239505 (medium)
95 percent confidence interval:
     inf      sup 
0.117339 1.125768 

Which is perfectly in line with the values reported in the paper (the paper used truncated values therefore a slight difference is observable for the upper limit)

Using the alternate version from psych package we have

ttr = apa::t_test(hf,nhf)
psych::cohen.d.ci(apa::cohens_d(ttr),
                   n1 = ttr$parameter + 1,
                   alpha = .05)

that returns

         lower    effect    upper
[1,] 0.3483285 0.6239505 0.895237

Which is clearly different from the expected result reported in the paper.


The function to generate a data set with given mean and stdev is:

generate_data <- function(n,m,stdev){
  x <- rnorm(n,m,stdev)
  sd.adj = stdev/sd(x)
  x <- x * sd.adj
  m.adj = m - mean(x)
  x <- x + m.adj
  return(x)
}
dfeistauer commented 6 years ago

Hi Marco,

I installed now devtools::install_github("mtorchiano/effsize") and have now version 0.7.2.

However, the result is still the same as reported first. Your example is not for paired data.

Thanks

Daniela

mtorchiano commented 6 years ago

I take another example from the same paper.

The paper reports a single sample case with an expected mean = 1, I turned into the equivalente paired two sample with one sample being all 1s:

moon.data = c(1.73, 1.06, 2.03, 1.40, 0.95, 1.13, 1.41, 1.73, 1.63, 1.56)
g1 = rep(1,length(moon.data))

The paper reports a 95% CI equal to 0.488 ≤ d ≤ 2.33

Using the psych::cohen.d.ci() function:

ttr = t.test(moon.data,g1,paired=TRUE)
psych::cohen.d.ci(apa::cohens_d(ttr),
                  n1 = ttr$parameter + 1,
                  alpha = .05)

I get:

         lower   effect    upper
[1,] 0.4655942 1.359018 2.216011

Using the effsize::cohen.d() function:

effsize::cohen.d(moon.data,g1,
                       paired = TRUE,
                       noncentral = TRUE,
                       conf.level = 0.95) 

I get

Cohen's d

d estimate: 1.359018 (large)
95 percent confidence interval:
    lower     upper 
0.4907785 2.3358769 

This latter result is clearly closer to the one in the paper, though the values in the paper have been computed using truncated results, this explains the slight difference.

dfeistauer commented 6 years ago

Thank you Marco. Now it seems to work.