Niveau et puissance de test (paramétrique)

Lors du dernier cours de statistique, nous avons vu comme se construisent des tests, et comment on compare leur puissance. Considérons pour cela des observations \{x_1,\cdots,x_n\} tirées suivant une loi uniforme \mathcal{U}([0,\theta]). Et on se pose la question de savoir si l’hypothèse H_0:\theta\leq\theta_0 est rejetée, ou acceptée, avec \theta_0=1 par exemple. L’hypothèse alternative sera ici le complémentaire, autrement dit H_1:\theta>\theta_0.

Comme on l’a vu, la construction d’un test se fait en utilisant une statistique de test (ou plus simplement un estimateur du paramètre inconnu), avec \psi(\mathbf{x})=\mathbf{1}(\widehat{\theta}(\mathbf{x}))\leq c. On prend la décision D_1 (consistant à rejeter H_0) si \psi(\mathbf{x})=1, et la décision D_0 (consistant à ‘accepter’ H_0) si \psi(\mathbf{x})=0.

Le risque d’erreur de première espèce est basé sur \mathbb{E}[\psi(\mathbf{X})]. Si l’hypothèse H_0 est une hypothèse simple H_0:\theta=\theta_0, le risque d’erreur de première espèce est https://latex.codecogs.com/gif.latex?\alpha=\mathbb{E}[\psi(\mathbf{X})\vert%20X_1,\cdots,X_n\sim%20F_{\theta_0}], soit \alpha=\mathbb{P}[\widehat{\theta}(\mathbf{X})>c \vert X_1,\cdots,X_n\sim F_{\theta_0}] pour un seuil c. On va déterminer le seuil c de manière à ce que le risque de première espèce soit 5%.

Pour cela, il nous fait un estimateur \widehat{\theta} dont on connaîtrait la loi. Commençons par l’estimateur du maximum de vraisemblance. Rappelons que l’estimateur du maximum de vraisemblance est \widehat{\theta}=\max\{x_i\}. Nous avons vu – dans la section sur les statistiques d’ordre – que \mathbb{P}[\max\{X_i\}]=\left(\frac{x}{\theta}\right)^n

pour x\leq\theta. Aussi, pour un seuil c, \mathbb{P}[\widehat{\theta}>c]=1-\left(\frac{c}{\theta}\right)^n

pour c\leq\theta. Le risque de première espèce \alpha=\underset{\theta\leq\theta_0}{\max}\left\lbrace\mathbb{P}\big[\widehat{\theta}\geq c\vert X\sim F_{\theta}\big]\right\rbrace

soit \alpha=\underset{\theta\leq\theta_0}{\max}\left\lbrace1-\left(\frac{c}{\theta}\right)^n\right\rbrace

qui atteindra son maximum en \theta_0. Aussi \alpha=1-\left(\frac{c}{\theta_0}\right)^n

soit c=\theta_0\cdot[1-\alpha]^{1/n}. Les régions d’acceptation et de rejet sont les suivantes,

Pour calculer l’erreur de seconde espèce, utilisons\beta=\mathbb{P}[D_0\vert H_1\text{vraie}] soit
\beta=\left(\frac{c}{\theta}\right)^n=(1-\alpha)\left(\frac{\theta_0}{\theta}\right)^n
\theta>\theta_0 puisque nous sommes sous H_1. La puissance du test, basé sur le maximum des observations (correspondant au maximum de vraisemblance), est alors 1-\left(\frac{c}{\theta}\right)^n=(1-\alpha)\left(\frac{\theta_0}{\theta}\right)^n

On obtient graphiquement

theta=seq(1,3,by=.1)
n=5
puissance=1-(1-.05)*(1/theta)^n
plot(theta,puissance,type="l",lwd=2,col="red")

et on note que la puissance vaut \alpha quand \theta\rightarrow\theta_0 et 1 lorsque \theta\rightarrow\infty.

Et que se passe-t-il si on avait utilisé l’estimateur de la méthode des moments, et non plus l’estimateur du maximum de vraisemblance ? Rappelons que si U\sim\mathcal{U}([0,\theta]) alors \mathbb{E}[U]=\theta/2. Aussi, l’estimateur de la méthode des moments est \widetilde{\theta}=2\overline{x}. Le risque de première espèce est \alpha=\mathbb{P}(2\overline{X}>c\vert X\sim F_{\theta_0})

Le soucis est que la loi de \overline{X} n’est pas simple à obtenir, analytiquement. Mais on peut l’avoir par simulations.

> ns = 1e5
> n = 5
> X = matrix(runif(ns*n),ns,n)
> Xbar = apply(X,1,mean)
> hist(Xbar,probability=TRUE)
> lines(density(Xbar),lwd=2,col="red")

La fonction de survie de \widehat{\theta} s’approche avec la fonction de survie empirique, obtenue sur 100,000 simulations.

> x = seq(0,2,by=.001)
> y = Vectorize(function(t) mean(2*Xbar>t))(x)
> plot(x,y,type="l")
> abline(h=.05,col="red",lty=2)

On peut alors déterminer le niveau du seuil
\alpha valant 5% en utilisant le quantile empirique

> (c = quantile(2*Xbar,.95))
95%
1.427396

Ce seuil permet de construire un test de risque de première espèce 5%. On peut alors calculer la puissance de ce test, numériquement,

> ns = 1e4
> Puissance = function(t){
+ X = matrix(runif(ns*n,0,t),ns,n)
+ Xbar = apply(X,1,mean)
+ mean(2*Xbar > c)}
> theta =s eq(1,3,by=.1)
> puissance = Vectorize(Puissance)(theta)
> lines(theta,puissance,type="l",lwd=2,col="blue")

On a ainsi construit un test de même niveau que le précédant (risque de première espèce), mais de puissance plus faible, tout le temps.


One thought on “Niveau et puissance de test (paramétrique)”

Leave a Reply

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