Tamanho da amostra e cáluclo do ‘poder’ do teste

Quantas observações eu preciso para estimar com precisão os parâmetros de meu modelo? Conforme Gelman e Hill (pg. 438) o tamanho da amostra nunca é grande o suficiente, pois, conforme o \(n\) cresce, estimamos mais interações ou, se os dados tiverem estruturados em grupos, mais subgrupos de nossos dados. O tamanho da amostra necessário depende do efeito que queremos estimar, o seu erro padrão e o nível de confiança desejado. Vamos supor que queremos estimar a proporção de pessoas que apoiam a descriminalização do aborto. Suponha também que suspeitamos que esse apoio é de 60%. Queremos descobrir a verdadeira proporção \(p\) com uma precisão (erro padrão) de não mais que 5% com uma amostra de tamanho \(n\). O erro padrão da média é dado por \(\sqrt{p(1-p)/n}\). Substituindo \(p\) por nosso palpite de 60% temos \(\sqrt{0.60 \times 0.40)/n} = 0.49/\sqrt{n}\). Queremos que essa quantia seja menor que 0.05, daí que \(0.49/\sqrt{n} \leq 0.05\) ou \(n \geq 96\).

Em geral não temos ideia de qual seja a proporção real, neste caso usamos a estimativa conservadora de 50% como ponto de comparação e temos \(\sqrt{(0.50 \times 0.50)/n} = 0.50/\sqrt{n} \leq 0.05\) e \(n \geq 100\). Suponha agora que queremos demonstrar que mais da metade das pessoas apoiam a descriminalização do aborto, isto é, que \(\hat{p} > 0.5\). Usando a estimativa conservadora para \(\hat{p} = 0.50/\sqrt{n}\), num intervalo de confiança de 95% em torno do verdadeiro \(p\) temos \([\hat{p} = \pm 1.96 \times 0.5/\sqrt{n}]\), assim teremos mais segurança em afirmar que \(p > 0.5\) se o intervalo estiver acima de \(\hat{p} > 0.5 + 1.96 \times 0.5/\sqrt{n}\), ou se a estimativa estiver, pelo menos, a 1.96 erros padrão do ponto de comparação de 0.5.

A fórmula utilizada para calcular o tamanho da amostra para proporções é dada por

\[n = \frac{Z^{2}\hat{p}(1-\hat{p})}{s^2}\]

Onde Z é o Z-score - o desvio com relação à média numa distribuição normal padronizada - que pode adotar valores como 1.96 para um intrevalo de 95%, 1.28 para 80% e 0.68 para 50% e \(s^2\) é a variância1

Uma ideia equivocada é que se nosso palpite sobre o valor de \(p\) estiver correto ele estará a 1.96 erros padrão de \(\hat{p}\) isto é \(0.60 > 0.5 + 1.96 \times 0.5/\sqrt{n}\) e bastaria \(n > (1.96 \times 0.5/0.1)^2 = 96\) observações para obter o erro desejado, conforme a figura abaixo.

plot(seq(0.2,0.7994,0.006),dnorm(seq(0.2,0.7994,0.006), 0.5,0.05), type = "n", xlab = "", ylab = "n")
lines(seq(0.2,0.7994,0.006),dnorm(seq(0.2,0.7994,0.006), 0.5,0.05))
abline(v = 0.5)
abline(v = 0.6, lty = 2)

No entanto se a proporção real \(p\) for 0.6, conforme nosso palpite, então \(\hat{p}\) depende da amostra e tem distribuição normal com média 0.6 e desvio padrão \(\sqrt{0.60 \times 0.40)/n} = 0.49/\sqrt{n}\). Com \(n = 96\) cerca da metade dos intervalos de confiança de 95% incluem o ponto de comparação 0.5, não sendo estatisiticamente diferentes dele, como mostra a figura abaixo

nsims <- 100
confins <- matrix(NA, nsims,2)
set.seed(2345)
for(i in 1:nsims){
  m.tmp <- rnorm(1, 0.6, 0.49/sqrt(96))
  s.tmp <- m.tmp + 1.96*(0.49/sqrt(96))
  i.tmp <- m.tmp - 1.96*(0.49/sqrt(96))
  confins[i,] <- c(i.tmp,s.tmp)
}
plot(seq(0.4,0.7996,0.004),1:100, type = "n", xaxt = "n", yaxt = "n", bty = "n", xlab = "", ylab = "")
axis(1, c(0.4,0.5,0.6,0.7,0.8), labels = c(0.4,0.5,0.6,0.7,0.8))

for(i in 1:100){
  lines(confins[i,], rep(i,2), col = "darkgray", lwd = 0.7)
}
abline(v = 0.5, lty = 2)

sum(confins[,1] > 0.5)/100
## [1] 0.51

O valor de 0.51 nos dá a probabilidade do intervalo de confiança de 95% estar acima do ponto de comparação de 0.5. Essa probabilidade é chamada de poder. O tamanho de amostra apropriada depende também do poder desejado. Por convenção o poder desejado de um teste é de 0.8.2 Quanto maior o \(n\) mais os valores estimados se aproximam do valor verdadeiro e a variabilidade diminui, diminuindo o intervalo de onfiança. No nosso exemplo, para garantir que 80% dos intervalos de confiança estejam acima do ponto de comparação de 0.5 precisamos que

\[ 0.5 + 1.96\text{e.p.} = 0.6 - 0.84\text{e.p.} \\ (1.96+0.84)\text{e.p.} = 0.1 \]

Substituindo \(\text{e.p.}\) por \(0.5/\sqrt{n}\) encontramos o tamanho da amostra necessário para detectar uma diferença de 0.1 com nível de significância de 0.05 e poder de 0.8

\[ 2.8\text{e.p.} = 0.1\\ 2.8 \times 0.5/\sqrt{n} = 0.1\\ 2.8 \times 0.5/0.1 = \sqrt{n}\\ 14^2 = 196 = n \]

nsims <- 100
confins <- matrix(NA, nsims,2)
set.seed(2345)
for(i in 1:nsims){
  m.tmp <- rnorm(1, 0.6, 0.49/sqrt(196))
  s.tmp <- m.tmp + 1.96*(0.49/sqrt(196))
  i.tmp <- m.tmp - 1.96*(0.49/sqrt(196))
  confins[i,] <- c(i.tmp,s.tmp)
}
plot(seq(0.4,0.7996,0.004),1:100, type = "n", xaxt = "n", yaxt = "n", bty = "n", xlab = "", ylab = "")
axis(1, c(0.4,0.5,0.6,0.7,0.8), labels = c(0.4,0.5,0.6,0.7,0.8))

for(i in 1:100){
  lines(confins[i,], rep(i,2), col = "darkgray", lwd = 0.7)
}
abline(v = 0.5, lty = 2)

sum(confins[,1] > 0.5)/100
## [1] 0.79

Quando estamos comparando duas amostras o erro padrão da diferença entre duas proporções é dado por \(\sqrt{p_1(1-p_1)/n_1 + p_2(1-p_2)/n_2}\) com um limite superior dados por \(0.5\sqrt{1/n_1 + 1/n_2}\) e, quando as amostras têm o mesmo tamanho esse limite se reduz para \(1/\sqrt{n}\). Um erro padrão desejado pode ser obtido com uma amostra de tamanho \(n=1/(\text{s.e.})^2\). Para obter um poder de 80% para distinguir duas proporções que diferem em 10% precisaríamos de um \(n = [2.8/(p_1-p_2)]^2\), sendo conservador, ou, mais precisamente, \(n = 2[p_1(1-p_1) + p_2(1-p_2)]\times [2.8/(p_1-p_2)]^2\).

# A vacina da Pfizer

# A WHO aceita vacinas com 50% de eficácia, ou seja, com metade dos infectados no tratamento.
# Uma amostra suficiente, com nivel de significância de 0.05 e poder de 0.8, seria dada por

2*(0.75*0.25)/(0.01^2)
## [1] 3750
# No estudo da Pfizer

n = 170
c = 162
t = 8

e <- 1 - t/c # eficácia de 95%

((162/170*(1-162/170)) + (8/170*(1-8/170)))/0.05^2
## [1] 35.87543

No caso de variáveis contínuas se queremos estimar uma certa média \(\theta\) das observações \(y_1,y_2,\dots,y_n\) com um determinado erro padrão esse erro é \(\sigma/\sqrt{n}\), onde \(\sigma\) é o desvio padrão de \(y\) na população, e \(n = (\sigma/\text{e.p.})^2\). Se queremos um poder de distinção de 80% então um tamanho conservador para a amostra seria de \(n = (2.8\sigma/(\theta - \theta_0))^2\)3.

Assim como no caso das proporções podemos estender estes resultados para comparação de médias, onde o erro padrão de \(y_1 - y_2\) é dado por \(\sqrt{\sigma_1^2/n_1 + \sigma_2^2/n_2}\). Com o mesmo número de casos e mesma variância nos grupos teríamos: \(\text{e.p.} = 2\sigma/\sqrt{n}\) e o tamanho da amostra seria dado por \(n = (2\sigma/\text{e.p.})^2\). Para um poder de distinção da diferença \(\Delta\) de 80% teríamos \(n = 2(\sigma_1^2 + \sigma_2^2)(2.8/\Delta)^2\) que se reduz a \((5.6\sigma/\text{e.p.})^2\) se o desvio padrão for o meso nos dois grupos.

Podemos utilizar esses cálculos para entender melhor os efeitos encontrados em nossos estudos. Suponha que você faça uma regressão de \(y\) em \(x\) e obtenha um \(\hat\beta = 0.7\%\) com um erro padrão de \(1.9\%\), portanto não significativo. Para que esta estimativa seja significativamente diferente de zero seu erro padrão deve ser de \(0.35\%\) o que implica um ajuste na amostra original de \((1.9\%/0.35\%)^2 = 29\). Suponha que nossa amostra tenha \(1000\) observações, para obterum \(\hat\beta\) com significância precisaríamos ter \(1000*29 = 29000\) observações. Para obter um poder de 80% dividiríamos \(0.7\%/2.8 = 0.25\%\) o que implicaria em uma amostra de \((1.9\%/0.35\%)^2\times 1000 = 57760\). Esses cálculos partem da premissa de que o verdadeiro \(\beta = 0.7\). Esse número pode ser menor o que exigiria uma amostra ainda maior. Isso quer dizer que não devemos sempre esperar significância estatística, mesmo quando o tamanho da amostra foi calculado corretamente.

Estimar interações é ainda mais difícil que estimar os efeitos principais. Para se estimar uma interação com o mesmo nível de precisão que um efeito principal você necessitaria de uma amostra quatro vezes maior. Se o efeito da interação for metade do efeito principal o tamanho da amostra precisaria ser 16 vezes maior4. Isso mostra que não é apropriado se espear significância estatística de interações e que não devemos considerar o efeito de interações que são teóricamente interessantes como sendo zero só porque não atingem o nível de significância.

Vamos ver isso no R (exemplo do livro GELMAN, Andrew; HILL, Jennifer; VEHTARI, Aki. Regression and other stories. Cambridge University Press, 2021)

N <- 1000
sigma <- 10
y <- rnorm(N, 0, sigma)
x1 <- sample(c(-0.5,0.5), N, replace=TRUE)
x2 <- sample(c(-0.5,0.5), N, replace=TRUE)
fake <- data.frame(c(y,x1,x2))

#' #### Fit models
fit_1a <- stan_glm(y ~ x1, data = fake, refresh = 0)
fit_1b <- stan_glm(y ~ x1 + x2 + x1:x2, data = fake, refresh = 0)
print(fit_1a)
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ x1
##  observations: 1000
##  predictors:   2
## ------
##             Median MAD_SD
## (Intercept) -0.3    0.3  
## x1           0.1    0.6  
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 9.9    0.2   
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg
print(fit_1b)
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ x1 + x2 + x1:x2
##  observations: 1000
##  predictors:   4
## ------
##             Median MAD_SD
## (Intercept) -0.3    0.3  
## x1           0.1    0.6  
## x2          -0.2    0.7  
## x1:x2       -2.2    1.2  
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 9.9    0.2   
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg

Como era de se esperar os erros padrões dos efeitos principais foram \(2\sigma/\sqrt{n} = 2*10/\sqrt{1000} = 0.6\) e o erro padrão da interação era duas vezez maior \(4\sigma/\sqrt{n} = 1.3\).

Uma maneira de se avaliar o desenho de um estudo é por meio de dados sintéticos. No cap. 16 de Gelman, Hill e Vehtari (2021) os autores usam a simulação dos efeitos de um tratamento para melhorar o desempenho de uma amostra de 100 alunos em um exame final.

n <- 100
y_if_control <- rnorm(n, 60, 20)
y_if_treated <- y_if_control + 5
z <- sample(rep(c(0,1), n/2))
y <- ifelse(z==1, y_if_treated, y_if_control)
fake <- data.frame(y, z)
diff <- mean(y[z==1]) - mean(y[z==0])
se_diff <- sqrt(sd(y[z==0])^2/sum(z==0) + sd(y[z==1])^2/sum(z==1))
print(c(diff, se_diff), digits=2)
## [1] 7.9 3.9
fit_1a <- stan_glm(y ~ z, data=fake, refresh=0)
print(fit_1a)
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ z
##  observations: 100
##  predictors:   2
## ------
##             Median MAD_SD
## (Intercept) 59.7    2.8  
## z            8.0    3.8  
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 19.7    1.4  
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg
# Devemos comparar o e.p. com o valor assumido do tratamento (5).
# Um e-.p. para atingir sgnificância seria de 2. Se obtemos em e.p. de 4
# com n = 100 o n necesário para obter um e.p. de 2 seria proporcional a
# 1/sqrt(n) = 2/(sqrt(100)*4) = 400

# Incluindo um preditor (pré-teste) para aumentar a eficiência sem aumentar o n

fake$x <- rnorm(n, 50, 20)
fit_1b <- stan_glm(y ~ z + x, data=fake, refresh=0)
print(fit_1b)
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ z + x
##  observations: 100
##  predictors:   3
## ------
##             Median MAD_SD
## (Intercept) 61.7    5.0  
## z            8.1    3.7  
## x            0.0    0.1  
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 19.8    1.4  
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg
# Simulation with realistic (não independente) pre-test
n <- 100
true_ability <- rnorm(n, 50, 16)
x <- true_ability + rnorm(n, 0, 12)
y_if_control <- true_ability + rnorm(n, 0, 12) + 10
y_if_treated <- y_if_control + 5
z <- sample(rep(c(0,1), n/2))
y <- ifelse(z==1, y_if_treated, y_if_control)
fake_2 <- data.frame(x, y, z)
fit_2a <- stan_glm(y ~ z, data=fake_2, refresh=0)
fit_2b <- stan_glm(y ~ z + x, data=fake_2, refresh=0)
print(fit_2a)
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ z
##  observations: 100
##  predictors:   2
## ------
##             Median MAD_SD
## (Intercept) 58.7    3.0  
## z            7.9    4.0  
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 20.6    1.5  
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg
print(fit_2b)
## stan_glm
##  family:       gaussian [identity]
##  formula:      y ~ z + x
##  observations: 100
##  predictors:   3
## ------
##             Median MAD_SD
## (Intercept) 27.6    4.8  
## z            7.5    3.3  
## x            0.6    0.1  
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## sigma 16.4    1.2  
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg

Gelman e Hill (2007:450) traz um exemplo de um desenho mais complexo onde os dados tem uma estrutra multinível. No exemplo o interesse está no efeito do tatamento com zinco em crianças com HIV. O HIV reduz a imunidade das pessoas ao atacar a proteina CD4 que reveste as células T. Níveis baixos de CD4 indicam o dano que o vírus causou ao sisitema imune. Acredita-se que alimentação rica em zinco ajudaria a evitar a queda da CD4 preservando o sistem aimune das crianças.

Para estimar o efeito do tratamento a ideia é medir periódicamente (K períodos) o nível de CD4 em J crianças. Dados preliminares mostraram que, em crianças sem tratamento, havia uma queda média de 0.5 a cada período. Isto é, \(\beta = -0.5\), \(\sigma_{\beta}^2 = 0.7\), \(\sigma_{y}^2 = 0.7\) e \(\sigma_{\alpha}^2 = 1.3\) estimados por \(y_j \sim N(\alpha_{j[k]} + \beta_{j[k]}t_{k}, \sigma_y)\).

Podemos modelar um efeito do tratamento a partir do seguinte modelo

\[ \left( \begin{array}{c} \alpha_j \\ \beta_j \end{array} \right) = N \left( \left( \begin{array}{c} \gamma_{0}^{\alpha} \\ \gamma_{0}^{\alpha} + \gamma_{1}^{\beta}z_j \end{array} \right), \left( \begin{array}{cc} \sigma_{\alpha}^2 & \rho\sigma_{\alpha}\sigma_{\beta} \\ \rho\sigma_{\alpha}\sigma_{\beta} & \sigma_{\beta} \end{array} \right) \right) \]

Onde \(z_j = 1\) se a criança recebeu o tratamento e 0 do contrário.

O próximo passo no exemplo de Gelman e Hill foi esepcificar o desenho do estudo. Assume-se que J crianças HIV-positivas seriam alocadas aleatoriamente a dois tratamentos com J/2 crianças recebendo os cuidados usuais e J/2 recebendo, além dos cuidados usuais, uma suplementação com zinco.Mede-se a porcentagem de CD4 a cada dois meses durante um ano (K = 7). Vamos então estabelecer o J necessário para, com um poder de 80%, obter um efeito de 0.5 (esse efeito significa zero queda de CD4 nos tratados).

No R

## Simulating the hypothetical data
CD4.fake <- function(J, K){
  time <- rep(seq(0,1,length=K), J)  # K measurements during the year
  person <- rep(1:J, each=K)         # person ID's
  treatment <- sample(rep(0:1, J/2))
  treatment1 <- treatment[person] 
#                                     # hyperparameters
  mu.a.true <- 4.8                    # more generally, these could
  g.0.true <- -.5                     # be specified as additional
  g.1.true <- .5                      # arguments to the function
  sigma.y.true <- .7
  sigma.a.true <- 1.3
  sigma.b.true <- .7
#                                     # personal-level parameters
  a.true <- rnorm(J, mu.a.true, sigma.a.true)
  b.true <- rnorm(J, g.0.true + g.1.true*treatment, sigma.b.true)
#                                     # data
  y <- rnorm(J*K, a.true[person] + b.true[person]*time, sigma.y.true)
  return (data.frame (y, time, person, treatment1))
}

## Fitting the model and checking the power
CD4.power <- function(J, K, n.sims=1000){
  signif <- rep(NA, n.sims)
  for (s in 1:n.sims){
    fake <- CD4.fake(J,K)
    lme.power <- lmer(y ~ time + time:treatment1 + (1 + time | person),
         data=fake)
    theta.hat <- fixef(lme.power)["time:treatment1"]
    theta.se <- se.fixef(lme.power)["time:treatment1"]
    signif[s] <- (theta.hat - 2*theta.se) > 0    # return TRUE or FALSE
  }
  power <- mean (signif)                         # proportion of TRUE
  return (power)
}

CD4.fake(100,7)[1:3,]
##          y      time person treatment1
## 1 2.632440 0.0000000      1          1
## 2 2.645487 0.1666667      1          1
## 3 4.064349 0.3333333      1          1
CD4.power(100,7,1000)
## boundary (singular) fit: see ?isSingular
## [1] 0.658

Para calcular analiticamente o tamanho da amostra em um desenho multinível partimos \(n = Jm\), onde temos m observações em J grupos. O erro padrão da média \(\overline{y}\) é dado por \(\sqrt{\sigma_y^2/n + \sigma_{\alpha}^2/J}\). Essa fórmula pode ser reescrita como

\[ \text{e.p. de }\overline{y} = \sqrt{\frac{\sigma_{total}^2}{Jm}[1+(m-1)\text{ICC}]}\]

Onde o ICC é a correlação intraclasse que vimos nas aulas anteriores (\(\sigma_{\alpha}^2/(\sigma_{\alpha}^2 + \sigma_{y}^2)\)). A ICC nos diz como o erro padrão é afetado pela estrutura multinível de nosso dados. Quanto mais semelhantes as observações no interior dos grupos, maior o impacto no erro padrão. Se não houver correlação intraclasse o erro padrão é simplesmente \(\sigma_{total}^2/\sqrt{n}\).

Sobre p-valor

Os exemlpos acima mostram que podemos manipular o p-valor. Aumentar a mostra, incluir preditores ou selecionar observações com o intuito de obter efeitos maiores (ex: testar vacina contra a Covid no Brasil) são práticas recomendadas para se desenhar um estudo com poder de detectar um certo efeito. Por outro lado o modo como codificamos as variáveis, a seleção de variáveis que incluímos no modelo, o tipo de modelo que usamos podem fazer um efeito significativo que não corresponde ao efeito verdadeiro.

A seguinte simulação mostra que com um número sufucente de variáveis podemos obter 1% de p-valor = 0.01, 5% de p-valor 0.05 etc., mesmo que não haja uma correlação entre as variáveis (ver https://www.r-bloggers.com/2018/04/p-values-sample-size-and-data-mining/).

set.seed(123456789)
vars = 100 # Número de variáveis explicativas
obs = 150 # número de observações

# Criando as variáveis explicativas (iid)
x = matrix(rnorm(obs*vars), obs, vars)
colnames(x) = paste0("x",1:vars)
# Uma variável resposta independente de qualquer X
y = rnorm(obs)

# Data frame
dat = cbind(data.frame(y=y),x)

# Regressão de y nas 100 variáveis explicativas
reg1 = lm(y~., data=dat)

# Coleta os resultados das regressões em um data frame
# e extrai os p-valores
library(broom)
res1 = broom::tidy(reg1)
p = res1$p.value[-1]

# Conta o número de variáveis significativas
# nos diferentes níveis de significância
sum(p <= 0.01)
## [1] 1
sum(p <= 0.05)
## [1] 5
sum(p <= 0.10)
## [1] 13

Esse resultado é esperado, pois 100 variáveis iid equivalem a 100 experimentos e um p-valor = 0.05 nos diz que 5% desses experimentos vão ser significativos a 0.05. O problema ocore quando o pesquisador seleciona variáveis com olho no p-valor, como mostra o restante do código onde um pesquisador trabalha apenas com variáveis com p < 0.25.

# Seleciona variáveis com p < 0.25

used.vars = which(res1$p.value[-1] <= 0.25)
length(used.vars) 
## [1] 37
dat2 = dat[,c(1,used.vars+1)]

# Regressão de y nas variáveis selecionadas

reg2 = lm(y~., data=dat2)

# Coleta os resultados das regressões em um data frame
# e extrai os p-valores
res2 = broom::tidy(reg2)
p = res2$p.value[-1]

# Conta o número de variáveis significativas
# nos diferentes níveis de significância
sum(p <= 0.01)
## [1] 7
sum(p <= 0.05)
## [1] 19
sum(p <= 0.10)
## [1] 25

Com a seleção das variáveis o número de p-valores com significância aumentou muito. A seguinte figura ilustra os efeitos do p-hacking.

  1. utilza-e uma correção para populações finitas divindo a equação acima por \(1 + (Z^{2}\hat{p}(1-\hat{p})/s^2N)\), onde \(N\) é o tamanho da população. Vemos que com \(N\) grande a fórmula se reduz à fórmula geral.↩︎

  2. Definimos antes o poder como 1 - \(\beta\), onde \(\beta\) é a especificiadade ou a probabilidade de falso positivo.↩︎

  3. Com amostras pequenas no lugar do Z-score de 2.8 deve-se utilizar o t de 3.1 para dar conta dos graus de liberdade perdidos↩︎

  4. Para ver como se chega a esses números ver Gelman, Hill e Vehtari (2021:301)↩︎