Puissance et tests statistiques (simples)

Pour comprendre la notion de puissance dans un test statistique, reprenons Statistical power analysis for the behavioral sciences de Jacob Cohen. Comme il le rappelle Fisher (dans The Design of Experiments) ne souhaite pas “prouver” qu’une hypothèse (appelée H_0) est valide, mais il espère, en réalité, rejeter l’hypothèse en question (c’est ce que rappellent d’ailleurs Hubbard et Bayarri en 2003). Autrement dit, supposons que l’on cherche à faire un test médical, et à détecter si une personne est malade (ou pas). Si une personne est généralement malade dès lors qu’un taux de quelque chose est faible (noté https://latex.codecogs.com/gif.latex?\theta), on va chercher à tester https://latex.codecogs.com/gif.latex?H_0:\theta_\star%3E\theta_0. Dans un test statistique, on espère rejeter l’hypothèse nulle https://latex.codecogs.com/gif.latex?H_0.

The power of a statistical test of a null hypothesis is the probability that it will lead to the rejection of the null hypothesis.

C’est (presque) la seule chose à retenir de cette section sur le test, le reste étant essentiellement du calcul. De manière générale, on fait un test de

https://latex.codecogs.com/gif.latex?H_0:\boldsymbol{\theta}_\star\in\Theta_0 contre https://latex.codecogs.com/gif.latex?H_1:\boldsymbol{\theta}_\star\in\Theta_1

où https://latex.codecogs.com/gif.latex?\boldsymbol{\theta} est un ensemble de paramètres d’intérêt. Bon, pour commencer, on va supposer que https://latex.codecogs.com/gif.latex?\boldsymbol{\theta} est un paramètre uni-dimensionnel (ça permet de faire des dessins, et des exemples simples, mais après l’idée est la même). Par exemple, faisons un test de moyenne, du genre

https://latex.codecogs.com/gif.latex?H_0:\mathbb{E}(X)\geq%20m contre https://latex.codecogs.com/gif.latex?H_1:\mathbb{E}(X)%3C%20m

Ou plus simple pour commencer,

https://latex.codecogs.com/gif.latex?H_0:\mathbb{E}(X)=%20m_0 contre https://latex.codecogs.com/gif.latex?H_1:\mathbb{E}(X)=%20m_1

(on fait, on verra qu’on va toujours de ramener à ce dernier cas… en gros, c’est le seul cas qu’on sait traiter). Et pour simplifier davantage, on va supposer que les observations sont Gaussiennes, de variance connue (que l’on peut fixer à 1, éventuellement). Bon, techniquement, ne pas connaître la variance ne change rien à la méthodologie. La seule chose est que la statistique de test ne suit plus une loi normale, mais une de Student à https://latex.codecogs.com/gif.latex?n-1 degrés de liberté. Bref, c’est un détail pour comprendre ce qu’on cherche à faire ici.

On a ici un test de deux hypothèses simples. Dans ce cas, on a le premier théorème de Neyman-Pearson, qui nous dit que pour un degré de significativité donné https://latex.codecogs.com/gif.latex?\alpha, on peut construire un test de puissance maximale avec une région critique de la forme

https://latex.codecogs.com/gif.latex?W=\left\{\boldsymbol{x}\in\mathbb{R}^n,\frac{\mathcal{L}(\theta_1,\boldsymbol{x})}{\mathcal{L}(\theta_0,\boldsymbol{x})}\geq%20k\right\}

Bref, la région de rejet peut s’exprimer en fonction de la moyenne de nos observations, https://latex.codecogs.com/gif.latex?\overline{y}, car le rapport de vraisemblance est une fonction croissante de la moyenne (on en reparlera bientôt…). Revenons un instant sur les deux quantités importantes dont on vient de parler : le niveau de significativité, et la puissance (on les avait défini en cours, en introduction aux tests statistiques). Le niveau de significativité est la probabilité de rejeter à tort l’hypothèse.

Dans le cas où on a deux hypothèses simples, et que , on aura tendance à rejeter https://latex.codecogs.com/gif.latex?H_0 si la moyenne est trop faible. La région de rejet est alors de la forme . On en reparlera un peu, mais le signe de  interviendra régulièrement, et aura la facheuse tendance à inverser le sens de nos inégalités. Bon, reste à comprendre ce que signifie ce conditionnement par https://latex.codecogs.com/gif.latex?H_0. Car si ça ne pose pas de soucis pour une hypothèse simple, pour une hypothèse composite, c’est plus tordu… mais on en reparlera…

Aussi, dans notre cas, on doit expliciter

https://latex.codecogs.com/gif.latex?\mathbb{P}(\overline{Y}%3Cc%20\vert%20H_0)=\alpha

soit

https://latex.codecogs.com/gif.latex?\mathbb{P}\left(\left.\sqrt{n}~\frac{\overline{Y}-m}{\sigma}%3C\sqrt{n}~\frac{c-m}{\sigma}%20\right\vert%20H_0\right)=\alpha

alors

https://latex.codecogs.com/gif.latex?\mathbb{P}\left(\left.Z%3C\sqrt{n}~\frac{c-m}{\sigma}%20\right\vert%20H_0\right)=\alpha

i.e.

https://latex.codecogs.com/gif.latex?c=m_0+\Phi^{-1}(\alpha)\frac{\sigma}{\sqrt{n}}

Notons d’ailleurs que c est fonction de

  • https://latex.codecogs.com/gif.latex?m_0 la moyenne sous https://latex.codecogs.com/gif.latex?H_0
  • https://latex.codecogs.com/gif.latex?n la taille de l’échantillon
  • https://latex.codecogs.com/gif.latex?\alpha le seuil de significativité

Mais le seuil ne dépend pas de l’hypothèse alternative, par exemple.

On détermine c à partir du niveau de significativité https://latex.codecogs.com/gif.latex?\alpha. Aussi, en pratique, on calcule https://latex.codecogs.com/gif.latex?\overline{y}, et le seuil. Si https://latex.codecogs.com/gif.latex?\overline{y}%3Cc, on va rejeter https://latex.codecogs.com/gif.latex?H_0.
Si https://latex.codecogs.com/gif.latex?H_0 est vraie, on risque de se tromper avec probabilité https://latex.codecogs.com/gif.latex?\alpha (c’est le risque de première espèce).

Maintenant, si https://latex.codecogs.com/gif.latex?H_1 était vraie ? Car pour le moment, on n’a jamais parlé de cette hypothèse alternative. La puissance est la probabilité de rejeter – avec raison – l’hypothèse https://latex.codecogs.com/gif.latex?H_0. On la note https://latex.codecogs.com/gif.latex?1-\beta (où https://latex.codecogs.com/gif.latex?\beta est la probabilité de commettre une erreur de seconde espèce).

https://latex.codecogs.com/gif.latex?\mathbb{P}_1(W)=\mathbb{P}(W\vert%20H_1)=\mathbb{P}(\overline{Y}%3Cc\vert%20H_1)

Pour conditionner par https://latex.codecogs.com/gif.latex?H_1, là encore, on peut le faire dans le cas d’une hypothèse simple, à savoir https://latex.codecogs.com/gif.latex?H_1:\mathbb{E}(X)=m_1. Sinon, on va conditionner par toutes les valeurs possibles, et tracer une fonction puissance.

On va introduire une notion appelée taille d’effet, qui est l’idée intuitive, ou a priori, que l’on se fait sur la distance entre la moyenne sous https://latex.codecogs.com/gif.latex?H_0 et la vraie moyenne https://latex.codecogs.com/gif.latex?\mathbb{E}(X). Dans le cas de test entre hypothèses simples, j’avais dit au début du billet qu’on cherche à rejeter . On pose alors https://latex.codecogs.com/gif.latex?H_0, donc https://latex.codecogs.com/gif.latex?m_1 est justement censé refléter l’idée qu’on se fait de https://latex.codecogs.com/gif.latex?\mathbb{E}(X).  On pose alors

https://latex.codecogs.com/gif.latex?d=\frac{\mathbb{E}(X)-\mu_0}{\sigma}

ou

https://latex.codecogs.com/gif.latex?d=\frac{m_1-\mu_0}{\sigma}

si on a deux hypothèses simples. La puissance est alors (je reviens sur notre comparaison entre deux moyennes)

https://latex.codecogs.com/gif.latex?1-\beta%20=%20\mathbb{P}(Z%3C\Phi^{-1}(\alpha)-d\sqrt{n})

(où classiquement https://latex.codecogs.com/gif.latex?Z\sim%20\mathcal{N}(0,1)). On peut faire le calcul à la main (sous R, quoi): si on a 10 observations, que l’on veut faire un test avec un écart de 1 entre les moyennes sous les deux hypothèses et l’on utilise un seuil de significativité de 5%,

> a=0.05
> d=-1
> n=10

la puissance est alors

> pnorm(qnorm(a)-d*sqrt(10))
[1] 0.9354202

mais on peut aussi utiliser une fonction de R

> library(pwr)
> pwr.norm.test(sig.level=a,d,n,
+ alternative="less")

     Mean power calculation for normal distribution with known variance 

              d = -1
              n = 10
      sig.level = 0.05
          power = 0.9354202
    alternative = less

On peut aussi faire un joli dessin,

> d=(seq(minf,msup,length=256)-m0)/sigma
> plot(seq(minf,msup,length=256),
+ pwr.norm.test(sig.level=a,d,n,alternative=
+ "less")$power,type="l",col="red",lwd=2)

Il n’est pas rare que l’on souhaite une puissance au moins égale à 80%. Comme on le voit sur le dessin, cela signifie qu’il faut avoir une moyenne différence entre les moyennes d’au plus -0.75.  On peut inverser la lecture en prenant

> d=(seq(minf,msup,length=256)-m0)/sigma
> plot(seq(minf,msup,length=256),
+ pwr.norm.test(sig.level=a,d,n,alternative=
+ "greater")$power,type="l",col="red",lwd=2)

Ce qui est intéressant, c’est que si d est trop faible, on peut utiliser un autre paramètre pour augmenter la puissance du test : la taille de l’échantillon. On avait vu que

https://latex.codecogs.com/gif.latex?1-\beta%20=%20\mathbb{P}(Z%3C\Phi^{-1}(\alpha)-d\sqrt{n})

aussi, il faut avoir

https://latex.codecogs.com/gif.latex?n=\left(\frac{\Phi^{-1}(1-\alpha)+\Phi^{-1}(1-\beta)}{d}\right)^2
Là encore, on peut utiliser la fonction de R pour visualiser l’impact de la taille de l’échantillon sur la puissance (par exemple en passant de 10 à 25 observations)

> lines(seq(minf,msup,length=256),
+ pwr.norm.test(sig.level=a,d,25,alternative=
+ "greater")$power,type="l",col="blue",lwd=2)

 

On peut considérer un autre test, avec une hypothèse bilatéral, par exemple, où https://latex.codecogs.com/gif.latex?H_1 est de la forme https://latex.codecogs.com/gif.latex?\mathbb{E}(X)\neq%20m.

> plot(seq(minf,msup,length=256),
+ pwr.norm.test(sig.level=a,d,10,alternative=
+ "two.sided")$power,type="l",col="red",lwd=2)

qui est assez intuitif quand on y pense : si la distance entre les deux moyennes est faible (différence proche de 0), le test n’est pas très puissance, et on risque de commettre beaucoup d’erreurs de seconde espèce, pour maintenant un bon niveau pour le taux d’erreur de première espèce (ici 5%).  En effet, si on accepte un plus fort taux d’erreur de première espèce (par exemple 20%), on augmente la puissance du test (en baisse le taux d’erreur de second espèce),

> lines(seq(minf,msup,length=256),
+ pwr.norm.test(sig.level=.2,d,10,alternative=
+ "two.sided")$power,type="l",col="blue",lwd=2)

On rejoint le point que je faisais il y a quelques semaines sur les tests : il y a une balance entre les erreurs de deux types. On ne peut pas gagner sur les deux tableaux !

Mercredi on fera un peu plus d’exercices, et on abordera les tests composites (je mettrais alors un billet en ligne, si j’ai le temps…)


Leave a Reply

Your email address will not be published. Required fields are marked *

This site uses Akismet to reduce spam. Learn how your comment data is processed.