Analyse de puissance pour les données binomiales lorsque l'hypothèse nulle est que

10

Je voudrais faire une analyse de puissance pour un seul échantillon à partir de données binomiales, avec , contre , où est la proportion de succès dans la population. Si , je pourrais utiliser soit l'approximation normale du binôme, soit -test, mais avec , les deux échouent. J'aimerais savoir s'il existe un moyen de faire cette analyse. J'apprécierais beaucoup toutes suggestions, commentaires ou références. Merci beaucoup!H 1 : p = 0,001 p 0 < p < 1 χ 2 p = 0H0:p=0H1:p=0.001p0<p<1χ2p=0

user765195
la source
Alors pourquoi n'utilisez-vous pas le test Clopper-Pearson exact?
Stéphane Laurent
2
J'espère que vous avez un très gros échantillon! Cela va être difficile à tester.
Peter Flom - Réintègre Monica

Réponses:

13

Vous avez une hypothèse alternative exacte unilatérale où et . p 1 = 0,001 p 0 = 0p1>p0p1=0.001p0=0

  • La première étape consiste à identifier un seuil pour le nombre de succès tel que la probabilité d'obtenir au moins succès dans un échantillon de taille soit très faible sous l'hypothèse nulle (classiquement ). Dans votre cas, , quel que soit votre choix particulier pour et .ccnα=0.05c=1n1α>0
  • La deuxième étape consiste à déterminer la probabilité d'obtenir au moins succès dans un échantillon de taille sous l'hypothèse alternative - c'est votre pouvoir. Ici, vous avez besoin d'un fixe tel que la distribution binomiale soit entièrement spécifiée.cnnB(n,p1)

La deuxième étape de R avec :n=500

> n  <- 500                 # sample size
> p1 <- 0.001               # success probability under alternative hypothesis
> cc <- 1                   # threshold
> sum(dbinom(cc:n, n, p1))  # power: probability for cc or more successes given p1
[1] 0.3936211

Pour avoir une idée de l'évolution de la puissance avec la taille de l'échantillon, vous pouvez dessiner une fonction de puissance: entrez la description de l'image ici

nn   <- 10:2000                 # sample sizes
pow  <- 1-pbinom(cc-1, nn, p1)  # corresponding power
tStr <- expression(paste("Power for ", X>0, " given ", p[1]==0.001))
plot(nn, pow, type="l", xaxs="i", xlab="sample size", ylab="power",
     lwd=2, col="blue", main=tStr, cex.lab=1.4, cex.main=1.4)

Si vous voulez savoir de quelle taille d'échantillon vous avez besoin pour obtenir au moins une puissance prédéfinie, vous pouvez utiliser les valeurs de puissance calculées ci-dessus. Disons que vous voulez une puissance d'au moins .0.5

> powMin <- 0.5
> idx    <- which.min(abs(pow-powMin))  # index for value closest to 0.5
> nn[idx]     # sample size for that index
[1] 693

> pow[idx]    # power for that sample size
[1] 0.5000998

Vous avez donc besoin d'un échantillon d'au moins pour atteindre une puissance de .6930.5

caracal
la source
Selon pwr.p.test, pour une puissance de 0,5, vous avez besoin d'au moins 677 observations. Mais la puissance = 0,5 est très faible!
Jessica
@caracal Utilisez-vous une approximation normale pour obtenir votre courbe de puissance? Une fonction de puissance binomiale exacte ne serait pas aussi fluide. Il est en fait en dents de scie que vous pouvez voir si l'axe de la taille de l'échantillon est agrandi. J'en discute dans mon article de 2002 dans le Statisticien américain co-écrit avec Christine Liu. De plus, le binôme est si fortement asymétrique à p très bas que n doit être grand pour que l'approximation normale fonctionne bien.
Michael R. Chernick
2
@MichaelChernick Non, cela provient des distributions binomiales, pas d'une approximation normale. Bien sûr, vous avez raison: en général, la puissance d'un test binomial est une fonction en dents de scie qui n'est pas monotone. Mais notez que nous avons un cas particulier ici avec . Cela signifie que la région d'acceptation de l'hypothèse alternative commence toujours à 1, quel que soit . Avec un seuil constant , une constante , la puissance est une fonction strictement croissante de . p0=0nc=1p1=0.001n
caracal
@Jessica Notez qu'il pwr.p.test()utilise une approximation normale, pas les distributions binomiales exactes. Tapez simplement pwr.p.testpour jeter un œil au code source. Vous trouverez les appels pour pnorm()indiquer qu'une approximation est utilisée.
caracal
1
@caracal Alors, je peux le voir de cette façon: sous l'hypothèse nulle, la probabilité de succès est 0, donc si jamais vous voyez un succès, vous pouvez rejeter l'hypothèse nulle. C'est pourquoi vous dites que le seuil est 1 car si la somme binomiale arrive à 1, vous pouvez la rejeter avec une erreur de type 2 de 0! Maintenant, selon l'alternative, la probabilité du premier succès au nième essai est (1-p) p. Cette probabilité passe à 0 lorsque n va à l'infini. Ainsi, une règle séquentielle passant par l'arrêt lorsque S = 1 aurait la puissance 1 pour tout p> 0. n1n
Michael R. Chernick
3

Vous pouvez répondre facilement à cette question avec le pwrpackage en R.

Vous devrez définir un niveau de signification, une puissance et une taille d'effet. En règle générale, le niveau de signification est défini sur 0,05 et la puissance est définie sur 0,8. Une puissance plus élevée nécessitera plus d'observations. Un niveau de signification inférieur diminuera la puissance.

La taille de l'effet pour les proportions utilisées dans ce package est le h de Cohen. Le seuil pour un petit h est souvent considéré comme 0,20. La coupure réelle varie selon l'application et peut être plus petite dans votre cas. Un h plus petit signifie que davantage d'observations seront nécessaires. Vous avez dit que votre alternative est . C'est très petitp=0.001

> ES.h(.001, 0)
[1] 0.0632561

Mais nous pouvons toujours continuer.

 > pwr.p.test(sig.level=0.05, power=.8, h = ES.h(.001, 0), alt="greater", n = NULL)

 proportion power calculation for binomial distribution (arcsine transformation) 

          h = 0.0632561
          n = 1545.124
  sig.level = 0.05
      power = 0.8
alternative = greater

En utilisant ces valeurs, vous avez besoin d'au moins 1546 observations.

Jessica
la source
1

Dans votre cas spécifique, il existe une solution simple et exacte:

Sous l'hypothèse nulle particulière vous ne devriez jamais observer de succès. Donc, dès que vous observez un succès, vous pouvez être sûr que .H0:p=0p0

Sous l'alternative Le nombre d'essais requis pour observer au moins 1 succès suit une distribution géométrique. Donc, afin d'obtenir la taille d'échantillon minimale pour atteindre une puissance de , vous devez trouver le plus petit k tel que,H1:p=0.0011β

1β1(1p)(k1)

Donc, avec pour obtenir puissance de , vous auriez besoin d'au moins 1610 échantillons.p=0.00180

flotte
la source
En lisant les commentaires de la solution 1, je me rends compte que c'est essentiellement la même solution que celle que vous obtenez si vous vous en tenez à y répondre. Néanmoins, cela ne nuit jamais à l'énoncé de certains résultats de base de la théorie des probabilités, sans qu'il soit nécessaire d'y arriver par intuition.
flotte le