From 156411965ddcc081950ae708b654460b41abbcdb Mon Sep 17 00:00:00 2001 From: Arthur DANJOU Date: Mon, 5 Jan 2026 16:53:31 +0100 Subject: [PATCH] =?UTF-8?q?Add=20TP1.Rmd=20for=20Statistiques=20Non=20Para?= =?UTF-8?q?m=C3=A9trique=20(M2)?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- M2/Statistiques Non Paramétrique/TP1.Rmd | 129 +++++++++++++++++++++++ 1 file changed, 129 insertions(+) create mode 100644 M2/Statistiques Non Paramétrique/TP1.Rmd diff --git a/M2/Statistiques Non Paramétrique/TP1.Rmd b/M2/Statistiques Non Paramétrique/TP1.Rmd new file mode 100644 index 0000000..67fa78e --- /dev/null +++ b/M2/Statistiques Non Paramétrique/TP1.Rmd @@ -0,0 +1,129 @@ + +```{r} +n = 100 +Z = sample.int(3, n, replace=T) +``` + + +```{r} +mu1 = -3 +mu2 = 0 +mu3 = 1 +sigma12 = 1 +sigma22 = 1 +sigma32 = 0.03 +X = rnorm(n, mu1, sqrt(sigma12)) +X[Z == 2] = rnorm(sum(Z == 2), mu2, sqrt(sigma22)) +X[Z == 3] = rnorm(sum(Z == 3), mu3, sqrt(sigma32)) +``` + + +```{r} +plot(density(X, bw = 0.1), ylim = c(0, 0.85), main = 'Estimation de la densité') +rug(X) +x = seq(-6, 2, length.out = 100) +points( + x, + (dnorm(x, mu1, sqrt(sigma12)) + + dnorm(x, mu2, sqrt(sigma22)) + + dnorm(x, mu3, sqrt(sigma32))) / + 3, + type = 'l', + col = 'darkgreen' +) +legend( + 'topleft', + c('Vraie densité', 'Estimation'), + lty = c(1, 1), + col = c('darkgreen', 1)) +``` + +La valeur 0. semble donner une meilleure estimation. + +```{r} +plot( + density(X, bw = 'SJ'), + ylim = c(0, 0.85), + main = 'Estimation de la densité' +) +lines(density(X, bw = 'ucv'), col = 2) +lines(density(X, bw = 'bcv'), col = 4) + +rug(X) +x = seq(-6, 2, length.out = 100) +points( + x, + (dnorm(x, mu1, sqrt(sigma12)) + + dnorm(x, mu2, sqrt(sigma22)) + + dnorm(x, mu3, sqrt(sigma32))) / + 3, + type = 'l', + col = 'darkgreen' +) +legend( + 'topleft', + c('Vraie densité', 'SJ', 'ucv', 'bcv'), + lty = rep(1, 4), + col = c('darkgreen', 1, 2, 4)) +``` + +La validation croisée non biaisée (UCV) semble donner les meilleurs résultats.. + + +```{r} +kmax = 50 +H = (1:kmax)^2 / (n * sqrt(2 * pi)) +hmin = min(H) +p = 5000 +fhat = matrix(NA, kmax, p) +fhatmin = density(X, bw = hmin, n = p) +x = fhatmin$x +fhat[1, ] = fhatmin$y +for (j in 2:kmax) { + fhat[j, ] = density(X, bw = H[j], n = p)$y +} +# Calcul du critere PCO +b = rep(NA, kmax) +v = rep(NA, kmax) +crit = rep(NA, kmax) +for (k in 1:kmax) { + b[k] = mean((fhat[k, ] - fhat[1, ])^2) + v[k] = 2 * mean(dnorm(x, sd = hmin) * dnorm(x, sd = H[k])) / n + crit[k] = b[k] + v[k] +} +khat = which.min(crit) +hhat = H[khat] +fhhat = fhat[khat, ] +plot( + x, + fhhat, + main = paste("Minimum atteint pour h=", signif(hhat, 2)), + xlab = 't', + ylab = expression(hat(f)[hat(h)]), + type = 'l' +) +points( + x, + (dnorm(x, mu1, sqrt(sigma12)) + + dnorm(x, mu2, sqrt(sigma22)) + + dnorm(x, mu3, sqrt(sigma32))) / + 3, + type = 'l', + col = 'darkgreen', + lwd = 2 +) +``` + +```{r} +plot(density(precip, bw = 'ucv')) +rug(precip) +``` + +Les données sont bimodales. On décèle deux groupes de villes : un groupe avec des précipitations situées autour de 15 inches avec une grande dispersion et un autre groupe avec des précipitations moyennes situées autours de 40 inches. + +```{r} +plot(density(faithful[, 1], bw = 'ucv')) +rug(faithful[, 1]) +``` + +On distingue là aussi deux groupes : un groupe de geysers ayant un temps d'éruption autour de 2 mins et un autre 4 mins/4 mins 30. La règle du pouce donne ds estimateurs très différents de la validation croisée et de la méthode de Sheather et Jones. \ No newline at end of file