Sobre Implicit Association Tests (IAT) – Olhe o README do repositório.

IAT: 0.15, 0.35, and 0.65 are considered small, medium, and large levels of bias for individual scores. Positive means bias towards arts / against Math.

Exemplo de análise de uma replicação

iat = read_csv(here::here(params$arquivo_dados), col_types = "cccdc")
iat = iat %>% 
    mutate(sex = factor(sex, levels = c("m", "f"), ordered = TRUE))
glimpse(iat)
## Rows: 113
## Columns: 5
## $ session_id  <chr> "2401243", "2401244", "2401246", "2401249", "2401250", "24…
## $ referrer    <chr> "brasilia", "brasilia", "brasilia", "brasilia", "brasilia"…
## $ sex         <ord> m, m, f, f, f, m, f, m, m, f, f, f, f, f, m, m, f, m, f, m…
## $ d_art       <dbl> 0.1480913, 0.6285349, 0.4977736, 0.3999447, 0.8314632, 1.1…
## $ iat_exclude <chr> "Include", "Include", "Include", "Include", "Include", "In…
iat %>%
    ggplot(aes(x = d_art, fill = sex, color = sex)) +
    geom_histogram(binwidth = .2, alpha = .4) +
    geom_rug() +
    facet_grid(sex ~ ., scales = "free_y") + 
    theme(legend.position = "None")

iat %>% 
    ggplot(aes(x = sex, y = d_art)) + 
    geom_quasirandom(width = .1)

iat %>% 
    ggplot(aes(x = sex, y = d_art)) + 
    geom_quasirandom(width = .1) + 
    stat_summary(geom = "point", fun.y = "mean", color = "red", size = 5)
## Warning: The `fun.y` argument of `stat_summary()` is deprecated as of ggplot2 3.3.0.
## ℹ Please use the `fun` argument instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Qual a diferença na amostra
resumo = iat %>%
    group_by(sex) %>%
    summarise(media = mean(d_art), dp = sd(d_art), n = n())
resumo
## # A tibble: 2 × 4
##   sex   media    dp     n
##   <ord> <dbl> <dbl> <int>
## 1 m     0.400 0.516    48
## 2 f     0.570 0.423    65
agrupado = iat %>%
        group_by(sex) %>%
        summarise(media = mean(d_art))
    m = agrupado %>% filter(sex == "m") %>% pull(media)
    f = agrupado %>% filter(sex == "f") %>% pull(media)
m - f
## [1] -0.1705546

Comparação via ICs

library(boot)
set.seed(123)

theta <- function(d, i) {
    agrupado = d %>% 
        slice(i) %>% 
        group_by(sex) %>% 
        summarise(media = mean(d_art))
    m = agrupado %>% filter(sex == "m") %>% pull(media)
    f = agrupado %>% filter(sex == "f") %>% pull(media)
    m - f
}

booted <- boot(data = iat, 
               statistic = theta, 
               R = 2000)

ci = tidy(booted, 
          conf.level = .95,
          conf.method = "bca",
          conf.int = TRUE)

glimpse(ci)
## Rows: 1
## Columns: 5
## $ statistic <dbl> -0.1705546
## $ bias      <dbl> 0.0007123825
## $ std.error <dbl> 0.08884168
## $ conf.low  <dbl> -0.3449753
## $ conf.high <dbl> 0.001041476
ci %>%
    ggplot(aes(
        x = "",
        y = statistic,
        ymin = conf.low,
        ymax = conf.high
    )) +
    geom_pointrange() +
    geom_point(size = 3) + 
    labs(x = "Diferença", 
         y = "IAT homens - mulheres")

p1 = iat %>% 
    ggplot(aes(x = sex, y = d_art)) +
    geom_quasirandom(width = .1) + 
    stat_summary(geom = "point", fun.y = "mean", color = "red", size = 5)

p2 = ci %>%
    ggplot(aes(
        x = "",
        y = statistic,
        ymin = conf.low,
        ymax = conf.high
    )) +
    geom_pointrange() +
    geom_point(size = 3) + 
    ylim(-1, 1) + 
    labs(x = "Diferença", 
         y = "IAT homens - mulheres")

grid.arrange(p1, p2, ncol = 2)

Atividade

Preencha os resultados e conclusões abaixo
media_f = resumo %>% filter(sex == "f") %>% pull(media)
dp_f = resumo %>% filter(sex == "f") %>% pull(dp)
n_f = resumo %>% filter(sex == "f") %>% pull(n)

media_m = resumo %>% filter(sex == "m") %>% pull(media)
dp_m = resumo %>% filter(sex == "m") %>% pull(dp)
n_m = resumo %>% filter(sex == "m") %>% pull(n)

diff_obs = m - f

Em média, as mulheres que participaram do experimento carregam uma associação implícita (medida pelo IAT) positiva com a matemática — ou seja, contra a matemática e a favor das artes — de magnitude média a forte (média 0.57, desv. padrão 0.423, N = 65). Os homens também tiveram uma associação positiva, mas menor que a das mulheres (média 0.4, desv. padrão 0.516, N = 48). A diferença entre os dois grupos é pequena a moderada, e ainda incerta (diferença das médias -0.171: os homens ficam cerca de 0.171 ponto abaixo das mulheres na escala IAT; 95% CI [-0.345, 0.001] pelo método BCa). Com os dados desta amostra, nossa melhor aposta é que as mulheres realmente têm uma rejeição implícita à matemática um pouco mais forte que a dos homens — mas o intervalo de confiança encosta em zero em uma das pontas, então não dá para descartar que, na população, essa diferença seja praticamente inexistente.

1. Bootstrap a partir de uma biblioteca (BCa)

Esse já foi feito lá em cima: usamos boot::boot() para gerar 2.000 réplicas bootstrap da diferença de médias, e broom::tidy(..., conf.method = "bca") para chegar no intervalo de confiança BCa (bias-corrected and accelerated), que corrige o viés e a assimetria da distribuição bootstrap.

2. Bootstraps implementados por nós: métodos percentil e normal

Agora sem depender do boot(): geramos as réplicas na mão, com replicate() e sample(..., replace = TRUE), reamostrando homens e mulheres separadamente e mantendo o tamanho de cada grupo fixo, do mesmo jeito que o exemplo da biblioteca faz por baixo dos panos.

m_scores = iat %>% filter(sex == "m") %>% pull(d_art)
f_scores = iat %>% filter(sex == "f") %>% pull(d_art)

set.seed(123)
R_manual = 2000
boot_diffs = replicate(R_manual, {
    m_i = mean(sample(m_scores, replace = TRUE))
    f_i = mean(sample(f_scores, replace = TRUE))
    m_i - f_i
})

erro_padrao_boot = sd(boot_diffs)
erro_padrao_boot
## [1] 0.0893725
media_boot <- mean(boot_diffs)
dens <- density(boot_diffs)
moda_boot <- dens$x[which.max(dens$y)]
tibble(diferenca = boot_diffs) %>%
  ggplot(aes(x = diferenca)) +
  geom_histogram( aes(y = after_stat(density)), bins = 40, fill = "steelblue", alpha = 0.6) +
  geom_density(color = "red", linewidth = 1.2 ) +
  geom_vline(xintercept = media_boot, color = "blue",linewidth = 1.0 ) +
  geom_vline(xintercept = moda_boot, color = "darkgreen", linewidth = 1.0 ) +
  annotate("text", x = media_boot,  y = Inf, label = paste("Média", round(media_boot, 3)), color = "blue", angle = 90, vjust = 1.0 ) +
  annotate("text", x = moda_boot,  y = Inf, label = paste("Moda ", round(moda_boot, 3)), color = "darkgreen", angle = 90, vjust = 3 ) +
  labs( x = "Diferença bootstrap (homens - mulheres)", y = "Densidade",    title = "Distribuição Bootstrap"
  ) +
  theme_minimal()

Olhando para o histograma, a distribuição bootstrap sai bem simétrica e parecida com uma normal — média e moda praticamente coincidem. Isso é um bom sinal de que métodos mais simples (percentil, normal) não devem se afastar muito do BCa.

Método percentil

A ideia aqui é a mais direta possível: olhamos para os 2000 valores que o bootstrap gerou e simplesmente pegamos os percentis 2.5% e 97.5% dessa distribuição como limites do intervalo.

ic_percentil = quantile(boot_diffs, probs = c(.025, .975))
ic_percentil
##          2.5%         97.5% 
## -0.3490261175 -0.0006090179
Método normal

Já este método parte do princípio de que a distribuição amostral da estatística é (aproximadamente) normal, e resume toda a variação observada no bootstrap em um único número: o erro padrão. O intervalo fica então centrado no valor observado, mais ou menos 1.96 * erro padrão bootstrap:

ic_normal = diff_obs + c(-1, 1) * qnorm(.975) * erro_padrao_boot
names(ic_normal) = c("2.5%", "97.5%")
ic_normal
##         2.5%        97.5% 
## -0.345721498  0.004612269

3. Comparação entre os três métodos

comparacao = tibble(
    metodo = c("BCa (biblioteca boot)", "Percentil (manual)", "Normal (manual)"),
    estimativa = c(ci$statistic, diff_obs, diff_obs),
    conf.low = c(ci$conf.low, ic_percentil[[1]], ic_normal[["2.5%"]]),
    conf.high = c(ci$conf.high, ic_percentil[[2]], ic_normal[["97.5%"]])
)
comparacao
## # A tibble: 3 × 4
##   metodo                estimativa conf.low conf.high
##   <chr>                      <dbl>    <dbl>     <dbl>
## 1 BCa (biblioteca boot)     -0.171   -0.345  0.00104 
## 2 Percentil (manual)        -0.171   -0.349 -0.000609
## 3 Normal (manual)           -0.171   -0.346  0.00461
comparacao %>%
    ggplot(aes(x = metodo, y = estimativa, ymin = conf.low, ymax = conf.high)) +
    geom_pointrange() +
    geom_hline(yintercept = 0, linetype = "dashed", color = "grey50") +
    labs(x = "Método", y = "IAT homens - mulheres (IC 95%)") +
    coord_flip()

Discussão. Na prática, os três métodos contam a mesma história. A estimativa pontual é sempre a mesma (afinal, ela vem direto dos dados, não do método de IC), e os intervalos de confiança ficam bem próximos entre si: todos apontam para uma diferença pequena a moderada, com a ponta superior do intervalo beirando o zero. Isso não é coincidência — como vimos no histograma, a distribuição bootstrap saiu bem simétrica (o que faz sentido, já que temos um número razoável de observações em cada grupo, 48 homens e 65 mulheres, e o Teorema Central do Limite trabalhando a nosso favor). E é justamente viés e assimetria que o BCa corrige e os outros dois métodos ignoram — se a distribuição já é simétrica, não há muito o que corrigir.

Por que escolher cada método:

  • Percentil é o caminho mais direto: pega os quantis da própria distribuição bootstrap, sem suposição nenhuma sobre sua forma. Funciona bem quando essa distribuição já é razoavelmente simétrica, como é o nosso caso, mas não corrige viés nem assimetria caso eles apareçam.
  • Normal é o mais “ingênuo” dos três: assume de cara que tudo se comporta como uma normal e resume a distribuição bootstrap inteira em um único número, o erro padrão. É o mais fácil de calcular na mão, mas também o que menos se sustenta se a distribuição real fugir muito de uma normal.
  • BCa, feito pela biblioteca boot, é o mais sofisticado dos três: corrige viés e assimetria automaticamente, e por isso costuma ser a recomendação padrão na literatura, principalmente quando a amostra é pequena ou a distribuição é enviesada.

Neste conjunto de dados específico, essa sofisticação extra do BCa não fez muita diferença na prática — os três métodos levam praticamente à mesma conclusão: a diferença entre sexos parece pequena a moderada, e a incerteza é grande o bastante para não podermos descartar uma diferença próxima de zero. Mas essa concordância não é garantida em geral: com amostras menores ou distribuições mais enviesadas, os métodos mais simples podem discordar bastante do BCa, e nesses cenários vale a pena confiar mais no método que corrige viés e assimetria.