Que fração das experiências repetidas terá um tamanho de efeito dentro do intervalo de confiança de 95% da primeira experiência?

12

Vamos nos ater a uma situação ideal com amostragem aleatória, populações gaussianas, variações iguais, sem hackers P etc.

Etapa 1. Você realiza uma experiência, digamos, comparando duas médias de amostra e calcula um intervalo de confiança de 95% para a diferença entre as duas médias da população.

Etapa 2. Você executa muito mais experimentos (milhares). A diferença entre médias varia de experimento para experimento devido a amostragem aleatória.

Pergunta: Que fração da diferença entre médias da coleta de experimentos na etapa 2 estará dentro do intervalo de confiança da etapa 1?

Isso não pode ser respondido. Tudo depende do que aconteceu na etapa 1. Se o experimento da etapa 1 foi muito atípico, a resposta para a pergunta pode ser muito baixa.

Então, imagine que os dois passos sejam repetidos várias vezes (com o passo 2 repetido muito mais vezes). Agora, acho que seria possível chegar a uma expectativa de que fração de experimentos repetidos, em média, tenha um tamanho de efeito dentro do intervalo de confiança de 95% do primeiro experimento.

Parece que a resposta a essas perguntas precisa ser entendida para avaliar a reprodutibilidade dos estudos, uma área muito quente atualmente.

Harvey Motulsky
fonte
Para cada original (passo 1) experimentar , definir x i como a fracção de subsequente (passo 2) resulta que as conclusões produto dentro de intervalo de confiança do resultado original. Deseja calcular a distribuição empírica de x ? ixix
Matthew Gunn
Sim, você entende o que eu estou pedindo
Harvey Motulsky
@MatthewGunn perguntou se você queria a distribuição empírica da "fração de captura" para observações futuras. Seu post perguntou "... acho que seria possível ter uma expectativa de que fração de experimentos repetidos, em média, tenha um tamanho de efeito dentro do intervalo de confiança de 95% do primeiro experimento" . Esta não é uma distribuição, mas um valor esperado (médio).
A análise de Whuber é ótima, mas se você precisar de uma citação, aqui está um artigo que discute exatamente essa questão em detalhes: Cumming & Maillardet, 2006, Intervalos de confiança e replicação: onde cairá a próxima média? . Eles chamam de porcentagem de captura de um intervalo de confiança.
Ameba diz Reinstate Monica

Respostas:

12

Análise

Por se tratar de uma questão conceitual, por simplicidade, vamos considerar a situação em que um intervalo de confiança [ ˉ x ( 1 ) + Z α / 2 s ( 1 ) / 1αé construído para uma médiaμusando uma amostra aleatóriax(1)do tamanhone uma segunda amostra aleatóriax(2)é coletada do tamanhom, todas da mesmadistribuiçãoNormal(μ,σ2). (Se desejar, você pode substituirZspor valores dadistribuiçãot deStudentden-1graus de liberdade; a análise a seguir não será alterada.)

[x¯(1)+Zα/2s(1)/n,x¯(1)+Z1α/2s(1)/n]
μx(1)nx(2)m(μ,σ2)Ztn1

A chance de a média da segunda amostra estar dentro do IC determinado pela primeira é

Pr(x¯(1)+Zα/2ns(1)x¯(2)x¯(1)+Z1α/2ns(1))=Pr(Zα/2ns(1)x¯(2)x¯(1)Z1α/2ns(1)).

Como a média da primeira amostra é independente do desvio padrão da primeira amostra s ( 1 ) (isso requer normalidade) e a segunda amostra é independente da primeira, a diferença na amostra significa U = ˉ x ( 2 ) - ˉ x ( 1 ) é independente do s ( 1 ) . Além disso, para este intervalo simétrico Z α / 2 = - Z 1 - α / 2x¯(1)s(1)U=x¯(2)x¯(1)s(1)Zα/2=Z1α/2. Portanto, ao escrever para a variável aleatória s ( 1 ) e ao quadrado das duas desigualdades, a probabilidade em questão é a mesma queSs(1)

Pr(U2(Z1α/2n)2S2)=Pr(U2S2(Z1α/2n)2).

As leis da expectativa implicam que tem uma média de 0 e uma variação deU0

Var(U)=Var(x¯(2)x¯(1))=σ2(1m+1n).

Como é uma combinação linear de variáveis ​​normais, também possui uma distribuição normal. Portanto U 2 é σ 2 ( 1UU2σ2(1n+1m) times a χ2(1) variable. We already knew that S2 is σ2/n times a χ2(n1) variable. Consequently, U2/S2 is 1/n+1/m times a variable with an F(1,n1) distribution. The required probability is given by the F distribution as

(1)F1,n1(Z1α/221+n/m).

Discussion

An interesting case is when the second sample is the same size as the first, so that n/m=1 and only n and α determine the probability. Here are the values of (1) plotted against α for n=2,5,20,50.

Figure

The graphs rise to a limiting value at each α as n increases. The traditional test size α=0.05 is marked by a vertical gray line. For largish values of n=m, the limiting chance for α=0.05 is around 85%.

By understanding this limit, we will peer past the details of small sample sizes and better understand the crux of the matter. As n=m grows large, the F distribution approaches a χ2(1) distribution. In terms of the standard Normal distribution Φ, the probability (1) then approximates

Φ(Z1α/22)Φ(Zα/22)=12Φ(Zα/22).

For instance, with α=0.05, Zα/2/21.96/1.411.386 and Φ(1.386)0.083. Consequently the limiting value attained by the curves at α=0.05 as n increases will be 12(0.083)=10.166=0.834. You can see it has almost been reached for n=50 (where the chance is 0.8383.)

For small α, the relationship between α and the complementary probability--the risk that the CI does not cover the second mean--is almost perfectly a power law. Another way to express this is that the log complementary probability is almost a linear function of logα. The limiting relationship is approximately

log(2Φ(Zα/22))1.79712+0.557203log(20α)+0.00657704(log(20α))2+

In other words, for large n=m and α anywhere near the traditional value of 0.05, (1) will be close to

10.166(20α)0.557.

(This reminds me very much of the analysis of overlapping confidence intervals I posted at /stats//a/18259/919. Indeed, the magic power there, 1.91, is very nearly the reciprocal of the magic power here, 0.557. At this point you should be able to re-interpret that analysis in terms of reproducibility of experiments.)


Experimental results

These results are confirmed with a straightforwward simulation. The following R code returns the frequency of coverage, the chance as computed with (1), and a Z-score to assess how much they differ. The Z-scores are typically less than 2 in size, regardless of n,m,μ,σ,α (or even whether a Z or t CI is computed), indicating the correctness of formula (1).

n <- 3      # First sample size
m <- 2      # Second sample size
sigma <- 2 
mu <- -4
alpha <- 0.05
n.sim <- 1e4
#
# Compute the multiplier.
#
Z <- qnorm(alpha/2)
#Z <- qt(alpha/2, df=n-1) # Use this for a Student t C.I. instead.
#
# Draw the first sample and compute the CI as [l.1, u.1].
#
x.1 <- matrix(rnorm(n*n.sim, mu, sigma), nrow=n)
x.1.bar <- colMeans(x.1)
s.1 <- apply(x.1, 2, sd)
l.1 <- x.1.bar + Z * s.1 / sqrt(n)
u.1 <- x.1.bar - Z * s.1 / sqrt(n)
#
# Draw the second sample and compute the mean as x.2.
#
x.2 <- colMeans(matrix(rnorm(m*n.sim, mu, sigma), nrow=m))
#
# Compare the second sample means to the CIs.
#
covers <- l.1 <= x.2 & x.2 <= u.1
#
# Compute the theoretical chance and compare it to the simulated frequency.
#
f <- pf(Z^2 / ((n * (1/n + 1/m))), 1, n-1)
m.covers <- mean(covers)
(c(Simulated=m.covers, Theoretical=f, Z=(m.covers - f)/sd(covers) * sqrt(length(covers))))
whuber
fonte
You say that using t instead of z won't make much difference. I believe you but haven't checked yet. With small sample size, the two critical values can differ a lot and the t distribution is the correct way to compute the CI. Why do you prefer to use z??
Harvey Motulsky
It's purely illustrative and Z is simpler. When you use t it is interesting that the curves in the figure start high and descend to their limit. In particular, the chance of reproducing a significant result is then much higher for small samples than for large! Note that there's nothing to check, because you are free to interpret Zα as a percentage point of the appropriate Student t distribution (or of any other distribution you might care to name). Nothing changes in the analysis. If you do want to see the particular effects, uncomment the qt line in the code.
whuber
1
+1. This is a great analysis (and your answer has way too few upvotes for what it is). I just came across a paper that discusses this very question in great detail and I thought you might be interested: Cumming & Maillardet, 2006, Confidence Intervals and Replication: Where Will the Next Mean Fall?. They call it capture percentage of a confidence interval.
amoeba says Reinstate Monica
@Amoeba Thank you for the reference. I especially appreciate one general conclusion therein: "Replication is central to the scientific method, and researchers should not turn a blind eye to it just because it makes salient the inherent uncertainty of a single study."
whuber
1
Update: Thanks to the ongoing discussion in the sister thread, I now believe my reasoning in the above comment was not correct. 95% CIs have 83% "replication-capture", but this is a statement about repeated sampling and cannot be interpreted as giving a probability conditioned on one particular confidence interval, at least not without further assumptions. (Perhaps both this and previous comments should better be deleted in order not to confuse further readers.)
amoeba says Reinstate Monica
4

[Editado para corrigir o bug apontado pelo WHuber.]

Alterei o código R do @ Whuber para usar a distribuição t e plotar a cobertura como uma função do tamanho da amostra. Os resultados estão abaixo. No tamanho da amostra alto, os resultados correspondem ao WHuber's, é claro.

enter image description here

And here is the adapted R code, run twice with alpha set to either 0.01 or 0.05.

sigma <- 2 
mu <- -4
alpha <- 0.01
n.sim <- 1e5
#
# Compute the multiplier.

for (n in c(3,5,7,10,15,20,30,50,100,250,500,1000))
{
   T <- qt(alpha/2, df=n-1)     
# Draw the first sample and compute the CI as [l.1, u.1].
#
x.1 <- matrix(rnorm(n*n.sim, mu, sigma), nrow=n)
x.1.bar <- colMeans(x.1)
s.1 <- apply(x.1, 2, sd)
l.1 <- x.1.bar + T * s.1 / sqrt(n)
u.1 <- x.1.bar - T * s.1 / sqrt(n)
#
# Draw the second sample and compute the mean as x.2.
#
x.2 <- colMeans(matrix(rnorm(n*n.sim, mu, sigma), nrow=n))
#
# Compare the second sample means to the CIs.
#
covers <- l.1 <= x.2 & x.2 <= u.1
#
Coverage=mean(covers)

print (Coverage)

}

And here is the GraphPad Prism file that made the graph.

Harvey Motulsky
fonte
Eu acredito que seus gráficos não usam a distribuição t , devido a um erro: você define o valor de Tfora do loop! Se você deseja ver as curvas corretas, plote-as diretamente usando o resultado teórico em minha resposta, conforme indicado no final do meu Rcódigo (em vez de confiar nos resultados simulados):curve(pf(qt(.975, x-1)^2 / ((x * (1/x + 1/x))), 1, x-1), 2, 1000, log="x", ylim=c(.8,1), col="Blue"); curve(pf(qt(.995, x-1)^2 / ((x * (1/x + 1/x))), 1, x-1), add=TRUE, col="Red")
whuber
1
@whuber. Caramba! Claro que voce esta certo. Embaraçoso. Eu consertei isso. Como você apontou, a cobertura é maior com amostras pequenas. (Corrigi as simulações e não tentei sua função teórica.)
Harvey Motulsky
I am glad you fixed it, because it is very interesting how high the coverage is for small sample sizes. We could also invert your question and use the formula to determine what value of Zα/2 to use if we wished to assure (before doing any experiments), with probability p=0.95 (say), that the mean of the second experiment would lie within the two-sided 1α confidence interval determined from the second. Doing so, as a routine practice, could be one intriguing way of addressing some criticism of NHST.
whuber
@whuber I think the next step is to look at the distribution of coverage. So far, we have the average coverage (average of many first experiments, with average of many second experiments each). But depending on what the first experiment is, in some cases the average coverage will be poor. It would be interesting to see the distribution. I'm trying to learn R well enough to find out.
Harvey Motulsky
Regarding the distributions, see the paper I linked to in the comments above.
amoeba says Reinstate Monica