Begin'R
Les statistiques avec R
Navigation
[Retour au sommaire]
# Exercices bilan :Objectifs Evaluer ses compétences sur la programmation en langage R. <style> div.pseudocode pre.r { background-color:#EEEEEE; } </style> :Exercice : Calcul de la somme des éléments d'un tableau {#exercice_sommeTableau, toggle=collapse, title-display=show} L'objectif de cet exercice est de créer un algorithme qui utilise une structure itérative. * Le pseudo-code suivant permet à partir d'un tableau `Tab` saisi par l'utilisateur de calculer et de renvoyer en sortie la somme des éléments de ce tableau. L'instruction `taille()` dans le pseudo-code permet de compter le nombre de cases d'un tableau. <div class = "pseudocode"> ```r Début Affectation de Tab Somme <- 0 t <- taille(Tab) i <- 1 Tant que (i <= t) faire Somme <- Somme + Tab[i] i <- i + 1 Fin du Tant que Afficher Somme Fin ``` </div> * Traduire cet algorithme en langage R en créant une fonction intitulée `'calcul_somme_tableau()'`. La fonction R nommée **`length()`** permet de compter le nombre de cases d'un tableau. * Valider cette fonction `'calcul_somme_tableau()'` en prenant `Tab = c(1,3,6)` et vérifier le résultat à l'aide de la fonction **`sum()`** qui calcule également la somme des éléments d'un tableau. :Corrigé {#corrige_sommeTableau, toggle=collapse, title-display=show} * Traduction en langage R : code de la fonction `'calcul_somme_tableau()'` à implémenter dans le fichier `'calcul_somme_tableau.r'` ```r calcul_somme_tableau <- function (Tab){ Somme <- 0 t <- length(Tab) i <- 1 while (i <= t) { Somme <- Somme + Tab[i] i <- i + 1 } return(Somme) } ``` * Utilisation de la fonction `'calcul_somme_tableau()'` (dans un autre script ou dans la console) et vérification du résultat avec la fonction `'sum()'` ```r source('calcul_somme_tableau.r') # déclaration de la fonction 'calcul_somme_tableau()' Tab<-c(1,3,6) Somme <- calcul_somme_tableau(Tab) print(Somme) ``` ``` ## [1] 10 ``` ```r # Vérification avec la fonction `sum()` sum(Tab) ``` ``` ## [1] 10 ``` :Exercice : Calcul de la somme des éléments d'un tableau supérieurs à un seuil {#exercice_sommeTableauseuil, toggle=collapse, title-display=show} L'objectif de cet exercice est de créer un algorithme qui associe une structure alternative et une structure itérative. * Le pseudo-code suivant permet à partir d'un tableau `Tab` et d'un seuil `S` saisis par l'utilisateur de calculer et de renvoyer en sortie la somme des éléments de ce tableau qui sont strictement supérieurs au seuil `S`. <div class = "pseudocode"> ```r Début Affectation de Tab Affectation de S Somme <- 0 t <- taille(Tab) i <- 1 Tant que (i <= t)faire Si (Tab[i]>S) alors Somme <- Somme + Tab[i] Fin du Si i <- i + 1 Fin du Tant que Afficher Somme Fin ``` </div> * Traduire cet algorithme en langage R en créant une fonction intitulée `'calcul_somme_tableau_seuil()'`. * Valider cette fonction `'calcul_somme_tableau_seuil()'` en prenant `Tab=c(1,3,6)` et `S=2`. :Corrigé {#corrige_sommeTableauseuil, toggle=collapse, title-display=show} * Traduction en langage R : code de la fonction `'calcul_somme_tableau_seuil()'` à implémenter dans le fichier `'calcul_somme_tableau_seuil.r'` ```r calcul_somme_tableau_seuil <- function (Tab, S){ Somme <- 0 t <- length(Tab) i <- 1 while (i <= t) { if (Tab[i]>S){ Somme <- Somme + Tab[i] } i <- i + 1 } return(Somme) } ``` * Utilisation de la fonction `'calcul_somme_tableau_seuil()'` (dans un autre script ou dans la console) ```r source('calcul_somme_tableau_seuil.r') # déclaration de la fonction 'calcul_somme_tableau_seuil()' Tab <- c(1,3,6) S <- 2 Somme <- calcul_somme_tableau_seuil(Tab, S) print(Somme) ``` ``` ## [1] 9 ``` * Vérification du résultat à l'aide de la fonction `sum()`. ```r # Vérification avec la fonction `sum()` sum(Tab[Tab>S]) ``` ``` ## [1] 9 ``` :Exercice : Comparaison de la moyenne de deux séries non appariées {#exercice_comparaisonmoyennes, toggle=collapse, title-display=show} L'objectif de cet exercice est de créer un algorithme qui fait appel à quelques notions de statistiques (test d'hypothèses). * Le pseudo-code suivant permet de comparer deux séries non appariées à l'aide d'un test de Student en vérifiant au préalable l'hypothèse d'homoscédasticité. L'utilisateur saisit en entrée les deux séries de données dans deux tableaux `Tab1` et `Tab2`, ainsi que la valeur `alpha` du risque de première espèce. L'algorithme renvoie en sortie une variable `diff` qui contient le résultat de l'interprétation statistique, à savoir `'Les moyennes ne sont pas significativement différentes'` ou `'Les moyennes sont significativement différentes'`. <div class = "pseudocode"> ```r Début Affectation de Tab1 Affectation de Tab2 Affectation de alpha Vérification de l égalité des variances (`var.test()`) avec récupération de la p-value p Si (p>alpha) alors Mise en place d un test de Student-t avec hypothèse d homosédasticite, et récupération de la p-value p Sinon Mise en place d un test de Student-t sans hypothèse d homosédasticite , et récupération de la p-value p Fin du Si Si(p>alpha) diff <- 'Les moyennes ne sont pas significativement différentes' Sinon diff <- 'Les moyennes sont significativement différentes' Fin du Si Afficher diff Fin ``` </div> * Traduire cet algorithme en langage R en créant une fonction intitulée `'comparaison_moyennes'`. * Valider cette fonction `'comparaison_moyennes()'` en prenant `Tab1 = c(1.5, 1.7, 1.4, 1.6, 1.4)`, `Tab2 = c(1.3, 1.9, 1.4, 1.5)` et `alpha = 0.05`. [Aide](#aide_comparaisonmoyennes) :Aide{#aide_comparaisonmoyennes, toggle=popup} Si l'objet `test` est le résultat d'un test statistique (par exemple `test <- t.test(...)`), on récupère la p-value de ce test statistique en utilisant l'instruction `p <- test$p.value`. :Corrigé {#corrige_comparaisonmoyennes, toggle=collapse, title-display=show} * Traduction en langage R : code de la fonction `'comparaison_moyennes()'` à implémenter dans le fichier `'comparaison_moyennes.r'` ```r comparaison_moyennes <- function (Tab1, Tab2, alpha){ test <- var.test(Tab1, Tab2, paired = FALSE) p <- test$p.value # récupération de la p-value du test de Fisher if (p>alpha){ test <- t.test(Tab1, Tab2, var.equal = TRUE) } else{ test <- t.test(Tab1, Tab2, var.equal = FALSE) } p <- test$p.value # récupération de la p-value du test de Student if(p>alpha){ resultat <- 'Les moyennes ne sont pas significativement différentes' } else{ resultat <- 'Les moyennes sont significativement différentes' } print(resultat) } ``` * Utilisation de la fonction `'comparaison_moyennes()'` (dans un autre script ou sur la console) ```r source('comparaison_moyennes.r') # déclaration de la fonction 'comparaison_moyennes()' Tab1 <- c(1.5, 1.7, 1.4, 1.6, 1.4) Tab2 <- c(1.3, 1.9, 1.4, 1.5) alpha <- 0.05 comparaison_moyennes(Tab1, Tab2, alpha) ``` ``` ## [1] "Les moyennes ne sont pas significativement différentes" ``` :Exercice : Comparaison de la moyenne de deux séries non appariées avec vérification de l'hypothèse de normalité {#exercice_shapiro, toggle=collapse, title-display=show} Cet exericie complète le précédent en vérifiant au préalable si l'hypothèse de normalité des observations est valide. Pour cela, un test de Shapiro sur chacune des séries de données est mis en place. Lorsque l'hypothèse de normalité est validée, la comparaison des moyennes s'effectue par le biais d'un test de Student. Dans le cas contraire, il faut s'orienter vers un test de Wilcoxon. Dans le cas où l'on s'oriente vers un test de Student, l'hypothèse d'homoscédasticité sera vérifiée au préalable. * Le pseudo-code suivant permet de comparer deux séries non appariées à l'aide d'un test de comparaison de moyennes. L'utilisateur saisit en entrée les deux séries de données dans deux tableaux `Tab1` et `Tab2`, ainsi que la valeur `alpha` du risque de première espèce. L'algorithme renvoie en sortie une variable `diff` qui contient le résultat de l'interprétation statistique, à savoir `'Les moyennes ne sont pas significativement différentes'` ou `'Les moyennes sont significativement différentes'`. <div class = "pseudocode"> ```r Début Affectation de Tab1 Affectation de Tab2 Affectation de alpha Vérification de la normalité (shapiro.test) avec récupération des p-value p1 et p2 pour les deux tableaux Si ((p1>alpha) && (p2>alpha)) alors vérification de l égalité de variance (var.test) et récupération de la p-value p Si (p>alpha) alors Test de Student avec hypothèse d homosédasticite et récupération de la p-value p Sinon Test de Student sans hypothèse d homosédasticite et récupération de la p-value p Fin du Si Sinon Test de Wilcoxon avec récupération de la p-value p Fin du Si Si(p>alpha) diff <- 'Les moyennes ne sont pas significativement différentes' Sinon diff <- 'Les moyennes sont significativement différentes' Fin du Si Afficher diff Fin ``` </div> * Traduire cet algorithme en langage R en créant une fonction intitulée `'compare_moyennes()'`. * Valider cette fonction `'compare_moyennes()'` en prenant `Tab1 = c(1.5, 1.7, 1.4, 1.6, 1.4)`, `Tab2 = c(1.3, 1.9, 1.4, 1.5)` et `alpha = 0.05`. [Aide](#aide_shapiro) :Aide{#aide_shapiro, toggle=popup} Sur R, le test de normalité de Shapiro est réalisé à l'aide de la fonction **`shapiro.test()`**. :Corrigé {#corrige_shapiro, toggle=collapse, title-display=show} * Traduction en langage R : code de la fonction `'compare_moyennes()'` à implémenter dans le fichier `'compare_moyennes.r'` ```r compare_moyennes <- function (Tab1, Tab2, alpha){ test_shapiro1 <- shapiro.test(Tab1) test_shapiro2 <- shapiro.test(Tab2) p1 <- test_shapiro1$p.value # récupération de la p-value du premier test de Shapiro p2 <- test_shapiro2$p.value # récupération de la p-value du deuxième test de Shapiro if ((p1>alpha) && (p2>alpha)){ test <- var.test(Tab1, Tab2, paired = FALSE) p <- test$p.value # récupération de la p-value du test de Fisher if (p>alpha){ test <- t.test(Tab1, Tab2, var.equal = TRUE) } else{ test <- t.test(Tab1, Tab2, var.equal = FALSE) } } else{ test <- wilcox.test(Tab1, Tab2) } p <- test$p.value # récupération de la p-value du test de comparaison de moyennes if(p>alpha){ resultat <- 'Les moyennes ne sont pas significativement différentes' } else{ resultat <- 'Les moyennes sont significativement différentes' } return(resultat) } ``` * Utilisation de la fonction `'compare_moyennes()'` dans un autre script ou sur la console ```r source('compare_moyennes.r') # déclaration de la fonction 'compare_moyennes()' Tab1 <- c(1.5, 1.7, 1.4, 1.6, 1.4) Tab2 <- c(1.3, 1.9, 1.4, 1.5) alpha <- 0.05 resultat <- compare_moyennes(Tab1, Tab2, alpha) print(resultat) ``` ``` ## [1] "Les moyennes ne sont pas significativement différentes" ``` :Exercice : Estimation d'une probabilité par tirage Monte Carlo {#exercice_calculprobaMonteCarlo, toggle=collapse, title-display=show} On suppose que `X` suit une loi normale centrée réduite et on souhaite calculer la probabilité que `X` soit strictement supérieur à `S`. Pour cela, on procèdera par une méthode de Monte Carlo qui consiste à générer un échantillon de taille `N` selon une loi normale centrée réduite et à estimer la probabilité `P(X>S)` en calculant le pourcentage de valeurs générées supérieures à `S`. Le pseudo-code suivant permet à partir des valeurs de `S` et `N` saisies par l'utilisateur de calculer `P(X>S)`. <div class = "pseudocode"> ```r Début Affectation de S Affectation de N X <- génération d'un échantillon de taille N selon une loi normale normale centrée réduite P <- 0 Pour i allant de 1 à N par pas de 1 faire Si (X[i]>S) Alors P <- P + 1 Fin du Si Fin du Pour P <- P / N Sortir P Fin ``` </div> Questions : * Traduire cet algorithme en langage R en créant une fonction intitulée `'estime_proba_monte_carlo()'`. * Valider cette fonction `'estime_proba_monte_carlo()'` en prenant `S=0.5` et `N=10^6`. [Aide](#aide_calculprobaMonteCarlo) :Aide{#aide_calculprobaMonteCarlo, toggle=popup} Sur R, la génération d'un échantillon selon une loi normale est réalisé à l'aide de la fonction `rnorm()`. :Corrigé {#corrige_calculprobaMonteCarlo, toggle=collapse, title-display=show} * Traduction en langage R : code de la fonction `'estime_proba_monte_carlo()'` à implémenter dans le fichier `'estime_proba_monte_carlo.r'` ```r estime_proba_monte_carlo <- function (S, N){ # génération d'un échantillon de taille N selon une loi normale centrée réduite X <- rnorm(N, mean = 0, sd = 1) P <- 0 for(i in seq(from=1,to=N,by=1)){ if(X[i]>S){ P <- P + 1 } } P <- P / N # ou plus rapidement # P <- length(X[X>S]) / N return(P) } ``` * Utilisation de la fonction `'estime_proba_monte_carlo()'` dans un autre script ou sur la console. ```r S <- 0.5 # Seuil pour le calcul de P(X>S) # déclaration de la fonction 'estime_proba_monte_carlo()' source('estime_proba_monte_carlo.r') N <- 10^6 # Taille de l'échantillon généré P <- estime_proba_monte_carlo (S, N) ``` ``` ## Error in estime_proba_monte_carlo(S, N): impossible de trouver la fonction "estime_proba_monte_carlo" ``` ```r P ``` ``` ## Error in eval(expr, envir, enclos): objet 'P' introuvable ``` Ce calcul de probabilité peut être vérifié à l'aide de la fonction `pnorm()` qui renvoie la fonction de répartition d'une loi normale. ```r # Vérification avec pnorm # P(X>S) <- 1 - P(X<=S) P <- 1 - pnorm(S, mean = 0, sd = 1) P ``` ``` ## [1] 0.3085375 ``` :Exercice : Estimation d'un intervalle de confiance par bootstrap {#exercice_bootstrap, toggle=collapse, title-display=show} Cet exercice permet d'estimer un intervalle de confiance sur la moyenne à partir d'une série `X` de `n` observations qui ne sont pas distribuées selon une loi Gaussienne. Pour cela, une méthode de bootstrap est utilisé. Le principe est le suivant. Tout d'abord, `m` ré-échantillons (ou échantillons bootstrap) sont simulés. Ces ré-échantillons de même taille que la série initiale `X` sont obtenus par échantillonnage avec remise dans la série. Pour chaque ré-échantillon, la valeur moyenne est calculée et stockée dans un tableau de taille `m`. L'intervalle de confiance `1-alpha` sur la moyenne est ensuite obtenu en calculant les quantiles `alpha/2` et `1-alpha/2` de ce tableau. `alpha` représente ici le risque que la moyenne se situe à l'extérieur de l'intervalle. Les valeurs `alpha=0.05` ou `alpha=0.01` sont communément choisies. Le pseudo-code suivant décrit le calcul de l'intervalle de confiance sur la moyenne à partir d'une série `X` de `n` observations, de la valeur `m` du nombre de ré-échantillons et du risque `alpha` saisis par l'utilisateur. <div class = "pseudocode"> ```r Début Affectation de X Affectation de m Affectation de alpha n <- taille(X) Pour i allant de 1 à m par pas de 1 faire bt <- tirage à partir de X d un échantillon boostrap avec remise de taille n moy[i] <- moyenne(bt) Fin du Pour IC_min <- quantile alpha/2 de moy IC_max <- quantile 1- alpha/2 de moy Afficher ICmin Afficher ICmax Fin ``` </div> Questions : * Traduire cet algorithme en langage R en créant une fonction intitulée `'bootstrap()'`. L'affichage de l'intervalle de confiance sera réalisé à l'intérieur de la fonction. * Valider cette fonction `'bootstrap()'` en considérant que `X` est la série de `n=15` observations définie par `c(0.07, 0.19, 0.41, 2.54, 2.27, 1.71, 0.65, 0.57, 0.08, 0.54, 0.06, 1.92, 2.11, 0.18, 1.03)`. Prendre un nombre de ré-échantillons `m=100000` et un risque `alpha` de première espèce égal à `0.05`. [Aide](#aide_bootstrap) :Aide{#aide_bootstrap, toggle=popup} * La fonction **`sample()`** permet d'effectuer une échantillonage avec remise. :Corrigé {#corrige_bootstrap, toggle=collapse, title-display=show} * Traduction en langage R : code de la fonction `'bootstrap()'` à implémenter dans le fichier `'bootstrap.r'` ```r bootstrap <- function (X, m, alpha){ n <- length(X) # Initialisation d'un tableau de m moyennes moy<-rep(0,m) # Simulation des m ré-échantillons et de leurs moyennes # Les ré-échantillons de même taille que la série initiale # sont obtenus par échantillonnage avec remise dans la série # Les ré-échantillons successifs ne sont pas préservés # Seule la moyenne est préservée for (i in seq(from=1,to=m,by=1)){ bt <- sample(X,size=n,replace=TRUE) moy[i] <- mean(bt) } # Extraction des quantiles alpha/2 et 1-alpha/2 qq <- quantile(moy,probs=c(alpha/2,1-alpha/2)) qq <- round(qq,2) # on arrondit à 2 chiffres après la virgule ICmin <- qq[1] ICmax <- qq[2] # Affichage chaine <- paste("Intervalle de confiance : [",ICmin,",",ICmax,"]", sep = "") print(chaine) } ``` * Utilisation de la fonction `'bootstrap()'` (dans un autre script ou sur la console) ```r # déclaration de la fonction 'bootstrap()' source('bootstrap.r') # échantillon X <- c(0.07, 0.19, 0.41, 2.54, 2.27, 1.71, 0.65, 0.57, 0.08, 0.54, 0.06, 1.92, 2.11, 0.18, 1.03) # Nombre de ré-échantillons (ou échantillons bootstrap) m <- 100000 # Risque alpha = 1 - confiance alpha <- 0.05 bootstrap(X, m, alpha) ``` ``` ## [1] "Intervalle de confiance : [0.53,1.41]" ``` :Suite Structure algorithmiques et fonctions {#algo, toggle=collapse, title-display=hidden} [Préambule : de l'algorithme au programme](caps_algo_1_de_l_algo_au_programme.html) [Les opérateurs](caps_algo_2_operateurs.html) [Structures algorithmiques](caps_algo_3_structures_algorithmiques.html) [Structures alternatives](caps_algo_4_structures_test.html) [Structures itératives](caps_algo_5_structures_iteratives.html) [Fonctions](caps_algo_6_fonctions.html) [Exercices bilan](caps_algo_7_exercice_bilan_programmation_R.html)
