tp4.utf8.md

Introduction

Pour commencer, lisez la partie du cours relative aux estimateurs et à l’estimation.

Illustration du biais et de la variance d’un estimateur:

biais et variance
source: http://scott.fortmann-roe.com/docs/BiasVariance.html

Exercice 1

1) Estimateurs de \(\mu\)

On souhaite tirer aléatoirement N échantillons de taille n d’une v.a \(X\) suivant une loi normale \(\mathcal N(\mu,1)\).

Choisir une valeur arbitraire de \(\mu\).

On choisit \(N=10000\) et \(n=100\).

On considère les estimateurs \(T_1\) des moyennes empiriques et \(T_2\) des médianes empiriques sur des échantillons de taille \(n\).

Les deux estimateurs estiment la valeur de \(\mu\) (rappelons que \(\mathbb E[X] = \mu\) et Mediane(\(X\)) = \(\mu\))

Evaluez \(N\) fois les estimateurs \(T_1\) et \(T_2\) puis estimer leurs biais, variances et erreurs quadratique moyenne pour les comparer.

# Helpers:
eval_mse <- function(x, theta) mean((x - theta)^2)
library(knitr) # Pour afficher des tableaux avec kable.

mu <- 1
n <- 100
N <- 10000
ech <- replicate(N, rnorm(n, mean=mu, sd=1))

# Estimateurs:
T1 <- apply(ech, 2, mean)
T2 <- apply(ech, 2, median)
Ts = cbind(T1, T2)

Biais = apply(Ts, 2, mean)- mu
Variance = apply(Ts, 2, var)
MSE <- apply(Ts, 2, eval_mse, theta=mu)
# Autrement (MSE = Biais^2 + Variance)
MSE_bis <- Biais^2 + Variance

kable(rbind(Biais, Variance, MSE, MSE_bis), digits=4)
T1 T2
Biais -3e-04 -0.0011
Variance 1e-02 0.0154
MSE 1e-02 0.0154
MSE_bis 1e-02 0.0154

L’erreur quadratique moyenne (EQM) ou encore MSE (de l’anglais mean squared error) permet de classer les estimateurs.

Vu que \(\rm{MSE}(T_1)\leq \rm{MSE}(T_2)\) alors l’estimateur \(T_1\) est meilleur que \(T_2\).

2) Estimateurs de \(\sigma^2\)

On choisit maintenant la loi \(\mathcal N(0,\sigma^2)\).

Choisir une valeur arbitraire de \(\sigma^2\).

On souhaite estimer l’écart type \(\sigma\) avec trois estimateurs :

  • L’estimateur \(T_1(X_1,\ldots, X_n) = \sqrt{\displaystyle \frac{1}{n-1} \displaystyle \sum_{k=1}^{n} (X_k-\bar X_n)^2}\).
  • La racine carrée de la variance empirique \(T_2(X_1,\ldots, X_n)=\sqrt{\displaystyle \frac{1}{n} \displaystyle \sum_{k=1}^{n} (X_k-\bar X_n)^2}\).
  • La distance interquartile (différence entre le troisième et le premier quartile) divisée par qnorm(0.75)-qnorm(0.25) : \(T_3(X_1,\ldots, X_n)=\displaystyle\frac{Q_3-Q_1}{1.34898}\).

Evaluez \(N\) fois les estimateurs \(T_1\), \(T_2\) et \(T_3\) puis comparer les. Visualiser des histogrammes (\(\color{blue}{\texttt{hist}}\)) ou des diagramme en boite (\(\color{blue}{\texttt{boxplot}}\)) des trois estimateurs.

sigma <- 1
n <- 100
N <- 10000
ech <- replicate(N, rnorm(n, mean=0, sd=sigma))

sd1 <- function(x) sqrt(sum((x - mean(x))^2) / (length(x) - 1))
# sd1 <- function(x) sqrt(var(x)) # Alternative
sd2 <- function(x) sqrt(mean((x - mean(x))^2))
# sd2 <- function(x) sqrt((n-1)/n * var(x))  # Alternative
sd3 <- function(x) (quantile(x,.75)-quantile(x,.25))/(qnorm(.75)-qnorm(.25))
# sd3 <- function(x) IQR(x) / (qnorm(.75)-qnorm(.25)) # Alternative

# Estimateurs
T1 <- apply(ech, 2, sd1)
T2 <- apply(ech, 2, sd2)
T3 <- apply(ech, 2, sd3)
Ts <- cbind(T1, T2, T3)

Biais = apply(Ts, 2, mean)- sigma
Variance = apply(Ts, 2, var)
MSE <- apply(Ts, 2, eval_mse, theta=sigma)
kable(rbind(Biais, Variance, MSE), digits=4)
T1 T2 T3
Biais -0.0033 -0.0083 -0.0145
Variance 0.0051 0.0050 0.0131
MSE 0.0051 0.0051 0.0133

Les deux estimateurs \(T_1\) et \(T_2\) sont presque equivalents, et les deux sont meilleurs que \(T_3\).

# Pour visualiser tous les plots sur la meme echelle on force des breaks identiques.
par(mfrow=c(3,1))
hist(T1, prob=T, breaks=seq(0.4,1.7, 0.05), col='red')
hist(T2, prob=T, breaks=seq(0.4,1.7, 0.05), col='blue')
hist(T3, prob=T, breaks=seq(0.4,1.7, 0.05), col='green')