2 min read

Compute effect sizes with R. A primer.

2017/11/21

A typical “cook book recipe” for doing data analysis is an applied stats course is:

  1. report descriptive statistics
  2. plot some nice diagrams
  3. test hypothesis
  4. report effect sizes

Let’s have a quick glance at these steps. We will use the dataset flights of the package nycflights13.

data(flights, package = "nycflights13")

This post will be tidyverse-driven.

library(tidyverse)
library(skimr)
library(mosaic)

Let’s compute some summaries:

flights %>% 
  select(arr_delay) %>% 
  skim
#> Skim summary statistics
#>  n obs: 336776 
#>  n variables: 1 
#> 
#> Variable type: numeric 
#>   variable missing complete      n mean    sd  p0 p25 p50 p75 p100
#>  arr_delay    9430   327346 336776  6.9 44.63 -86 -17  -5  14 1272
#>      hist
#>  ▇▁▁▁▁▁▁▁

Alternatively, using mosaic:

mosaic::favstats(~arr_delay, data = flights)
#>  min  Q1 median Q3  max     mean       sd      n missing
#>  -86 -17     -5 14 1272 6.895377 44.63329 327346    9430

Subgroup statistics

Differentiating between origin levels:

flights %>% 
  select(arr_delay, origin) %>%
  group_by(origin) %>% 
  skim
#> Skim summary statistics
#>  n obs: 336776 
#>  n variables: 2 
#>  group variables: origin 
#> 
#> Variable type: numeric 
#>  origin  variable missing complete      n mean    sd  p0 p25 p50 p75 p100
#>     EWR arr_delay    3708   117127 120835 9.11 45.53 -86 -16  -4  16 1109
#>     JFK arr_delay    2200   109079 111279 5.55 44.28 -79 -18  -6  13 1272
#>     LGA arr_delay    3522   101140 104662 5.78 43.86 -68 -17  -5  12  915
#>      hist
#>  ▇▁▁▁▁▁▁▁
#>  ▇▁▁▁▁▁▁▁
#>  ▇▁▁▁▁▁▁▁

Alternatively, using mosaic:

favstats(arr_delay~origin, data = flights)
#>   origin min  Q1 median Q3  max     mean       sd      n missing
#> 1    EWR -86 -16     -4 16 1109 9.107055 45.52918 117127    3708
#> 2    JFK -79 -18     -6 13 1272 5.551481 44.27745 109079    2200
#> 3    LGA -68 -17     -5 12  915 5.783488 43.86227 101140    3522

Effect sizes

Cohen’s d

library(effsize)

We need two groups not three:

flights2 <-
  filter(flights, origin != "JFK") %>%
  sample_n(1000) %>% 
  na.omit
cohen.d(d = flights2$arr_delay,
        f = flights2$origin)
#> 
#> Cohen's d
#> 
#> d estimate: -0.01059389 (negligible)
#> 95 percent confidence interval:
#>        inf        sup 
#> -0.1372892  0.1161014

Plot mean difference

ggplot(flights2) +
  aes(x = origin, y = arr_delay) +
  geom_point(color = "grey80", position = "jitter") +
  stat_summary(fun.y = mean, geom = "point", color = "red", size = 5)

Other effect sizes

Other effect sizes can quite conveniently be derived from the package compute.es.