6.2 Les boucles avec R

6.2.1 Les boucles

Les boucles dans R permettent de réaliser des tâches répétitives et ennuyeuses. Elles sont plus rapides que les tâches exécutées à la main et réduisent les occasions de commettre des erreurs. Par exemple, on peut utiliser une boucle pour effectuer une manipulation d’un jeu de données, telle que d’importer une série de fichiers stockés dans le même répertoire afin d’y extraire les données d’intérêt.

6.2.1.1 Structure

Il existe plusieurs façons de créer une boucle. Celle que nous verrons ici utilise la fonction for( ) pour mettre en place la boucle. Le contenu de la boucle sera inclus entre des accolades { }. Comme premier exemple, nous allons générer 100 échantillons contenant chacun 30 observations qui proviennent d’une distribution normale N(\(\mu = 10\), \(\sigma = 1\)), puis nous calculerons la moyenne de chaque échantillon.

Premièrement, on crée un objet pour stocker les moyennes calculées :

##objet pour stocker résultats
output <- rep(x = NA, times = 100)
output

##   [1] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
##  [46] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
##  [91] NA NA NA NA NA NA NA NA NA NA

La boucle comme telle sera initiée par la commande :

##on crée la boucle
for(i in 1:100) {
  ##on génère l'échantillon
  simdata <- rnorm(n = 30, mean = 10, sd = 1)
  ##on remplace le NA par la moyenne
  output[i] <- mean(simdata)
}

Décortiquons le code présenté ci-dessus:

La commande for( ) indique que la boucle sera effectuée de i = 1 jusqu’à i = 100. À noter que la boucle nécessite un nom de variable quelconque, par exemple, i, j, iter, sim, ou tout autre nom de notre choix. Toutefois, la boucle doit se référer à cette même variable plus tard. Immédiatement après la commande for( ), on trouvera une accolade { qui indique le début du corps de la boucle.
À la deuxième ligne, on crée un objet temporaire simdata qui contiendra 30 observations provenant d’une distribution normale ayant une moyenne de 10 et un écart-type de 1.
La prochaine ligne calcule la moyenne des observations de l’objet simdata et stocke cette moyenne dans l’objet output. Ainsi, le premier NA sera remplacé par la première moyenne.
La dernière ligne contient une accolade } qui indique la fin de la première itération (i = 1).
À l’itération i = 2, l’objet simdata stockera 30 nouvelles observations, et cette moyenne sera stockée dans output[i = 2]. La boucle continue jusqu’à l’itération i = 100}`.

Dans notre cas, la boucle a donné le résultat, qui sera différent sur votre machine puisque nous n’avons pas spécifié set.seed( ) :

output

##   [1]  9.843143  9.843070  9.905685 10.194206  9.792114  9.897694  9.962324  9.898893 10.186835  9.955964  9.932570 10.047830  9.575321
##  [14]  9.712321 10.054829 10.092733  9.989643  9.982827 10.043200 10.228745 10.197045 10.086565 10.249376 10.177044 10.298926  9.943895
##  [27] 10.046831  9.395522  9.990081  9.933483  9.918902 10.033769 10.107036 10.085682 10.126150 10.112908  9.953102  9.725974 10.070861
##  [40] 10.270107 10.063887 10.098440 10.061172 10.437588 10.166183 10.194932  9.820259 10.192904  9.777376 10.090552  9.763093  9.824848
##  [53] 10.127299 10.249598  9.790525 10.253496  9.876486  9.753068  9.975168  9.916710  9.890595  9.959431 10.028949 10.053114  9.968057
##  [66] 10.011333 10.006457  9.764378 10.126387 10.181128 10.020601  9.862828 10.048828  9.955702  9.872555 10.035752 10.225581  9.989991
##  [79] 10.107747 10.055392 10.146808  9.739223  9.797985  9.995382 10.134123 10.402918  9.991563  9.951723  9.947842 10.018451  9.989671
##  [92]  9.820996 10.023143 10.155746 10.350405 10.155261  9.895186  9.974955  9.903979  9.901210

Nous avons 100 moyennes calculées à partir de 100 échantillons de 30 observations provenant d’une population N(\(\mu = 10\), \(\sigma = 1\)).

On peut utiliser les boucles pour effectuer des manipulations plus complexes. D’ailleurs, les boucles sont utiles pour effectuer des tests de randomisation ou des simulations. Comme deuxième exemple, illustrons comment estimer la probabilité de commettre une erreur au niveau de l’expérience et au niveau des comparaisons telles que discutées dans le texte Analyse de variance à un critère. Dans l’exemple en question, nous avons considéré les probabilités de commettre des erreurs au niveau de l’expérience et au niveau des comparaisons avec trois groupes.

##illustration des erreurs au niveau de l'expérience et 
##au niveau des comparaisons
##3 groupes de n = 50 qui ne diffèrent pas les uns des autres
##tous provenant d'une même population normale avec mu=5 et sd=3
##alpha = 0.05
#############################
##on crée une matrice pour stocker les résultats
set.seed(seed = 221)
nsims <- 1000   #nombre de simulations
out <- matrix(NA, nrow = nsims, ncol = 4)
colnames(out) <- c("1 vs 2", "1 vs 3", "2 vs 3", "Erreur liée à l'expérience")
##débuter boucle
for (i in 1:nsims) {
  ##on crée 1er groupe
  ##tous les groupes sont égaux (H0 est vraie)
  groupe1 <- rnorm(n = 50, mean = 5, sd = 3)
  ##on crée 2ième groupe
  groupe2 <- rnorm(n = 50, mean = 5, sd = 3)
  ##on crée 3ième groupe
  groupe3 <- rnorm(n = 50, mean = 5, sd = 3)
  ##1vs2
  comp1 <- t.test(groupe1, groupe2)
  ##1vs3
  comp2 <- t.test(groupe1, groupe3)
  ##2vs3
  comp3 <- t.test(groupe2, groupe3)
  out[i, 1] <- ifelse(comp1$p.value<=0.05, "*", "NS")
  out[i, 2] <- ifelse(comp2$p.value<=0.05, "*", "NS")
  out[i, 3] <- ifelse(comp3$p.value<=0.05, "*", "NS")
  out[i, 4] <- ifelse(any(c(comp1$p.value, 
                            comp2$p.value,
                            comp3$p.value) <= 0.05), 1, 0)
}

Le code ci-dessus commence en spécifiant une valeur initiale afin d’alimenter l’algorithme du générateur de nombres aléatoires, ce qui permet de reproduire les mêmes résultats sur tous les ordinateurs. La prochaine ligne crée l’objet nsims pour spécifier le nombre d’itérations désirées pour la boucle (ici, 1000 itérations). La création d’un objet pour spécifier les itérations permet de se référer à cet objet tout le long du code. Ainsi, si on désire modifier le nombre d’itérations, nous n’aurons qu’à modifier l’objet nsims et le reste du code demeurera inchangé. La ligne suivante crée la matrice out de 4 colonnes et nsims rangées. À noter qu’ici, les résultats des trois comparaisons (groupe 1 vs groupe 2, groupe 1 vs groupe 3 et groupe 2 vs groupe 3) seront sauvegardés dans les trois premières colonnes alors que la dernière colonne sera réservée pour déterminer l’erreur au niveau de l’expérience.

La fonction colnames permet d’ajouter des étiquettes aux colonnes pour clairement identifier chaque colonne. La ligne débutant par for(i in 1:nsims) indique qu’on veut une boucle de i = 1 jusqu’à i = nsims. Ensuite, nous avons trois lignes de codes pour générer les données de chaque groupe conforme à \(H_0\), c’est-à-dire qu’il n’y a pas de différences entre les trois groupes (ils proviennent de la même population). Nous effectuons ensuite les comparaisons entre chaque groupe à l’aide d’un test \(t\). Pour chacun des tests, si la valeur de \(P \leq \alpha\), nous avons décidé d’ajouter un * dans la colonne, autrement nous avons inséré les caractères NS. À noter que si on rejette \(H_0\) pour un test \(t\) sur ces données, nous avons commis une erreur de type I puisque les données simulées sont réellement tirées de la même population (la boucle a été construite pour générer des données conformes à \(H_0\)).

Pour une itération i donnée, à chaque fois que la colonne 1 contient *, c’est une erreur au niveau de la comparaison. Il en va de même avec les colonnes 2 et 3. La quatrième colonne de la matrice out prend la valeur de 1 si au moins une des colonnes à l’itération i a rejeté \(H_0\) ce qui correspond à l’erreur au niveau de l’expérience, autrement la colonne prend la valeur de 0. Ces étapes sont répétées i fois.

Lorsque toutes les itérations de la boucle sont complétées, nous obtenons les résultats suivants :

##20 premières itérations
out[1:20, ]

##       1 vs 2 1 vs 3 2 vs 3 Erreur liée à l'expérience
##  [1,] "NS"   "NS"   "NS"   "0"                       
##  [2,] "NS"   "NS"   "NS"   "0"                       
##  [3,] "NS"   "NS"   "NS"   "0"                       
##  [4,] "NS"   "NS"   "NS"   "0"                       
##  [5,] "NS"   "NS"   "NS"   "0"                       
##  [6,] "NS"   "NS"   "NS"   "0"                       
##  [7,] "NS"   "NS"   "NS"   "0"                       
##  [8,] "NS"   "*"    "*"    "1"                       
##  [9,] "NS"   "NS"   "NS"   "0"                       
## [10,] "NS"   "*"    "NS"   "1"                       
## [11,] "NS"   "NS"   "NS"   "0"                       
## [12,] "NS"   "NS"   "NS"   "0"                       
## [13,] "NS"   "NS"   "NS"   "0"                       
## [14,] "NS"   "NS"   "NS"   "0"                       
## [15,] "NS"   "*"    "NS"   "1"                       
## [16,] "NS"   "NS"   "NS"   "0"                       
## [17,] "NS"   "NS"   "NS"   "0"                       
## [18,] "NS"   "NS"   "NS"   "0"                       
## [19,] "NS"   "NS"   "NS"   "0"                       
## [20,] "NS"   "NS"   "NS"   "0"

##erreurs au niveau de la comparaison 1 vs 2
err.comp1 <- sum(out[, 1] == "*")/nsims
err.comp2 <- sum(out[, 2] == "*")/nsims
err.comp3 <- sum(out[, 3] == "*")/nsims
err.comp1

## [1] 0.045

err.comp2

## [1] 0.05

err.comp3

## [1] 0.057

##erreur au niveau de l'expérience
err.exp <- sum(out[, 4] == "1")/nsims

Afin de déterminer la probabilité de commettre une erreur au niveau des comparaisons, on peut calculer la proportion du nombre de résultats significatifs (\(P < \alpha\)) par rapport au nombre total de simulations. On obtient les valeurs de 0.045, 0.05 et 0.057, pour les comparaions 1 vs 2, 1 vs 3 et 2 vs 3, respectivement. Ces valeurs oscillent très près du seuil de signification car elles dépendent directement de cette valeur (ici, 0.05). L’erreur au niveau de l’expérience est de 0.124.

Si on avait utilisé un seuil de 0.10, nous aurions obtenu les valeurs 0.109, 0.104 et 0.101 pour les comparaisons 1 vs 2, 1 vs 3 et 2 vs 3, respectivement. L’erreur au niveau de l’expérience dans ce cas serait de 0.245, obtenu en divisant le nombre de lignes avec au moins un résultat significatif par le nombre total d’itérations. Les calculs pour les erreurs au niveau de l’expérience et au niveau des comparaisons ne sont pas montrés ici pour un seuil \(\alpha = 0.10\). Nous vous suggérons de les faire comme exercice.

6.2.2 Alternatives aux boucles

L’utilisation de boucles dans R permet beaucoup de flexibilité pour réaliser des tâches très variées. Les boucles peuvent être imbriquées les unes dans les autres. Nous avons illustré l’utilisation de for( ), bien qu’il est aussi possible de réaliser des boucles avec des fonctions comme while( ) ou break( ). Pour de très gros jeux de données ou certains types de calculs, les boucles peuvent être lentes. Certaines fonctions peuvent aussi agir un peu comme des boucles. C’est le cas de la fonction apply( ) que nous verrons dans la section qui suit.

6.2.2.1 `apply( )` et fonctions apparentées

La fonction apply( ) et ses fonctions apparentées telles que lapply( ) et tapply( ) (pour plus d’options, faites help.search(apropos = "apply") permettent d’appliquer des calculs ou des manipulations à une matrice, à une liste ou à un jeu de données, entre autres. Les fonctions de la famille apply( ) utilisent un processus qu’on appelle la vectorisation, c’est-à-dire que les calculs sont effectués simultanément sur des vecteurs au lieu d’une valeur à la fois comme avec la boucle. Les fonctions apply( ) sont souvent beaucoup plus rapides que les boucles. On peut utiliser tapply( ) pour calculer une statistique (FUN) d’une variable numérique X = Ozone pour différents groupes définis par une variable catégorique INDEX = Jardin dans le jeu de données jardins.txt.

##importer le jeu de données jardins.txt
ozone <- read.table("Module_6/data/jardins.txt", header = TRUE)
str(ozone)

## 'data.frame':    45 obs. of  2 variables:
##  $ Ozone : num  14.82 20.59 6.18 7.61 31.35 ...
##  $ Jardin: chr  "A" "A" "A" "A" ...

##on calcule la moyenne de Ozone pour chacun des groupes
tapply(X = ozone$Ozone, INDEX = ozone$Jardin, FUN = mean)

##        A        B        C 
## 16.34200 18.20133 18.07800

##on calcule la SD de Ozone pour chacun des groupes
tapply(X = ozone$Ozone, INDEX = ozone$Jardin, FUN = sd)

##        A        B        C 
## 8.550944 4.506536 9.037485

##on calcule le n pour chacun des groupes
tapply(X = ozone$Ozone, INDEX = ozone$Jardin, FUN = length)

##  A  B  C 
## 15 15 15

##on détermine si chaque groupe contient des valeurs numériques
tapply(X = ozone$Ozone, INDEX = ozone$Jardin, FUN = is.matrix)

##     A     B     C 
## FALSE FALSE FALSE

La fonction tapply( ) est destinée aux jeux de données. De façon plus générale, on peut appliquer une fonction (FUN) sur les rangées ou colonnes (MARGIN) d’une matrice (X) à l’aide de apply( ) :

##on crée une matrice
mat.sim <- matrix(data = 1:24, nrow = 6, ncol = 4)
mat.sim

##      [,1] [,2] [,3] [,4]
## [1,]    1    7   13   19
## [2,]    2    8   14   20
## [3,]    3    9   15   21
## [4,]    4   10   16   22
## [5,]    5   11   17   23
## [6,]    6   12   18   24

##on calcule la somme de chaque rangée
apply(X = mat.sim, MARGIN = 1, FUN = sum)

## [1] 40 44 48 52 56 60

##on compare à rowSums - identique
rowSums(mat.sim)

## [1] 40 44 48 52 56 60

##on calcule le produit de chaque rangée
apply(X = mat.sim, MARGIN = 1, FUN = prod)

## [1]  1729  4480  8505 14080 21505 31104

##on calcule la somme de chaque colonne 
apply(X = mat.sim, MARGIN = 2, FUN = sum)

## [1]  21  57  93 129

##on compare à colSums - identique
colSums(mat.sim)

## [1]  21  57  93 129

##on calcule le produit de chaque colonne
apply(X = mat.sim, MARGIN = 2, FUN = prod)

## [1]      720   665280 13366080 96909120

On comprend que lorsque l’argument MARGIN prend la valeur de 1, les calculs sont effectués sur chacune des rangées. Si MARGIN a la valeur de 2, les calculs se réalisent sur chaque colonne.