Warum erhöht das Sammeln von Daten bis zum Erreichen eines signifikanten Ergebnisses die Fehlerrate von Typ I?

Ich habe mich genau gefragt, warum das Sammeln von Daten bis zu einem signifikanten Ergebnis (z. B. ) die Typ-I-Fehlerrate erhöht.$p \lt .05$

Ich würde mich auch sehr über eine RDemonstration dieses Phänomens freuen.

Das Problem ist, dass Sie sich zu viele Chancen geben, um den Test zu bestehen. Es ist nur eine schicke Version dieses Dialogs:

Ich drehe dich um, um zu sehen, wer für das Abendessen bezahlt.

OK, ich rufe Köpfe.

Ratten, du hast gewonnen. Beste zwei von drei?

---

Um dies besser zu verstehen, betrachten Sie ein vereinfachtes, aber realistisches Modell dieser sequentiellen Prozedur . Angenommen, Sie beginnen mit einem "Probelauf" einer bestimmten Anzahl von Beobachtungen, sind jedoch bereit, länger zu experimentieren, um einen p-Wert von weniger als . Die Nullhypothese ist, dass jede Beobachtung (unabhängig) von einer Standardnormalverteilung stammt. Die Alternative besteht darin, dass unabhängig von einer Normalverteilung mit Einheitsvarianz mit einem Mittelwert ungleich Null vorliegt. Die Teststatistik ist der Mittelwert aller Beobachtungen, , dividiert durch ihren Standardfehler, . Bei einem zweiseitigen Test sind die kritischen Werte die$0.05$$X_i$$X_i$$n$$\bar X$$1/\sqrt{n}$$0.025$ und Prozentpunkte der Standardnormalverteilung, ungefähr.$0.975$$Z_\alpha=\pm 1.96$

Dies ist ein guter Test - für ein einzelnes Experiment mit einer festen Stichprobengröße . Es hat genau eine Wahrscheinlichkeit von , die Nullhypothese abzulehnen, egal was mag.$n$$5\%$$n$

Lassen Sie uns dies algebraisch in einen äquivalenten Test umwandeln, der auf der Summe aller Werte basiert :S n = X 1 + X 2 + ⋯ + X n = n ˉ X$n$

$$S_n=X_1+X_2+\cdots+X_n = n\bar X.$$

Somit sind die Daten "signifikant", wenn

$$\left| Z_\alpha\right| \le \left| \frac{\bar X}{1/\sqrt{n}} \right| = \left| \frac{S_n}{n/\sqrt{n}} \right| = \left| S_n \right| / \sqrt{n};$$

das ist,

$$\left| Z_\alpha\right| \sqrt{n} \le \left| S_n \right| .\tag{1}$$

Wenn wir klug sind, werden wir unsere Verluste reduzieren und aufgeben, sobald sehr groß wird und die Daten immer noch nicht in den kritischen Bereich gelangen.$n$

Dies beschreibt einen Zufallsspaziergang . Die Formel läuft darauf hinaus, einen gekrümmten parabolischen "Zaun" oder eine Barriere um das Grundstück des zufälligen Weges herum zu errichten : Das Ergebnis ist "signifikant", wenn irgendein Punkt des zufälligen Weges auf den Zaun trifft.$S_n$$(1)$$(n, S_n)$

Es ist eine Eigenschaft von zufälligen Spaziergängen, dass, wenn wir lange genug warten, es sehr wahrscheinlich ist, dass das Ergebnis irgendwann signifikant aussehen wird.

Hier sind 20 unabhängige Simulationen bis zu einer Grenze von Proben. Sie alle beginnen mit der Prüfung bei Stichproben. An diesem Punkt prüfen wir, ob jeder Punkt außerhalb der nach Formel eingezeichneten Barrieren liegt . Ab dem Punkt, an dem der statistische Test zuerst "signifikant" ist, werden die simulierten Daten rot gefärbt.$n=5000$$n=30$$(1)$

Sie können sehen, was los ist: Der zufällige Gang peitscht mehr und mehr auf und ab, wenn zunimmt. Die Barrieren breiten sich in etwa gleichem Maße aus - aber nicht schnell genug, um das zufällige Laufen zu vermeiden.$n$

In 20% dieser Simulationen wurde ein "signifikanter" Unterschied festgestellt - normalerweise ziemlich früh - obwohl in jedem von ihnen die Nullhypothese absolut korrekt ist! Wenn Sie weitere Simulationen dieses Typs ausführen, liegt die tatsächliche Testgröße eher bei als bei dem beabsichtigten Wert von Dies bedeutet, dass Sie bereit sind, bis zu einer Stichprobengröße von nach "Signifikanz" zu suchen gibt Ihnen eine Chance von , die Null zurückzuweisen, auch wenn die Null wahr ist.$25\%$$\alpha=5\%$$5000$$25\%$

Beachten Sie, dass in allen vier "signifikanten" Fällen die Daten bei fortgesetzten Tests an einigen Punkten nicht mehr signifikant aussahen. Im wirklichen Leben verliert ein Experimentator, der früh aufhört, die Chance, solche "Umkehrungen" zu beobachten. Diese Selektivität durch optionales Stoppen wirkt sich negativ auf die Ergebnisse aus.

Bei sequentiellen Tests, bei denen es sich um ehrliche Tests handelt, sind die Barrieren Linien. Sie breiten sich schneller aus als die hier gezeigten gekrümmten Barrieren.

library(data.table)
library(ggplot2)

alpha <- 0.05   # Test size
n.sim <- 20     # Number of simulated experiments
n.buffer <- 5e3 # Maximum experiment length
i.min <- 30     # Initial number of observations
#
# Generate data.
#
set.seed(17)
X <- data.table(
  n = rep(0:n.buffer, n.sim),
  Iteration = rep(1:n.sim, each=n.buffer+1),
  X = rnorm((1+n.buffer)*n.sim)
)
#
# Perform the testing.
#
Z.alpha <- -qnorm(alpha/2)
X[, Z := Z.alpha * sqrt(n)]
X[, S := c(0, cumsum(X))[-(n.buffer+1)], by=Iteration]
X[, Trigger := abs(S) >= Z & n >= i.min]
X[, Significant := cumsum(Trigger) > 0, by=Iteration]
#
# Plot the results.
#
ggplot(X, aes(n, S, group=Iteration)) +
  geom_path(aes(n,Z)) + geom_path(aes(n,-Z)) +
  geom_point(aes(color=!Significant), size=1/2) +
  facet_wrap(~ Iteration)

Menschen, die noch keine Erfahrung mit Hypothesentests haben, sind der Meinung, dass ein p-Wert unter 0,05 liegt, wenn mehr Teilnehmer hinzugefügt werden, der p-Wert nur noch weiter sinkt. Das stimmt aber nicht. Unter der Nullhypothese ist der p-Wert gleichmäßig zwischen 0 und 1 verteilt und kann in diesem Bereich ziemlich stark schwanken.

Ich habe einige Daten in R simuliert (meine R-Kenntnisse sind ziemlich einfach). In dieser Simulation sammle ich 5 Datenpunkte - jeweils mit einer zufällig ausgewählten Gruppenzugehörigkeit (0 oder 1) und jeweils mit einem zufällig ausgewählten Ergebnismaß ~ N (0,1). Ab Teilnehmer 6 führe ich bei jeder Iteration einen T-Test durch.

for (i in 6:150) {
  df[i,1] = round(runif(1))
  df[i,2] = rnorm(1)
  p = t.test(df[ , 2] ~ df[ , 1], data = df)$p.value
  df[i,3] = p
}

Die p-Werte sind in dieser Figur angegeben. Beachten Sie, dass ich signifikante Ergebnisse finde, wenn die Stichprobengröße zwischen 70 und 75 liegt. Wenn ich dort aufhöre, werde ich am Ende glauben, dass meine Ergebnisse bedeutend sind, weil ich die Tatsache übersehen habe, dass meine p-Werte mit einer größeren Stichprobe wieder angestiegen sind (dies ist mir tatsächlich einmal mit echten Daten passiert). Da ich weiß, dass beide Populationen einen Mittelwert von 0 haben, muss dies ein falsch positiver Wert sein. Dies ist das Problem beim Hinzufügen von Daten bis p <.05. Wenn Sie genügend Tests hinzufügen, wird p möglicherweise den Schwellenwert von 0,05 überschreiten, und Sie können feststellen, dass jeder Datensatz eine signifikante Auswirkung hat.

Diese Antwort betrifft nur die Wahrscheinlichkeit, dass letztendlich ein "signifikantes" Ergebnis erzielt wird, und die Verteilung der Zeit auf dieses Ereignis nach dem @ whuber-Modell.

Wie im Modell von @whuber bezeichne den Wert der Teststatistik, nachdem Beobachtungen gesammelt wurden, und nehme an, dass die Beobachtungen normal sind . Dann ist so dass sich wie eine zeitkontinuierliche Brownsche Standardbewegung verhält, wenn wir dies für den Moment ignorieren die Tatsache, dass wir einen zeitdiskreten Prozess haben (linkes Diagramm unten).$S(t)=X_1 + X_2 + \dots + X_t$$t$$X_1,X_2,\dots$

$$S(t+h)|S(t)=s_0 \sim N(s_0, h), \tag{1}$$ $S(t)$

Es sei die erste Durchgangszeit von über die zeitabhängigen Barrieren (die Anzahl der Beobachtungen, die vor dem Test erforderlich sind, wird signifikant).$T$$S(t)$$\pm z_{\alpha/2}\sqrt{t}$

Betrachten Sie den transformierten Prozess , der durch Skalieren von mit seiner Standardabweichung zum Zeitpunkt und durch Zulassen der neuen Zeitskala so dass Aus (1) und (2) folgt, dass normalerweise mit und $Y(\tau)$$S(t)$$t$$\tau=\ln t$

$$Y(\tau)=\frac{S(t(\tau))}{\sqrt{t(\tau)}}=e^{-\tau/2}S(e^\tau). \tag{2}$$ $Y(\tau+\delta)$ $$\begin{align} E(Y(\tau+\delta)|Y(\tau)=y_0) &=E(e^{-(\tau+\delta)/2}S(e^{\tau+\delta})|S(e^\tau)=y_0e^{\tau/2}) \\&=y_0e^{-\delta/2} \tag{3} \end{align}$$ $$\begin{align} \operatorname{Var}(Y(\tau+\delta)|Y(\tau)=y_0) &=\operatorname{Var}(e^{(\tau+\delta)/2}S(e^{\tau+\delta})|S(e^\tau)=y_0e^{\tau/2}) \\&=1-e^{-\delta}, \tag{4} \end{align}$$ das heißt, ist ein Ornstein-Uhlenbeck (OU) -Prozess mit einem Mittelwert von Null mit einer stationären Varianz von 1 und einer Rückkehrzeit von 2 (rechte Darstellung unten).$Y(\tau)$

Für das transformierte Modell werden die Barrieren zu zeitunabhängigen Konstanten, die . Es ist dann bekannt ( Nobile et al. 1985 ; Ricciardi & Sato, 1988 ), dass die erste Durchlaufzeit des OU-Prozesses über diese Barrieren ungefähr exponentiell mit einem Parameter (abhängig von den Barrieren bei ) (geschätzt auf für unten). Es gibt auch eine zusätzliche Punktmasse in der Größe in . "Ablehnung" von$\pm z_{\alpha/2}$$\mathcal{T}$$Y(\tau)$$\lambda$$\pm z_{\alpha/2}$$\hat\lambda=0.125$$\alpha=0.05$$\alpha$$\tau=0$$H_0$passiert schließlich mit Wahrscheinlichkeit 1. Daher folgt (die Anzahl der Beobachtungen, die gesammelt werden müssen, bevor ein "signifikantes" Ergebnis erhalten wird) ungefähr einer logarithmischen Exponentialverteilung mit dem erwarteten Wert Somit hat nur dann eine endliche Erwartung, wenn (für ausreichend große Signifikanzstufen ).$T=e^\mathcal{T}$

$$ET\approx 1+(1-\alpha)\int_0^\infty e^\tau \lambda e^{-\lambda \tau}d\tau.\tag{5}$$ $T$$\lambda>1$$\alpha$

Das Obige ignoriert die Tatsache, dass für das reale Modell diskret ist und dass der reale Prozess diskret und nicht zeitkontinuierlich ist. Daher überschätzt das obige Modell die Wahrscheinlichkeit, dass die Barriere überschritten wurde (und unterschätzt ), weil der zeitkontinuierliche Abtastpfad die Barriere möglicherweise nur vorübergehend zwischen zwei benachbarten diskreten Zeitpunkten und überschreitet . Solche Ereignisse sollten jedoch eine vernachlässigbare Wahrscheinlichkeit für große . $T$$ET$$t$$t+1$$t$

Die folgende Abbildung zeigt eine Kaplan-Meier-Schätzung von auf der Log-Log-Skala zusammen mit der Überlebenskurve für die exponentielle zeitkontinuierliche Approximation (rote Linie).$P(T>t)$

R-Code:

# Fig 1
par(mfrow=c(1,2),mar=c(4,4,.5,.5))
set.seed(16)
n <- 20
npoints <- n*100 + 1
t <- seq(1,n,len=npoints)
subset <- 1:n*100-99
deltat <- c(1,diff(t))
z <- qnorm(.975)
s <- cumsum(rnorm(npoints,sd=sqrt(deltat)))
plot(t,s,type="l",ylim=c(-1,1)*z*sqrt(n),ylab="S(t)",col="grey")
points(t[subset],s[subset],pch="+")
curve(sqrt(t)*z,xname="t",add=TRUE)
curve(-sqrt(t)*z,xname="t",add=TRUE)
tau <- log(t)
y <- s/sqrt(t)
plot(tau,y,type="l",ylim=c(-2.5,2.5),col="grey",xlab=expression(tau),ylab=expression(Y(tau)))
points(tau[subset],y[subset],pch="+")
abline(h=c(-z,z))

# Fig 2
nmax <- 1e+3
nsim <- 1e+5
alpha <- .05
t <- numeric(nsim)
n <- 1:nmax
for (i in 1:nsim) {
  s <- cumsum(rnorm(nmax))
  t[i] <- which(abs(s) > qnorm(1-alpha/2)*sqrt(n))[1]
}
delta <- ifelse(is.na(t),0,1)
t[delta==0] <- nmax + 1
library(survival)
par(mfrow=c(1,1),mar=c(4,4,.5,.5))
plot(survfit(Surv(t,delta)~1),log="xy",xlab="t",ylab="P(T>t)",conf.int=FALSE)
curve((1-alpha)*exp(-.125*(log(x))),add=TRUE,col="red",from=1,to=nmax)

Es muss gesagt werden, dass die obige Diskussion für eine frequentistische Weltanschauung ist, für die die Vielfältigkeit von der Wahrscheinlichkeit herrührt, dass Sie Daten als extremer einstufen, und nicht von der Wahrscheinlichkeit, dass Sie einen Effekt erzielen. Die Hauptursache für das Problem ist, dass p-Werte und Typ-I-Fehler eine Rückwärts-Zeit-Rückwärts-Informationsfluss-Konditionierung verwenden, was es wichtig macht, "wie Sie hierher gekommen sind" und was stattdessen passiert sein könnte . Andererseits verschlüsselt das Bayes'sche Paradigma die Skepsis gegenüber einer Auswirkung auf den Parameter selbst, nicht auf die Daten. Dadurch wird jede hintere Wahrscheinlichkeit gleich interpretiert, unabhängig davon, ob Sie vor 5 Minuten eine andere hintere Wahrscheinlichkeit für einen Effekt berechnet haben oder nicht. Weitere Details und eine einfache Simulation finden Sie unter http://www.fharrell.com/2017/10/continuous-learning-from-data-no.

Wir betrachten einen Forscher eine Stichprobe der Größe zu sammeln , , eine Hypothese zu testen . Er lehnt ab, wenn eine geeignete Teststatistik ihren Level- kritischen Wert überschreitet . Ist dies nicht der Fall, sammelt er eine weitere Stichprobe der Größe , , und lehnt ab, wenn der Test für die kombinierte Stichprobe . Erhält er immer noch keine Ablehnung, so verfährt er auf diese Weise bis zu mal insgesamt.x 1 & thgr ; = & thgr ; 0 t & agr; c n x 2 ( x 1 , x 2 ) $n$$x_1$$\theta=\theta_0$$t$$\alpha$$c$$n$$x_2$$(x_1,x_2)$$K$

Dieses Problem scheint bereits von P. Armitage, CK McPherson und BC Rowe (1969), Journal of the Royal Statistical Society, angesprochen worden zu sein . Serie A (132), 2, 235-244: "Wiederholte Signifikanztests an akkumulierten Daten" .

Der Bayes'sche Standpunkt zu diesem Thema, der auch hier erörtert wird, wird übrigens in Berger und Wolpert (1988), "The Likelihood Principle" , Abschnitt 4.2 , erörtert .

Hier ist eine unvollständige Replikation von Armitage et al die Ergebnisse (Code unten), die zeigt , wie Signifikanzniveaus aufzublasen , wenn , sowie mögliche Korrekturfaktoren zur Wiederherstellung Ebene- kritische Werte. Beachten Sie, dass die Rastersuche eine Weile dauert - die Implementierung ist möglicherweise ziemlich ineffizient.$K>1$$\alpha$

Größe der Standard-Ablehnungsregel in Abhängigkeit von der Anzahl der Versuche$K$

Größe als Funktion der Erhöhung der kritischen Werte für verschiedene$K$

Angepasste kritische Werte zur Wiederherstellung von 5% -Tests als Funktion von$K$

reps <- 50000

K <- c(1:5, seq(10,50,5), seq(60,100,10)) # the number of attempts a researcher gives herself
alpha <- 0.05
cv <- qnorm(1-alpha/2)

grid.scale.cv <- cv*seq(1,1.5,by=.01) # scaled critical values over which we check rejection rates
max.g <- length(grid.scale.cv)
results <- matrix(NA, nrow = length(K), ncol=max.g)

for (kk in 1:length(K)){
  g <- 1
  dev <- 0
  K.act <- K[kk]
  while (dev > -0.01 & g <= max.g){
    rej <- rep(NA,reps)
    for (i in 1:reps){
      k <- 1
      accept <- 1
      x <- rnorm(K.act)
      while(k <= K.act & accept==1){
        # each of our test statistics for "samples" of size n are N(0,1) under H0, so just scaling their sum by sqrt(k) gives another N(0,1) test statistic
        rej[i] <- abs(1/sqrt(k)*sum(x[1:k])) > grid.scale.cv[g] 
        accept <- accept - rej[i]
        k <- k+1
      }
    }
    rej.rate <- mean(rej)
    dev <- rej.rate-alpha
    results[kk,g] <- rej.rate
    g <- g+1
  }
}
plot(K,results[,1], type="l")
matplot(grid.scale.cv,t(results), type="l")
abline(h=0.05)

cv.a <- data.frame(K,adjusted.cv=grid.scale.cv[apply(abs(results-alpha),1,which.min)])
plot(K,cv.a$adjusted.cv, type="l")