Cours complet pour débutants pour apprendre la programmation en R

| Vincent Goulet | Laurent Caron |

Il s'agit d'un cours complet d'introduction à la programmation. En partant des concepts mathématiques, vous allez apprendre, de bout en bout, à programmer avec R, pour de meilleures analyses et modélisations. Plusieurs exercices vous sont proposés pour une prise en main pratique.

Un espace vous est proposé sur le forum pour recevoir vos avis. Commentez

- Titre : Cours complet pour débutants pour apprendre la programmation en R

- Auteurs : Vincent Goulet, Laurent Caron

- Parution : 10 mars 2018

- Public visé :

- Licence : Le contenu de cet article est rédigé par Vincent Goulet et est mis à disposition selon les termes de la Licence Creative Commons Attribution - Pas d'Utilisation Commerciale - Partage dans les Mêmes Conditions 3.0 non transposé. Les logos Developpez.com, en-tête, pied de page, css, et look & feel de l'article sont Copyright ® 2013 Developpez.com.

Sommaire

- Introduction

- 1. Présentation du langage R

- 1-1. Bref historique

- 1-2. Description sommaire de R

- 1-3. Interfaces

- 1-4. Stratégies de travail

- 1-5. Éditeurs de texte et environnements intégrés

- 1-6. Anatomie d'une session de travail

- 1-7. Répertoire de travail

- 1-8. Consulter l'aide en ligne

- 1-9. Où trouver de la documentation

- 1-10. Exemples

- 1-11. Exercices

- 2. Bases du langage R

- 2-1. Commandes R

- 2-2. Conventions pour les noms d'objets

- 2-3. Les objets R

- 2-3-1. Modes et types de données

- 2-3-2. Longueur

- 2-3-3. Objet spécial NULL

- 2-3-4. Valeurs manquantes, indéterminées et infinies

- 2-3-5. Attributs

- 2-4. Vecteurs

- 2-5. Matrices et tableaux

- 2-6. Listes

- 2-7. Data frames

- 2-8. Indiçage

- 2-9. Exemples

- 2-10. Exercices

- 3. Opérateurs et fonctions

- 3-1. Opérations arithmétiques

- 3-2. Opérateurs

- 3-3. Appels de fonctions

- 3-4. Quelques fonctions utiles

- 3-4-1. Manipulation de vecteurs

- 3-4-2. Recherche d'éléments dans un vecteur

- 3-4-3. Arrondi

- 3-4-4. Sommaires et statistiques descriptives

- 3-4-5. Sommaires cumulatifs et comparaisons élément par élément

- 3-4-6. Opérations sur les matrices

- 3-4-7. Produit extérieur

- 3-5. Structures de contrôle

- 3-5-1. Exécution conditionnelle

- 3-5-2. Boucles

- 3-6. Fonctions additionnelles

- 3-7. Exemples

- 3-8. Exercices

- 4. Exemples résolus

- 4-1. Calcul de valeurs actuelles

- 4-2. Fonctions de masse de probabilité

- 4-3. Fonction de répartition de la loi gamma

- 4-4. Algorithme du point fixe

- 4-5. Suite de Fibonacci

- 4-6. Exercices

- 5. Fonctions définies par l'usager

- 5-1. Définition d'une fonction

- 5-2. Retourner des résultats

- 5-3. Variables locales et globales

- 5-4. Exemple de fonction

- 5-5. Fonctions anonymes

- 5-6. Débogage de fonctions

- 5-7. Styles de codage

- 5-8. Exemples

- 5-9. Exercices

- 6. Concepts avancés

- 6-1. Argument '...'

- 6-2. Fonction apply

- 6-3. Fonctions lapply et sapply

- 6-4. Fonction mapply

- 6-5. Fonction replicate

- 6-6. Classes et fonctions génériques

- 6-7. Exemples

- 6-8. Exercices

- 7. Fonctions d'optimisation

- 7-1. Fonctions d'optimisation et de calcul de racines

- 7-1-1. Fonction uniroot

- 7-1-2. Fonction optimize

- 7-1-3. Fonction nlm

- 7-1-4. Fonction nlminb

- 7-1-5. Fonction optim

- 7-1-6. polyroot

- 7-2. Astuce Ripley

- 7-3. Pour en savoir plus

- 7-4. Exemples

- 7-5. Exercices

- 8. Générateurs de nombres aléatoires

- 8-1. Générateurs de nombres aléatoires

- 8-2. Fonctions de simulation de variables aléatoires non uniformes

- 8-3. Exemples

- 8-4. Exercices

- A. GNU Emacs et ESS : la base

- A-1. Mise en contexte

- A-2. Installation

- A-3. Description sommaire

- A-4. Emacs-ismes et Unix-ismes

- A-5. Commandes de base

- A-5-1. Les essentielles

- A-5-2. Manipulation de fichiers

- A-5-3. Déplacements simples du curseur

- A-5-4. Sélection de texte, copier, coller, couper

- A-5-5. Manipulation de fenêtres

- A-5-6. Manipulation de fihiers de script dans le mode ESS

- A-5-7. Interaction avec l'invite de commande R

- A-5-8. Consultation des rubriques d'aide de R

- A-6. Anatomie d'une session de travail (bis)

- A-7. Configuration de l'éditeur

- A-8. Aide et documentation

- B. RStudio : une introduction

- B-1. Installation

- B-2. Description sommaire

- B-3. Projets

- B-4. Commandes de base

- B-5. Anatomie d'une session de travail (ter)

- B-6. Configuration de l'éditeur

- B-7. Aide et documentation

- C. Planification d'une simulation en R

- C-1. Contexte

- C-2. Première approche : avec une boucle

- C-3. Seconde approche : avec sapply

- C-4. Variante de la seconde approche

- C-5. Gestion des fichiers

- C-6. Exécution en lot

- C-7. Conclusion

- D. Installation de packages dans R

- Réponses des exercices

- Chapitre 2

- Chapitre 3

- Chapitre 4

- Chapitre 5

- Chapitre 6

- Chapitre 7

- Chapitre 8

- Code source

- Bibliographie

- Note de la rédaction Developpez.com

- Index

Introduction

Le système R connaît depuis plus d'une décennie une progression remarquable dans ses fonctionnalités, dans la variété de ses domaines d'application ou, plus simplement, dans le nombre de ses utilisateurs. La documentation disponible suit la même tendance : plusieurs maisons d'édition proposent dans leur catalogue des ouvrages — voire des collections complètes — dédiés spécifiquement aux utilisations que l'on fait de R en sciences naturelles, en sciences sociales, en finance, etc. Néanmoins, peu d'ouvrages se concentrent sur l'apprentissage de R en tant que langage de programmation sous-jacent aux fonctionnalités statistiques. C'est la niche que nous tâchons d'occuper.

Le présent ouvrage se base sur des notes de cours et des exercices utilisés à l'École d'actuariat de l'Université Laval. L'enseignement du langage R est axé sur l'exposition à un maximum de code — que nous avons la prétention de croire bien écrit — et sur la pratique de la programmation. C'est pourquoi les chapitres sont rédigés de manière synthétique et qu'ils comportent peu d'exemples au fil du texte. En revanche, le lecteur est appelé à lire et à exécuter le code informatique se trouvant dans les sections d'exemples à la fin de chacun des chapitres. Ce code et les commentaires qui l'accompagnent reviennent sur l'essentiel des concepts du chapitre et les complémentent souvent. Nous considérons l'exercice d'« étude active » consistant à exécuter du code et à voir ses effets comme essentiel à l'apprentissage du langage R. Afin d'ancrer la présentation dans un contexte concret, plusieurs chapitres proposent également d'entrée de jeu un problème à résoudre. Nous fournissons des indications en cours de chapitre et la solution complète à la fin. Afin d'être facilement identifiables, ces éléments de contenu se présentent dans des encadrés de couleur contrastante.

Le texte des sections d'exemples est disponible en format électronique sous la rubrique de la documentation par des tiers (Contributed) du site Comprehensive R Archive Network. On peut obtenir directement l'archive en suivant le lien fourni à la page précédente.

Un symbole de lecture vidéo indique qu'une capsule vidéo est disponible dans la chaîne YouTube de l'auteur sur le sujet en hyperlien.

Certains exemples et exercices trahissent le premier public de ce document : on y fait à l'occasion référence à des concepts de base de la théorie des probabilités et des mathématiques financières. Les contextes actuariels demeurent néanmoins peu nombreux et ne devraient généralement pas dérouter le lecteur pour qui ces notions sont moins familières. Les réponses de tous les exercices se trouvent en annexe. En consultation électronique, le numéro d'un exercice est un hyperlien vers sa réponse, et vice versa.

On trouvera également en annexe de brèves introductions à l'éditeur de texte GNU Emacs et à l'environnement de développement intégré RStudio, un court exposé sur la planification d'une simulation en R, ainsi que des conseils sur l'administration d'une bibliothèque de packages R. Nous tenons à remercier M. Mathieu Boudreault pour sa collaboration dans la rédaction des exercices et Mme Mireille Côté pour la révision linguistique de la seconde édition.

Vincent Goulet
Québec, avril 2016

1 - Présentation du langage R

Objectifs du chapitre

- Comprendre ce qu'est un langage de programmation interprété.

- Connaître la provenance du langage R et les principes ayant guidé son développement.

- Mettre en place sur son poste de travail un environnement de développement en R.

- Démarrer une session R et exécuter des commandes simples.

- Utiliser des fichiers de script R de manière interactive.

- Créer, modifier et sauvegarder ses propres fichiers de script R.

1-1 - Bref historique

À l'origine fut le S, un langage pour « programmer avec des données » développé chez Bell Laboratories à partir du milieu des années 1970 par une équipe de chercheurs menée par John M. Chambers. Au fil du temps, le S a connu quatre principales versions communément identifiées par la couleur du livre dans lequel elles étaient présentées : version « originale » (Brown Book ; Becker et Chambers, 1984), version 2 (Blue Book ; Becker et collab., 1988), version 3 (White Book ; Chambers et Hastie, 1992) et version 4 (Green Book ; Chambers, 1998) ; voir aussi Chambers (2000) et Becker (1994) pour plus de détails.

Dès la fin des années 1980 et pendant près de vingt ans, le S a principalement été popularisé par une mise en œuvre commerciale nommée S-PLUS. En 2008, Lucent Technologies a vendu le langage S à Insightful Corporation, ce qui a effectivement stoppé le développement du langage par ses auteurs originaux. Aujourd'hui, le S est commercialisé de manière relativement confidentielle sous le nom Spotfire S+ par TIBCO Software.

Ce qui a fortement contribué à la perte d'influence de S-PLUS, c'est une nouvelle mise en œuvre du langage développée au milieu des années 1990. Inspirés à la fois par le S et par Scheme (un dérivé du Lisp), Ross Ihaka et Robert Gentleman proposent un langage pour l'analyse de données et les graphiques qu'ils nomment R (Ihaka et Gentleman, 1996). À la suggestion de Martin Maechler de l'ETH de Zurich, les auteurs décident d'intégrer leur nouveau langage au projet GNU(1), faisant de R un logiciel libre.

Ainsi disponible gratuitement et ouvert aux contributions de tous, R gagne rapidement en popularité là même où S-PLUS avait acquis ses lettres de noblesse, soit dans les milieux académiques. De simple dérivé « not unlike S », R devient un concurrent sérieux à S-PLUS, puis le surpasse lorsque les efforts de développement se rangent massivement derrière le projet libre. D'ailleurs John Chambers place aujourd'hui ses efforts de réflexion et de développement dans le projet R (Chambers, 2008).

1-2 - Description sommaire de R

R est un environnement intégré de manipulation de données, de calcul et de préparation de graphiques. Toutefois, ce n'est pas seulement un « autre » environnement statistique (comme SPSS ou SAS, par exemple), mais aussi un langage de programmation complet et autonome.

Tel que mentionné précédemment, le R est un langage principalement inspiré du S et de Scheme (Abelson et collab., 1996). Le S était à son tour inspiré de plusieurs langages, dont l'APL (autrefois un langage très prisé par les actuaires) et le Lisp. Comme tous ces langages, le R est interprété, c'est-à-dire qu'il requiert un autre programme — l'interprète — pour que ses commandes soient exécutées. Par opposition, les programmes de langages compilés, comme le C ou le C++, sont d'abord convertis en code machine par le compilateur puis directement exécutés par l'ordinateur.

Cela signifie donc que lorsque l'on programme en R, il n'est pas possible de plaider l'attente de la fin de la phase de compilation pour perdre son temps au travail. Désolé !

Le programme que''on lance lorsque l'on exécute R est en fait l'interprète. Celui-ci attend qu'on lui soumette des commandes dans le langage R, commandes qu'il exécutera immédiatement, une à une et en séquence. Par analogie, Excel est certes un logiciel de manipulation de données, de mise en forme et de préparation de graphiques, mais c'est aussi au sens large un langage de programmation interprété. On utilise le langage de programmation lorsque l'on entre des commandes dans une cellule d'une feuille de calcul. L'interprète exécute les commandes et affiche les résultats dans la cellule.

Le R est un langage particulièrement puissant pour les applications mathématiques et statistiques (et donc actuarielles) puisque précisément développé dans ce but. Parmi ses caractéristiques particulièrement intéressantes, on note :

- langage basé sur la notion de vecteur, ce qui simplifie les calculs mathématiques et réduit considérablement le recours aux structures itératives (boucles for, while, etc.) ;

- pas de typage ni de déclaration obligatoire des variables ;

- programmes courts, en général quelques lignes de code seulement ;

- temps de développement très court.

1-3 - Interfaces

R est d'abord et avant tout une application n'offrant qu'une invite de commande du type de celle présentée à la figure 1.1. En soi, cela n'est pas si différent d'un tableur tel que Excel : la zone d'entrée de texte dans une cellule n'est rien d'autre qu'une invite de commande(2), par ailleurs aux capacités d'édition plutôt réduites.

Fig. 1.1 - Fenêtre de la console sous OS X au démarrage de R

- Sous Windows, une interface graphique plutôt rudimentaire est disponible. Elle facilite certaines opérations telles que l'installation de packages externes, mais elle offre autrement peu de fonctionnalités additionnelles pour l'édition de code R.

- L'interface graphique de R sous Mac OS X est la plus élaborée. Outre la console présentée à la figure 1.1, l'application R.app comporte de nombreuses fonctionnalités, dont un éditeur de code assez complet.

- Sous Unix et Linux, R n'est accessible que depuis la ligne de commande du système d'exploitation (terminal). Aucune interface graphique n'est offerte avec la distribution de base de R.

Peu importe la plateforme utilisée — quoique dans une moindre mesure sous OS X — nous recommandons d'interagir avec R par le biais d'un éditeur de texte pour programmeur ou d'un environnement de développement intégré ; voir la section 1.5.

1-4 - Stratégies de travail

Dans la mesure où R se présente essentiellement sous forme d'une invite de commande, il existe deux grandes stratégies de travail avec cet environnement statistique.

1. On entre des expressions à la ligne de commande pour les évaluer immédiatement :

> 2 + 3
[1] 5

On peut également créer des objets contenant le résultat d'un calcul. Ces objets sont stockés en mémoire dans l'espace de travail de R :

> x <- exp(2)
> x
[1] 7.389056

Lorsque la session de travail est terminée, on sauvegarde une image de l'espace de travail sur le disque dur de l'ordinateur afin de pouvoir conserver les objets pour une future séance de travail :

> save.image()

Par défaut, l'image est sauvegardée dans un fichier nommé .RData dans le dossier de travail actif (voir la section 1.7) et cette image est automatiquement chargée en mémoire au prochain lancement de R, tel qu'indiqué à la fin du message d'accueil :

[Sauvegarde de la session précédente restaurée]

Cette approche, dite de « code virtuel et objets réels » a un gros inconvénient : le code utilisé pour créer les objets n'est pas sauvegardé entre les sessions de travail. Or, celui-ci est souvent bien plus compliqué que l'exemple ci-dessus. De plus, sans accès au code qui a servi à créer l'objet x, comment savoir ce que la valeur 7.389056 représente au juste ?

2. L'approche dite de « code réel et objets virtuels » considère que ce qu'il importe de conserver d'une session de travail à l'autre n'est pas tant les objets que le code qui a servi à les créer. Ainsi, on sauvegardera dans ce que l'on nommera des fichiers de script nos expressions R et le code de nos fonctions personnelles. Par convention, on donne aux fichiers de script un nom se terminant avec l'extension .R.

Avec cette approche, les objets sont créés au besoin en exécutant le code des fichiers de script. Comment ? Simplement en copiant le code du fichier de script et en le collant dans l'invite de commande de R. La figure 1.2 illustre schématiquement ce que le programmeur R a constamment sous les yeux : d'un côté son fichier de script et de l'autre, l'invite de commande R dans laquelle son code a été exécuté.

La méthode d'apprentissage préconisée dans cet ouvrage suppose que le lecteur utilisera cette seconde approche d'interaction avec R.

1-5 - Éditeurs de texte et environnements intégrés

Dans la mesure où l'on a recours à des fichiers de script comme expliqué à la section précédente, l'édition de code R est rendue beaucoup plus aisée avec un bon éditeur de texte pour programmeur ou un environnement de développement intégré (integrated development environment, IDE).

- Un éditeur de texte est différent d'un traitement de texte en ce qu'il s'agit d'un logiciel destiné à la création, l'édition et la sauvegarde de fichiers texte purs, c'est-à-dire dépourvus d'information de présentation et de mise en forme. Les applications Bloc-notes sous Windows ou TextEdit sous OS X sont deux exemples d'éditeurs de texte simples.

1. ## Fichier de script simple contenant des expressions R pour

2. ## faire des calculs et créer des objets.

3. 2 + 3

4. ## Probabilité d'une loi de Poisson(10)

5. x <- 7

6. 10^x * exp(-10) / factorial(x)

7. ## Petite fonction qui fait un calcul trivial

8. f <- function(x) x^2

9. ## Évaluation de la fonction

10. f(2)

1. R version 3.2.3 (2015-12-10) -- "Wooden Christmas-Tree"

2. Copyright (C) 2015 The R Foundation for Statistical Computing

3. Platform: x86_64-apple-darwin13.4.0 (64-bit)

4. [...]

5. > ## Fichier de script simple contenant des expressions R pour

6. > ## faire des calculs et créer des objets.

7. > 2 + 3

8. [1] 5

9. >

10. > ## Probabilité d'une loi de Poisson(10)

11. > x <- 7

12. > 10^x * exp(-10) / factorial(x)

13. [1] 0.09007923

14. >

15. > ## Petite fonction qui fait un calcul trivial

16. > f <- function(x) x^2

17. >

18. > ## Évaluation de la fonction

19. > f(2)

20. [1] 4

Fig. 1.2 - Fichier de script (en haut) et invite de commande R dans laquelle les expressions R ont été exécutées (en bas). Les lignes débutant par # dans le fichier de script sont des commentaires ignorés par l'interprète de commandes.

Toute analogie est boiteuse, mais celle-ci peut néanmoins s'avérer utile pour illustrer l'approche « code réel, objets virtuels ». Avec un tableur comme Excel, que choisiriez-vous si l'on vous proposait les options suivantes :

1. Ne sauvegarder que les valeurs calculées d'une feuille de calcul ;

2. Ne sauvegarder que les formules d'une feuille de calcul, quitte à les réévaluer par la suite pour obtenir les valeurs correspondantes.

La seconde, sûrement… Or, sauvegarder l'espace de travail R dans un fichier .RData (l'approche « objets réels, code virtuel ») correspond à la première option, alors que sauvegarder le code source des fichiers de script correspond à la seconde.

Là où l'analogie est boiteuse, la raison pour laquelle on n'a jamais à se poser la question avec un tableur, c'est que le progiciel évalue constamment et en temps réel les formules pour afficher les résultats dans les cellules.

- Un éditeur de texte pour programmeur saura en plus reconnaître la syntaxe d'un langage de programmation et assister à sa mise en forme : indentation automatique du code, marquage des mots-clés, manipulation d'objets, etc.

- Un éditeur compatible avec R réduira, entre autres, l'opération de copier-coller à un simple raccourci clavier.

Le lecteur peut utiliser le logiciel de son choix pour l'édition de code R. Certains offrent simplement plus de fonctionnalités que d'autres.

- Dans la catégorie des éditeurs de texte, nous recommandons le vénérable et très puissant éditeur pour programmeur GNU Emacs. À la question 6.2 de la foire aux questions de R (Hornik, 2013), « Devrais-je utiliser R à l'intérieur de Emacs ? », la réponse est : « Oui, absolument ».

En effet, combiné avec le mode ESS (Emacs Speaks Statistics), Emacs offre un environnement de développement aussi riche qu'efficace. Entre autres fonctionnalités uniques à cet éditeur, le fichier de script et l'invite de commandes R sont regroupés dans la même fenêtre, comme on peut le voir à la figure 1.3.

Fig. 1.3 - Fenêtre de GNU Emacs sous OS X en mode d'édition de code R. Dans la partie du haut, on retrouve le fichier de script de la figure 1.2 et dans la partie du bas, l'invite de commandes R.

Emblème du logiciel libre, Emacs est disponible gratuitement et à l'identique sur toutes les plateformes supportées par R, dont Windows, OS X et Linux.

Consulter l'annexe A pour en savoir plus sur GNU Emacs et apprendre les commandes essentielles pour y faire ses premiers pas.

- Emacs est toutefois un logiciel difficile à apprivoiser, surtout pour les personnes moins à l'aise avec l'informatique.

Fig. 1.4 - Fenêtre de RStudio sous OS X dans sa configuration par défaut. Dans le sens des aiguilles d'une montre, on y retrouve : le fichier de script de la figure 1.2 ; la liste des objets de l'environnement (vide, ici) ; une interface pour charger des packages ; l'invite de commandes R.

- RStudio est un environnement de développement intégré (IDE) créé spécifiquement pour travailler avec R. Sa popularité connaît une progression foudroyante depuis 2014. Il permet de consulter dans une interface conviviale ses fichiers de script, la ligne de commande R, les rubriques d'aide, les graphiques, etc. ; voir la figure 1.4.
RStudio est disponible sur les plateformes Windows, OS X et Linux.

- Il existe plusieurs autres options pour éditer efficacement du code R — et le Bloc-notes de Windows n'en fait pas partie ! Nous recommandons plutôt :

- sous Windows, l'éditeur Notepad++ additionné de l'extension NppToR (Redd, 2010), tous deux des logiciels libres, ou le partagiciel WinEdt muni de l'extension libre R-WinEdt (Ligges, 2003) ;

- sous OS X, tout simplement l'éditeur de texte très complet intégré à l'application R.app, ou alors l'éditeur de texte commercial TextMate (essai gratuit de 30 jours) ;

- sous Linux, Vim et Kate semblent les choix les plus populaires dans la communauté R, après Emacs et RStudio.

1-6 - Anatomie d'une session de travail

Dans ses grandes lignes, une session de travail avec R réunit les étapes ci-dessous.

1. Démarrer une session R en cliquant sur l'icône de l'application si l'on utilise une interface graphique ou RStudio, ou alors en suivant la procédure appropriée pour travailler avec son éditeur de texte.

2. Ouvrir un fichier de script existant ou en créer un nouveau à l'aide de l'éditeur de texte de son choix.

3. Bâtir graduellement un fichier de script en y consignant le code R que l'on souhaite sauvegarder et les commentaires qui permettront de s'y retrouver plus tard. Tester les commandes à la ligne de commande. Au cours de la phase de développement, on fera généralement de nombreux allers-retours entre la ligne de commande et le fichier de script.

4. Sauvegarder son fichier de script et quitter l'éditeur ou l'environnement de développement.

5. Si nécessaire — et c'est rarement le cas — sauvegarder l'espace de travail de la session R avec save.image(). En fait, on ne voudra sauvegarder nos objets R que lorsque ceux-ci sont très longs à créer, par exemple, les résultats d'une simulation.

6. Quitter R en tapant q() à la ligne de commande ou en fermant l'interface graphique par la procédure usuelle. Encore ici, la manière de procéder est quelque peu différente dans GNU Emacs ; voir l'annexe A.

Les étapes 1 et 2 sont interchangeables, tout comme les étapes 4, 5 et 6.

L'annexe A explique plus en détail la procédure à suivre lorsque l'on utilise GNU Emacs. L'annexe B fait de même lorsque l'on utilise l'environnement intégré RStudio.

1-7 - Répertoire de travail

Le répertoire de travail (workspace) de R est le dossier par défaut dans lequel le logiciel : 1) va rechercher des fichiers de script ou de données ; et 2) va sauvegarder l'espace de travail dans le fichier .RData. Le répertoire de travail est déterminé au lancement de R.

- Les interfaces graphiques de R démarrent avec un répertoire de travail par défaut. Pour le changer, utiliser l'entrée appropriée dans le menu Fichier (Windows) ou Divers (OS X). Consulter aussi les foires aux questions spécifiques aux interfaces graphiques (Ripley et Murdoch, 2013 ; Iacus et collab., 2013) pour des détails additionnels sur la gestion des répertoires de travail.

- Dans RStudio, on change le répertoire de travail via le menu Session.

- Avec GNU Emacs, la situation est un peu plus simple puisque l'on doit spécifier un répertoire de travail chaque fois que l'on démarre un processus R ; voir l'annexe A.

1-8 - Consulter l'aide en ligne

Les rubriques d'aide des diverses fonctions de R contiennent une foule d'informations ainsi que des exemples d'utilisation. Leur consultation est tout à fait essentielle.

- Pour consulter la rubrique d'aide de la fonction foo, on peut entrer à la ligne de commande

> ?foo

ou

> help(foo)

1-9 - Où trouver de la documentation

La documentation officielle de R se compose de six guides accessibles depuis le menu Aide des interfaces graphiques ou encore en ligne dans le site du projet R(3). Pour le débutant, seuls An Introduction to R et, possiblement, R Data Import/Export peuvent s'avérer des ressources utiles à court terme.

Plusieurs livres — en versions papier ou électronique, gratuits ou non — ont été publiés sur R. On en trouvera une liste exhaustive dans la section Documentation du site du projet R.

Depuis plusieurs années maintenant, les ouvrages de Venables et Ripley (2000, 2002) demeurent des références standards de facto sur les langages S et R. Plus récent, Braun et Murdoch (2007) participe du même effort que le présent ouvrage en se concentrant sur la programmation en R plutôt que sur ses applications statistiques.

1-10 - Exemples

1. ### Générer deux vecteurs de nombres pseudo-aléatoires issus

2. ### d'une loi normale centrée réduite.

3. x <- rnorm(50)

4. y <- rnorm(x)

6. ### Graphique des couples (x, y).

7. plot(x, y)

9. ### Graphique d'une approximation de la densité du vecteur x.

10. plot(density(x))

12. ### Générer la suite 1, 2, ... , 10.

13. 1:10

15. ### La fonction 'seq' sert à générer des suites plus générales.

16. seq(from = -5, to = 10, by = 3)

17. seq(from = -5, length = 10)

19. ### La fonction 'rep' sert à répéter des valeurs.

20. rep(1, 5)          # répéter 1 cinq fois

21. rep(1:5, 5)        # répéter le vecteur 1, ... ,5 cinq fois

22. rep(1:5, each = 5) # répéter chaque élément du vecteur cinq fois

24. ### Arithmétique vectorielle.

25. v <- 1:12       # initialisation d'un vecteur

26. v + 2              # additionner 2 à chaque élément de v

27. v * -12:-1         # produit élément par élément

28. v + 1:3            # le vecteur le plus court est recyclé

30. ### Vecteur de nombres uniformes sur l'intervalle [1, 10].

31. v <- runif(12, min = 1, max = 10)

32. v

34. ### Pour afficher le résultat d'une affectation, placer la

35. ### commande entre parenthèses.

36. ( v <- runif(12, min = 1, max = 10) )

38. ### Arrondi des valeurs de v à l'entier près.

39. ( v <- round(v) )

41. ### Créer une matrice 3 x 4 à partir des valeurs de

42. ### v. Remarquer que la matrice est remplie par colonne.

43. ( m <- matrix(v, nrow = 3, ncol = 4) )

45. ### Les opérateurs arithmétiques de base s'appliquent aux

46. ### matrices comme aux vecteurs.

47. m + 2

48. m * 3

49. m ^ 2

51. ### Éliminer la quatrième colonne afin d'obtenir une matrice

52. ### carrée.

53. ( m <- m[,-4] )

55. ### Transposée et inverse de la matrice m.

56. t(m)

57. solve(m)

59. ### Produit matriciel.

60. m %*% m               # produit de m avec elle-même

61. m %*% solve(m)        # produit de m avec son inverse

62. round(m %*% solve(m)) # l'arrondi donne la matrice identité

64. ### Consulter la rubrique d'aide de la fonction 'solve'.

65. ?solve

67. ### Liste des objets dans l'espace de travail.

68. ls()

70. ### Nettoyage.

71. rm(x, y, v, m)

1-11 - Exercices

1.1 Démarrer une session R et entrer une à une les expressions ci-dessous à la ligne de commande. Observer les résultats.

1. > ls()

2. > pi

3. > (v <- c(1, 5, 8))

4. > v * 2

5. > x <- v + c(2, 1, 7)

6. > x

7. > ls()

8. > q()

1.2 Ouvrir dans un éditeur de texte le fichier de script contenant le code de la section précédente. Exécuter le code ligne par ligne et observer les résultats. Répéter l'exercice avec un ou deux autres éditeurs de texte afin de les comparer et de vous permettre d'en choisir un pour la suite.

1.3 Consulter les rubriques d'aide d'une ou plusieurs des fonctions rencontrées lors de l'exercice précédent. Observer d'abord comment les rubriques d'aide sont structurées — elles sont toutes identiques — puis exécuter quelques expressions tirées des sections d'exemples.

1.4 Exécuter le code de l'exemple de session de travail R que l'on trouve à l'annexe A de Venables et collab. (2013). En plus d'aider à se familiariser avec R, cet exercice permet de découvrir les fonctionnalités du logiciel en tant qu'outil statistique.

2 - Bases du langage R

Objectifs du chapitre

- Écrire et interpréter la syntaxe et la sémantique du langage R.

- Identifier les principaux types de données disponibles dans R.

- Utiliser les divers modes d'objets (en particulier numeric, character et logical) et la conversion automatique de l'un à l'autre.

- Créer et manipuler des vecteurs, matrices, tableaux, listes et data frames.

- Extraire des données d'un objet ou y affecter de nouvelles valeurs à l'aide des diverses méthodes d'indiçage.

Pour utiliser un langage de programmation, il faut en connaître la syntaxe et la sémantique, du moins dans leurs grandes lignes. C'est dans cet esprit que ce chapitre introduit des notions de base du langage R telles que l'expression, l'affectation et l'objet. Le concept de vecteur se trouvant au cœur du langage, le chapitre fait une large place à la création et à la manipulation des vecteurs et autres types d'objets de stockage couramment employés en programmation en R.

Énoncé du problème

Une ligue de hockey compte huit équipes. Le classement de la ligue est disponible quotidiennement dans le journal dans le format habituel ; voir le tableau 2.1.

Équipe

MJ

V

D

DP

PTS

Washington

55

36

16

3

75

Dallas

56

32

19

5

69

Chicago

57

30

21

6

66

Los Angeles

58

30

22

6

66

St-Louis

56

25

19

12

62

Détroit

57

25

21

11

61

Montréal

56

22

27

7

51

Boston

57

24

31

2

50

Tableau 2.1 : Classement de la ligue de hockey pour le problème à résoudre du chapitre. MJ - matchs joués ; V - victoires ; D - défaites ; DP : défaites en prolongation ; PTS - points.

Afin d'effectuer différentes analyses statistiques, on désire intégrer ces données dans un espace de travail R. On doit donc déterminer le type d'objet R approprié pour stocker le classement de la ligue. Ensuite, on souhaite extraire les valeurs suivantes de l'objet précédemment créé.

1. Le nombre d'équipes de la ligue.

2. La fiche complète de l'équipe de Montréal.

3. La fiche complète de l'équipe à la septième position du classement.

2-1 - Commandes R

Comme vu au chapitre précédent, l'utilisateur de R interagit avec l'interprète R en entrant des commandes à l'invite de commande. Toute commande R est soit une expression, soit une affectation.

- Normalement, une expression est immédiatement évaluée et le résultat est affiché à l'écran :

1. > 2 + 3

2. [1] 5

3. > pi

4. [1] 3.141593

5. > cos(pi/4)

6. [1] 0.7071068

Dans les anciennes versions de S et R, on pouvait affecter avec le caractère de soulignement « _ ». Cet emploi n'est plus permis, mais la pratique subsiste dans le mode ESS de Emacs. Ainsi, taper le caractère « _ » hors d'une chaîne de caractères dans Emacs génère automatiquement ␣<-␣. Si on souhaite véritablement obtenir le caractère de soulignement, il suffit d'appuyer deux fois successives sur « _ ».

- Lors d'une affectation, une expression est évaluée, mais le résultat est stocké dans un objet (variable) et rien n'est affiché à l'écran. Le symbole d'affectation est « <- », c'est-à-dire les deux caractères « < » et « - » placés obligatoirement l'un à la suite de l'autre :

1. > a <- 5

2. > a

3. [1] 5

4. > b <- a

5. > b

6. [1] 5

- Pour affecter le résultat d'un calcul dans un objet et simultanément afficher ce résultat, il suffit de placer l'affectation entre parenthèses pour ainsi créer une nouvelle expression(4) :

> (a <- 2 + 3)
[1] 5

- Le symbole d'affectation inversé « -> » existe aussi, mais il est rarement utilisé.

- Éviter d'utiliser l'opérateur « = » pour affecter une valeur à une variable puisque cette pratique est susceptible d'engendrer de la confusion avec les constructions nom = valeur dans les appels de fonction.

Que ce soit dans les fichiers de script ou à la ligne de commande, on sépare les commandes R les unes des autres par un point-virgule « ; » ou par un retour à la ligne.

- On considère généralement comme du mauvais style d'employer les deux, c'est-à-dire de placer des points-virgules à la fin de chaque ligne de code, surtout dans les fichiers de script.

- Le point-virgule peut être utile pour séparer deux courtes expressions ou plus sur une même ligne :

> a <- 5; a + 2
[1] 7

C'est le seul emploi du point-virgule que l'on rencontrera dans cet ouvrage.

On peut regrouper plusieurs commandes en une seule expression en les entourant d'accolades { }.

- Le résultat du regroupement est la valeur de la dernière commande :

1. > {

2. +      a <- 2 + 3

3. +      b <- a

4. +      b

5. + }

6. [1] 5

- Par conséquent, si le regroupement se termine par une assignation, aucune valeur n'est retournée ni affichée à l'écran :

1. > {

2. +      a <- 2 + 3

3. +      b <- a

4. + }

- Les règles ci-dessus joueront un rôle important dans la composition de fonctions ; voir le chapitre 5.

- Comme on peut le voir ci-dessus, lorsqu'une commande n'est pas complète à la fin de la ligne, l'invite de commande de R change de >␣ à +␣ pour nous inciter à compléter notre commande.

2-2 - Conventions pour les noms d'objets

Les caractères permis pour les noms d'objets sont les lettres minuscules a-z et majuscules A-Z, les chiffres 0-9, le point « . » et le caractère de soulignement « _ ». Selon l'environnement linguistique de l'ordinateur, il peut être permis d'utiliser des lettres accentuées, mais cette pratique est fortement découragée puisqu'elle risque de nuire à la portabilité du code.

- Les noms d'objets ne peuvent commencer par un chiffre. S'ils commencent par un point, le second caractère ne peut être un chiffre.

- Le R est sensible à la casse, ce qui signifie que foo, Foo et FOO sont trois objets distincts. Un moyen simple d'éviter des erreurs liées à la casse consiste à n'employer que des lettres minuscules.

- Certains noms sont utilisés par le système R, aussi vaut-il mieux éviter de les utiliser. En particulier, éviter d'utiliser

c, q, t, C, D, I, diff, length, mean, pi , range, var

- Certains mots sont réservés et il est interdit de les utiliser comme nom d'objet. Les mots réservés pour le système sont :

1. break, else, for, function, if, in, next, repeat, return, while,

2. TRUE, FALSE,

3. Inf, NA, NaN, NULL,

4. NA_integer_, NA_real_, NA_complex_, NA_character_,

5. ... , ..1, ..2, etc.

Oui, ... (point-point-point) est véritablement un nom d'objet dans R ! Son usage est expliqué à la section 6.1.

- Les variables T et F prennent par défaut les valeurs TRUE et FALSE, respectivement, mais peuvent être réaffectées :

1. > T

2. [1] TRUE

3. > F

4. [1] FALSE

5. > TRUE <- 3

6. Error in TRUE <- 3 : membre gauche de l′assignation

7. (do_set) incorrect

8. > (T <- 3)

9. [1] 3

- Nous recommandons de toujours écrire les valeurs booléennes TRUE et FALSE au long pour éviter des bogues difficiles à détecter.

Tab. 2.2 - Modes disponibles et contenus correspondants

Mode

Contenu de l'objet

numeric

nombres réels

complex

nombres complexes

logical

valeurs booléennes (vrai/faux)

character

chaînes de caractères

function

fonction

list

données quelconques

expression

expressions non évaluées

2-3 - Les objets R

Tout dans le langage R est un objet : les variables contenant des données, les fonctions, les opérateurs, même le symbole représentant le nom d'un objet est lui-même un objet. Les objets possèdent au minimum un mode et une longueur et certains peuvent être dotés d'un ou plusieurs attributs.

- Le mode d'un objet est obtenu avec la fonction mode :

> v <- c(1, 2, 5, 9)
> mode(v)
[1] "numeric"

- La longueur d'un objet est obtenue avec la fonction length :

> length(v)
[1] 4

2-3-1 - Modes et types de données

Le mode prescrit ce qu'un objet peut contenir. À ce titre, un objet ne peut avoir qu'un seul mode. Le tableau 2.2 contient la liste des principaux modes disponibles en R. À chacun de ces modes correspond une fonction du même nom servant à créer un objet de ce mode.

- Les objets de mode « numeric », « complex », « logical » et « character » sont des objets simples (atomic en anglais) qui ne peuvent contenir que des données d'un seul type.

- En revanche, les objets de mode « list » ou « expression » sont des objets récursifs qui peuvent contenir d'autres objets. Par exemple, une liste peut contenir une ou plusieurs autres listes ; voir la section 2.6 pour plus de détails.

- La fonction typeof permet d'obtenir une description plus précise de la représentation interne d'un objet (c'est-à-dire au niveau de la mise en œuvre en C). Le mode et le type d'un objet sont souvent identiques.

2-3-2 - Longueur

La longueur d'un objet est égale au nombre d'éléments qu'il contient.

- La longueur, au sens R du terme, d'une chaîne de caractères est toujours 1. Un objet de mode character doit contenir plusieurs chaînes de caractères pour que sa longueur soit supérieure à 1 :

1. > v1 <- "actuariat"

2. > length(v1)

3. [1] 1

4. > v2 <- c("a", "c", "t", "u", "a", "r", "i", "a", "t")

5. > length(v2)

6. [1] 9

- Il faut utiliser la fonction nchar pour obtenir le nombre de caractères dans une chaîne :

1. > nchar(v1)

2. [1] 9

3. > nchar(v2)

4. [1] 1 1 1 1 1 1 1 1 1

- Un objet peut être de longueur 0 et doit alors être interprété comme un contenant qui existe, mais qui est vide :

1. > v <- numeric(0)

2. > length(v)

3. [1] 0

2-3-3 - Objet spécial NULL

L'objet spécial NULL représente « rien », ou le vide.

- Son mode est NULL.

- Sa longueur est 0.

- Toutefois différent d'un objet vide :

- un objet de longueur 0 est un contenant vide ;

- NULL est « pas de contenant ».

- La fonction is.null teste si un objet est NULL ou non.

2-3-4 - Valeurs manquantes, indéterminées et infinies

Dans les applications statistiques, il est souvent utile de pouvoir représenter des données manquantes. Dans R, l'objet spécial NA remplit ce rôle.

- Par défaut, le mode de NA est logical, mais NA ne peut être considéré ni comme TRUE, ni comme FALSE.

- Toute opération impliquant une donnée NA a comme résultat NA.

- Certaines fonctions (sum, mean, par exemple) ont par conséquent un argument na.rm qui, lorsque TRUE, élimine les données manquantes avant de faire un calcul.

- La valeur NA n'est égale à aucune autre, pas même elle-même (selon la règle ci-dessus, le résultat de la comparaison est NA) :

> NA == NA
[1] NA

- Par conséquent, pour tester si les éléments d'un objet sont NA ou non il faut utiliser la fonction is.na :

> is.na(NA)
[1] TRUE

La norme IEEE 754 régissant la représentation interne des nombres dans un ordinateur (IEEE, 2003) prévoit les valeurs mathématiques spéciales +∞ et −∞ ainsi que les formes indéterminées du type ou . R dispose d'objets spéciaux pour représenter ces valeurs.

- Inf Représente .

- -Inf représente .

- NaN (Not a Number) représente une forme indéterminée.

- Ces valeurs sont testées avec les fonctions is.infinite, is.finite et is.nan.

Tab. 2.3 - Attributs les plus usuels d'un objet

Attribut

Utilisation

class

affecte le comportement d'un objet

dim

dimensions des matrices et tableaux

dimnames

étiquettes des dimensions des matrices et tableaux

names

étiquettes des éléments d'un objet

2-3-5 - Attributs

Les attributs d'un objet sont des éléments d'information additionnels liés à cet objet. La liste des attributs les plus fréquemment rencontrés se trouve au tableau 2.3. Pour chaque attribut, il existe une fonction du même nom servant à extraire l'attribut correspondant d'un objet.

- Plus généralement, la fonction attributes permet d'extraire ou de modifier la liste des attributs d'un objet. On peut aussi travailler sur un seul attribut à la fois avec la fonction attr.

- On peut ajouter à peu près ce que l'on veut à la liste des attributs d'un objet. Par exemple, on pourrait vouloir attacher au résultat d'un calcul la méthode de calcul utilisée :

1. > x <- 3

2. > attr(x, "methode") <- "au pif"

3. > attributes(x)

4. $methode

5. [1] "au pif"

- Extraire un attribut qui n'existe pas retourne NULL :

> dim(x)
NULL

- À l'inverse, donner à un attribut la valeur NULL efface cet attribut :

1. > attr(x, "methode") <- NULL

2. > attributes(x)

3. NULL

2-4 - Vecteurs

En R, à toutes fins pratiques, tout est un vecteur. Contrairement à certains autres langages de programmation, il n'y a pas de notion de scalaire en R ; un scalaire est simplement un vecteur de longueur 1. Comme nous le verrons au chapitre 3, le vecteur est l'unité de base dans les calculs.

- Dans un vecteur simple, tous les éléments doivent être du même mode. Nous nous restreignons à ce type de vecteurs pour le moment.

- Les fonctions de base pour créer des vecteurs sont :

- c (concaténation) ;

- numeric (vecteur de mode numeric) ;

- logical (vecteur de mode logical) ;

- character (vecteur de mode character).

- Il est possible (et souvent souhaitable) de donner une étiquette à chacun des éléments d'un vecteur.

1. > (v <- c(a = 1, b = 2, c = 5))

2. a b c

3. 1 2 5

4. > v <- c(1, 2, 5)

5. > names(v) <- c("a", "b", "c")

6. > v

7. a b c

8. 1 2 5

Ces étiquettes font alors partie des attributs du vecteur.

- L'indiçage dans un vecteur se fait avec les crochets [ ]. On peut extraire un élément d'un vecteur par sa position ou par son étiquette, si elle existe (auquel cas cette approche est beaucoup plus sûre).

1. > v[3]

2. c5

3. > v["c"]

4. c5

La section 2.8 traite plus en détail de l'indiçage des vecteurs et des matrices.

Dans un vecteur simple, tous les éléments doivent être du même mode. Or, les informations du classement de la ligue comportent à la fois des chaînes de caractères et des nombres. De plus, le classement se présente sous forme d'un tableau à deux dimensions, alors qu'un vecteur n'en compte qu'une seule. Le vecteur simple n'est donc pas le type d'objet approprié pour stocker le classement de la ligue.

2-5 - Matrices et tableaux

Le R étant un langage spécialisé pour les calculs mathématiques, il supporte tout naturellement et de manière intuitive — à une exception près, comme nous le verrons — les matrices et, plus généralement, les tableaux à plusieurs dimensions.

Les matrices et tableaux ne sont rien d'autre que des vecteurs dotés d'un attribut dim. Ces objets sont donc stockés, et peuvent être manipulés, exactement comme des vecteurs simples.

- Une matrice est un vecteur avec un attribut dim de longueur 2. Cela change implicitement la classe de l'objet pour « matrix » et, de ce fait, le mode d'affichage de l'objet ainsi que son interaction avec plusieurs opérateurs et fonctions.

- La fonction de base pour créer des matrices est matrix :

1. > matrix(1:6, nrow = 2, ncol = 3)

2. [,1] [,2] [,3]

3. [1,]    1    3    5

4. [2,]    2    4    6

5. > matrix(1:6, nrow = 2, ncol = 3, byrow = TRUE)

6. [,1] [,2] [,3]

7. [1,]    1    2    3

8. [2,]    4    5    6

- La généralisation d'une matrice à plus de deux dimensions est un tableau (array). Le nombre de dimensions du tableau est toujours égal à la longueur de l'attribut dim. La classe implicite d'un tableau est « array ».

- La fonction de base pour créer des tableaux est array :

1. > array(1:24, dim = c(3, 4, 2))

2. , , 1

3. [,1] [,2] [,3] [,4]

4. [1,]    1    4    7   10

5. [2,]    2    5    8   11

6. [3,]    3    6    9   12

7. , , 2

8. [,1] [,2] [,3] [,4]

9. [1,]   13   16   19   22

10. [2,]   14   17   20   23

11. [3,]   15   18   21   24

On remarquera ci-dessus que les matrices et tableaux sont remplis en faisant d'abord varier la première dimension, puis la seconde, etc. Pour les matrices, cela revient à remplir par colonne. On conviendra que cette convention, héritée du Fortran, n'est pas des plus intuitives.

La fonction matrix a un argument byrow qui permet d'inverser l'ordre de remplissage. Cependant, il vaut mieux s'habituer à la convention de R que d'essayer constamment de la contourner.

L'ordre de remplissage inhabituel des tableaux rend leur manipulation difficile si on ne les visualise pas correctement. Imaginons un tableau de dimensions 3 × 4 × 5.

- Il faut voir le tableau comme cinq matrices 3 × 4 (remplies par colonne !) les unes derrière les autres.

- Autrement dit, le tableau est un prisme rectangulaire haut de 3 unités, large de 4 et profond de 5.

- Si l'on ajoute une quatrième dimension, cela revient à aligner des prismes les uns derrière les autres, et ainsi de suite.

La figure 2.1 fournit une représentation schématique des tableaux à trois et quatre dimensions.

Comme pour les vecteurs, l'indiçage des matrices et tableaux se fait avec les crochets [ ].

Représentation schématique de tableaux. Les chiffres encerclés identifient l'ordre de remplissage.

- On extrait un élément d'une matrice en précisant sa position dans chaque dimension de celle-ci, séparées par des virgules :

1. > (m <- matrix(c(40, 80, 45, 21, 55, 32),

2. + nrow = 2, ncol = 3))

3. [,1] [,2] [,3]

4. [1,]   40   45   55

5. [2,]   80   21   32

6. > m[1, 2]

7. [1] 45

- On peut aussi ne donner que la position de l'élément dans le vecteur sous-jacent :

> m[3]
[1] 45

- Lorsqu'une dimension est omise dans les crochets, tous les éléments de cette dimension sont extraits :

> m[2, ]
[1] 80 21 32

- Les idées sont les mêmes pour les tableaux.

- Pour le reste, les règles d'indiçage de vecteurs exposées à la section 2.8 s'appliquent à chaque position de l'indice d'une matrice ou d'un tableau.

Des fonctions permettent de fusionner des matrices et des tableaux ayant au moins une dimension identique.

- La fonction rbind permet de fusionner verticalement deux matrices (ou plus) ayant le même nombre de colonnes.

1. > n <- matrix(1:9, nrow = 3)

2. > rbind(m, n)

3. [,1] [,2] [,3]

4. [1,]   40   45   55

5. [2,]   80   21   32

6. [3,]    1    4    7

7. [4,]    2    5    8

8. [5,]    3    6    9

- La fonction cbind permet de fusionner horizontalement deux matrices (ou plus) ayant le même nombre de lignes.

1. > n <- matrix(1:4, nrow = 2)

2. > cbind(m, n)

3. [,1] [,2] [,3] [,4] [,5]

4. [1,]   40   45   55    1    3

5. [2,]   80   21   32    2    4

Une matrice convient bien pour stocker un tableau de données. Toutefois, puisque la matrice est en fait un vecteur avec un attribut dim de longueur 2, tous les éléments doivent être du même mode, comme c'était le cas avec les vecteurs simples. Impossible dans ce cas d'y stocker le nom des équipes. La matrice n'est toujours pas le type d'objet approprié.

2-6 - Listes

La liste est le mode de stockage le plus général et polyvalent du langage R. Il s'agit d'un type de vecteur spécial dont les éléments peuvent être de n'importe quel mode, y compris le mode list. Cela permet donc d'emboîter des listes, d'où le qualificatif de récursif pour ce type d'objet.

- La fonction de base pour créer des listes est list :

1. > (x <- list(size = c(1, 5, 2), user = "Joe", new = TRUE))

2. $size

3. [1] 1 5 2

4. $user

5. [1] "Joe"

6. $new

7. [1] TRUE

Ci-dessus, le premier élément de la liste est de mode « numeric », le second de mode « character » et le troisième de mode « logical ».

- Nous recommandons de nommer les éléments d'une liste. En effet, les listes contiennent souvent des données de types différents et il peut s'avérer difficile d'identifier les éléments s'ils ne sont pas nommés. De plus, comme nous le verrons ci-dessous, il est très simple d'extraire les éléments d'une liste par leur étiquette.

- La liste demeure un vecteur. On peut donc l'indicer avec l'opérateur [ ]. Cependant, cela retourne une liste contenant le ou les éléments indicés. C'est rarement ce que l'on souhaite.

- Pour indicer un élément d'une liste et n'obtenir que cet élément, et non une liste contenant l'élément, il faut utiliser l'opérateur d'indiçage [[ ]].

Comparer

> x[1]
$size
[1] 1 5 2

et

> x[[1]]
[1] 1 5 2

Évidemment, on ne peut extraire qu'un seul élément à la fois avec les crochets doubles [[ ]] .

- Petite subtilité peu employée, mais élégante. Si l'indice utilisé dans [[ ]] est un vecteur, il est utilisé récursivement pour indicer la liste : cela sélectionnera la composante de la liste correspondant au premier élément du vecteur, puis l'élément de la composante correspondant au second élément du vecteur, et ainsi de suite.

- Une autre — la meilleure, en fait — façon d'indicer un seul élément d'une liste est par son étiquette avec l'opérateur $ :

> x$size
[1] 1 5 2

- La fonction unlist convertit une liste en un vecteur simple. Elle est surtout utile pour concaténer les éléments d'une liste lorsque ceux-ci sont des scalaires. Attention, cette fonction peut être destructrice si la structure interne de la liste est importante.

2-7 - Data frames

Les vecteurs, les matrices, les tableaux et les listes sont les types d'objets les plus fréquemment utilisés en programmation en R. Toutefois, un grand nombre de procédures statistiques — pensons à la régression linéaire, par exemple — reposent davantage sur les data frames pour le stockage des données.

- Un est une liste de classe « data.frame » dont tous les éléments sont de la même longueur (ou comptent le même nombre de lignes si les éléments sont des matrices).

- Il est généralement représenté sous la forme d'un tableau à deux dimensions. Chaque élément de la liste sous-jacente correspond à une colonne.

- Bien que visuellement similaire à une matrice un data frame est plus général puisque les colonnes peuvent être de modes différents ; pensons à un tableau avec des noms (mode character) dans une colonne et des notes (mode numeric) dans une autre.

- On crée un data frame avec la fonction data.frame ou, pour convertir un autre type d'objet en data frame, avec as.data.frame.

- Le data frame peut être indicé à la fois comme une liste et comme une matrice.

- Les fonctions rbind et cbind peuvent être utilisées pour ajouter des lignes ou des colonnes à un data frame.

- On peut rendre les colonnes d'un data frame (ou d'une liste) visibles dans l'espace de travail avec la fonction attach, puis les masquer avec detach.

La liste permettrait de stocker à la fois le nom des équipes et leurs statistiques puisqu'elle peut contenir des objets de modes différents.

On crée d'abord des vecteurs simples contenant les données de chaque colonne du classement des équipes.

1. > Equipe <- c("Washington", "Dallas", "Chicago",

2. + "Los Angeles", "St-Louis", "Détroit",

3. + "Montréal", "Boston")

4. > MJ  <- c(55, 56, 57, 58, 56, 57, 56, 57)

5. > V   <- c(36, 32, 30, 30, 25, 25, 22, 24)

6. > D   <- c(16, 19, 21, 22, 19, 21, 27, 31)

7. > DP  <- c( 3,  5,  6,  6, 12, 11,  7,  2)

8. > PTS <- c(75, 69, 66, 66, 62, 61, 51, 50)

On les combine ensuite sous forme de liste nommée.

1. > list(Equipe = Equipe, MJ = MJ, V = V, D = D,

2. + DP = DP, PTS = PTS)

3. $Equipe

4. [1]  "Washington"   "Dallas" "Chicago"

5. [4] "Los Angeles" "St-Louis" "Détroit"

6. [7]    "Montréal"   "Boston"

7. $MJ

8. [1] 55 56 57 58 56 57 56 57

9. $V

10. [1] 36 32 30 30 25 25 22 24

11. $D

12. [1] 16 19 21 22 19 21 27 31

13. $DP

14. [1] 3 5 6 6 12 11 7 2

15. $PTS

16. [1] 75 69 66 66 62 61 51 50

On constate que R ne présente pas le contenu de la liste sous forme d'un tableau. Ce n'est donc pas le type d'objet le mieux approprié pour stocker un classement. En fait, la liste est un mode de stockage trop général pour le type de données dont nous disposons.

L'élément distinctif entre un data frame et une liste générale, c'est que tous les éléments du premier doivent être de la même longueur et que, par conséquent, R les dispose en colonnes. Nous avons donc ici le type d'objet tout désigné pour stocker des données de modes différents, mais qui se présentent sous forme de tableau à deux dimensions.

2-8 - Indiçage

L'indiçage des vecteurs et matrices a déjà été brièvement présenté aux sections 2.4 et 2.5. La présente section contient plus de détails sur cette procédure des plus communes lors de l'utilisation du langage R. On se concentre toutefois sur le traitement des vecteurs.

- L'indiçage sert principalement à deux choses : soit extraire des éléments d'un objet avec la construction x[i], ou les remplacer avec la construction x[i] <- y.

- Il est utile de savoir que ces opérations sont en fait traduites par l'interprète R en des appels à des fonctions nommées [ et [<-, dans l'ordre.

- De même, les opérations d'extraction et de remplacement d'un élément d'une liste de la forme x$etiquette et x$etiquette <- y correspondent à des appels aux fonctions $ et $<-.

Il existe cinq façons d'indicer un vecteur dans le langage R. Dans tous les cas, l'indiçage se fait à l'intérieur de crochets [ ].

1. Avec un vecteur d'entiers positifs. Les éléments se trouvant aux positions correspondant aux entiers sont extraits du vecteur, dans l'ordre. C'est la technique la plus courante :

1. > x <- c(A = 2, B = 4, C = -1, D = -5, E = 8)

2. > x[c(1, 3)]

3. A  C

4. 2 -1

2. Avec un vecteur d'entiers négatifs. Les éléments se trouvant aux positions correspondant aux entiers négatifs sont alors éliminés du vecteur :

1. > x[c(-2, -3)]

2. A  D  E

3. 2 -5  8

3. Avec un vecteur booléen. Le vecteur d'indiçage doit alors être de la même longueur que le vecteur indicé. Les éléments correspondant à une valeur TRUE sont extraits du vecteur, alors que ceux correspondant à FALSE sont éliminés :

1. > x > 0

2. A    B     C     D    E

3. TRUE TRUE FALSE FALSE TRUE

4. > x[x > 0]

5. A B E

6. 2 4 8

4. Avec un vecteur de chaînes de caractères. Utile pour extraire les éléments d'un vecteur à condition que ceux-ci soient nommés :

1. > x[c("B", "D")]

2. B  D

3. 4 -5

5. L'indice est laissé vide. Tous les éléments du vecteur sont alors sélectionnés :

1. > x[]

2. A  B  C  D  E

3. 2  4 -1 -5  8

Cette méthode est essentiellement utilisée avec les matrices et tableaux pour sélectionner tous les éléments d'une dimension (voir l'exemple). Laisser l'indice vide est différent d'indicer avec un vecteur vide ; cette dernière opération retourne un vecteur vide.

Solution du problème

Nous avons déjà créé dans cette astuce des vecteurs contenant les données des différentes colonnes du classement de la ligue.

1. > Equipe

2. [1]  "Washington"   "Dallas" "Chicago"

3. [4] "Los Angeles" "St-Louis" "Détroit"

4. [7]    "Montréal"   "Boston"

5. > MJ

6. [1] 55 56 57 58 56 57 56 57

7. > V

8. [1] 36 32 30 30 25 25 22 24

9. > D

10. [1] 16 19 21 22 19 21 27 31

11. > DP

12. [1] 3 5 6 6 12 11 7 2

13. > PTS

14. [1] 75 69 66 66 62 61 51 50

On crée l'objet classement qui contiendra le classement de la ligue avec la fonction data.frame. Celle-ci prend en arguments les différents vecteurs de données.

1. > (classement <- data.frame(Equipe, MJ, V, D, DP, PTS))

2. Equipe MJ  V  D DP PTS

3. 1  Washington 55 36 16  3  75

4. 2      Dallas 56 32 19  5  69

5. 3     Chicago 57 30 21  6  66

6. 4 Los Angeles 58 30 22  6  66

7. 5    St-Louis 56 25 19 12  62

8. 6     Détroit 57 25 21 11  61

9. 7    Montréal 56 22 27  7  51

10. 8      Boston 57 24 31  2  50

On répond ensuite aux questions de traitement des données.

1. Le nombre d'équipes dans la ligue correspond au nombre de lignes de l'objet :

> dim(classement)[1]
[1] 8

ou, comme nous le verrons dans le code informatique de la section suivante,

> nrow(classement)
[1] 8

2. Pour extraire la fiche complète de l'équipe de Montréal, il suffit d'extraire cette ligne du tableau de données :

1. > ## indiçage par vecteur booléen

2. > classement[classement$Equipe == "Montréal", ]

3. Equipe MJ  V  D DP PTS

4. 7 Montréal 56 22 27  7  51

5. > ## extraction par une fonction pratique et élégante

6. > subset(classement, Equipe == "Montréal")

7. Equipe MJ  V  D DP PTS

8. 7 Montréal 56 22 27  7  51

3. Pour extraire la fiche complète de l'équipe à la septième position du classement, on peut indicer le tableau par position :

1. > classement[7, ]

2. Equipe MJ  V  D DP PTS

3. 7 Montréal 56 22 27  7  51

2-9 - Exemples

1. ###

2. ### COMMANDES R

3. ###

4. ## Les expressions entrées à la ligne de commande sont

5. ## immédiatement évaluées et le résultat est affiché à

6. ## l'écran, comme avec une grosse calculatrice.

7. 1                               # une constante

8. (2 + 3 * 5)/7                   # priorité des opérations

9. 3^5                             # puissance

10. exp(3)                          # fonction exponentielle

11. sin(pi/2) + cos(pi/2)           # fonctions trigonométriques

12. gamma(5)                        # fonction gamma

14. ## Lorsqu'une expression est syntaxiquement incomplète,

15. ## l'invite de commande change de '> ' à '+ '.

16. 2 -                             # expression incomplète

17. 5 *                             # toujours incomplète

18. 3                               # complétée

20. ## Taper le nom d'un objet affiche son contenu. Pour une

21. ## fonction, c'est son code source qui est affiché.

22. pi                              # constante numérique intégrée

23. letters                         # chaîne de caractères intégrée

24. LETTERS                         # version en majuscules

25. matrix                          # fonction

27. ## Ne pas utiliser '=' pour l'affectation. Les opérateurs

28. ## d'affectation standard en R sont '<-' et '->'.

29. x <- 5                       # affecter 5 à l'objet 'x'

30. 5 -> x                          # idem, mais peu usité

31. x                               # voir le contenu

32. (x <- 5)                     # affecter et afficher

33. y <- x                       # affecter la valeur de 'x' à 'y'

34. x <- y <- 5               # idem, en une seule expression

35. y                               # 5

36. x <- 0                       # changer la valeur de 'x'...

37. y                               # ... ne change pas celle de 'y'

39. ## Pour regrouper plusieurs expressions en une seule commande,

40. ## il faut soit les séparer par un point-virgule ';', soit les

41. ## regrouper à l'intérieur d'accolades { } et les séparer par

42. ## des retours à la ligne.

43. x <- 5; y <- 2; x + y     # compact; éviter dans les scripts

44. x <- 5;                      # éviter les ';' superflus

45. {                               # début d'un groupe

46. x <- 5                   # première expression du groupe

47. y <- 2                   # seconde expression du groupe

48. x + y                       # résultat du groupe

49. }                               # fin du groupe et résultat

50. {x <- 5; y <- 2; x + y}   # valide, mais redondant

52. ###

53. ### NOMS D'OBJETS

54. ###

56. ## Quelques exemples de noms valides et invalides.

57. foo <- 5                     # valide

58. foo.123 <- 5                 # valide

59. foo_123 <- 5                 # valide

60. 123foo <- 5                  # invalide; commence par un chiffre

61. .foo <- 5                    # valide

62. .123foo <- 5                 # invalide; point suivi d'un chiffre

64. ## Liste des objets dans l'espace de travail. Les objets dont

65. ## le nom commence par un point sont considérés cachés.

66. ls()                            # l'objet '.foo' n'est pas affiché

67. ls(all.names = TRUE)            # objets cachés aussi affichés

69. ## R est sensible à la casse

70. foo <- 1

71. Foo

72. FOO

74. ###

75. ### LES OBJETS R

76. ###

78. ## MODES ET TYPES DE DONNÉES

80. ## Le mode d'un objet détermine ce qu'il peut contenir. Les

81. ## vecteurs simples ("atomic") contiennent des données d'un

82. ## seul type.

83. mode(c(1, 4.1, pi))             # nombres réels

84. mode(c(2, 1 + 5i))              # nombres complexes

85. mode(c(TRUE, FALSE, TRUE))      # valeurs booléennes

86. mode("foobar")                  # chaînes de caractères

88. ## Si l'on mélange dans un même vecteur des objets de modes

89. ## différents, il y a conversion automatique vers le mode pour

90. ## lequel il y a le moins de perte d'information, c'est-à-dire

91. ## vers le mode qui permet le mieux de retrouver la valeur

92. ## originale des éléments.

93. c(5, TRUE, FALSE)               # conversion en mode 'numeric'

94. c(5, "z")                       # conversion en mode 'character'

95. c(TRUE, "z")                    # conversion en mode 'character'

96. c(5, TRUE, "z")                 # conversion en mode 'character'

98. ## La plupart des autres types d'objets sont récursifs. Voici

99. ## quelques autres modes.

100. mode(seq)                       # une fonction

101. mode(list(5, "foo", TRUE))      # une liste

102. mode(expression(x <- 5))     # une expression non évaluée

104. ## LONGUEUR

106. ## La longueur d'un vecteur est égale au nombre d'éléments

107. ## dans le vecteur.

108. (x <- 1:4)

109. length(x)

111. ## Une chaîne de caractères ne compte que pour un seul

112. ## élément.

113. (x <- "foobar")

114. length(x)

116. ## Pour obtenir la longueur de la chaîne, il faut utiliser

117. ## nchar().

118. nchar(x)

120. ## Un objet peut néanmoins contenir plusieurs chaînes de

121. ## caractères.

122. (x <- c("f", "o", "o", "b", "a", "r"))

123. length(x)

125. ## La longueur peut être 0, auquel cas on a un objet vide,

126. ## mais qui existe.

127. (x <- numeric(0))            # création du contenant

128. length(x)                       # l'objet 'x' existe...

129. x[1] <- 1                    # possible, 'x' existe

130. X[1] <- 1                    # impossible, 'X' n'existe pas

132. ## L'OBJET SPECIAL 'NULL'

133. mode(NULL)                      # le mode de 'NULL' est NULL

134. length(NULL)                    # longueur nulle

135. x <- c(NULL, NULL)           # s'utilise comme un objet normal

136. x; length(x); mode(x)           # mais donne toujours le vide

138. ## L'OBJET SPÉCIAL 'NA'

139. x <- c(65, NA, 72, 88)       # traité comme une valeur

140. x + 2                           # tout calcul avec 'NA' donne NA

141. mean(x)                         # voilà qui est pire

142. mean(x, na.rm = TRUE)           # éliminer les 'NA' avant le calcul

143. is.na(x)                        # tester si les données sont 'NA'

145. ## VALEURS INFINIES ET INDÉTERMINÉES

146. 1/0                             # +infini

147. -1/0                            # -infini

148. 0/0                             # indétermination

149. x <- c(65, Inf, NaN, 88)     # s'utilisent comme des valeurs

150. is.finite(x)                    # quels sont les nombres réels ?

151. is.nan(x)                       # lesquels ne sont « pas un nombre » ?

153. ## ATTRIBUTS

155. ## Les objets peuvent être dotés d'un ou plusieurs attributs.

156. data(cars)                      # jeu de données intégré

157. attributes(cars)                # liste de tous les attributs

158. attr(cars, "class")             # extraction d'un seul attribut

160. ## Attribut 'class'. Selon la classe d'un objet, certaines

161. ## fonctions (dites « fonctions génériques ») vont se comporter

162. ## différemment.

163. x <- sample(1:100, 10)       # échantillon aléatoire de 10

164. # nombres entre 1 et 100

165. class(x)                        # classe de l'objet

166. plot(x)                         # graphique pour cette classe

167. class(x) <- "ts"             # 'x' est maintenant une série

168. # chronologique

169. plot(x)                         # graphique pour les séries

170. # chronologiques

171. class(x) <- NULL; x          # suppression de l'attribut 'class'

173. ## Attribut 'dim'. Si l'attribut 'dim' compte deux valeurs,

174. ## l'objet est traité comme une matrice. S'il en compte plus

175. ## de deux, l'objet est traité comme un tableau (array).

176. x <- 1:24                    # un vecteur

177. dim(x) <- c(4, 6)            # ajoute un attribut 'dim'

178. x                               # l'objet est une matrice

179. dim(x) <- c(4, 2, 3)         # change les dimensions

180. x                               # l'objet est maintenant un tableau

182. ## Attribut 'dimnames'. Permet d'assigner des étiquettes (ou

183. ## noms) aux dimensions d'une matrice ou d'un tableau.

184. dimnames(x) <- list(1:4, c("a", "b"), c("A", "B", "C"))

185. dimnames(x)                     # remarquer la conversion

186. x                               # affichage avec étiquettes

187. attributes(x)                   # tous les attributs de 'x'

188. attributes(x) <- NULL; x     # supprimer les attributs

190. ## Attributs 'names'. Similaire à 'dimnames', mais pour les

191. ## éléments d'un vecteur ou d'une liste.

192. names(x) <- letters[1:24]    # attribution d'étiquettes

193. x                               # identification facilitée

195. ###

196. ### VECTEURS

197. ###

199. ## La fonction de base pour créer des vecteurs est 'c'. Il

200. ## peut s'avérer utile de donner des étiquettes aux éléments

201. ## d'un vecteur.

202. x <- c(a = -1, b = 2, c = 8, d = 10) # création d'un vecteur

203. names(x)                                # extraire les étiquettes

204. names(x) <- letters[1:length(x)]     # changer les étiquettes

205. x[1]                            # extraction par position

206. x["c"]                          # extraction par étiquette

207. x[-2]                           # élimination d'un élément

209. ## La fonction 'vector' sert à initialiser des vecteurs avec

210. ## des valeurs prédéterminées. Elle compte deux arguments : le

211. ## mode du vecteur et sa longueur. Les fonctions 'numeric',

212. ## 'logical', 'complex' et 'character' constituent des

213. ## raccourcis pour des appels à 'vector'.

214. vector("numeric", 5)            # vecteur initialisé avec des 0

215. numeric(5)                      # équivalent

216. numeric                         # en effet, voici la fonction

217. logical(5)                      # initialisé avec FALSE

218. complex(5)                      # initialisé avec 0 + 0i

219. character(5)                    # initialisé avec chaînes vides

221. ###

222. ### MATRICES ET TABLEAUX

223. ###

225. ## Une matrice est un vecteur avec un attribut 'dim' de

226. ## longueur 2 une classe implicite "matrix". La manière

227. ## naturelle de créer une matrice est avec la fonction

228. ## 'matrix'.

229. (x <- matrix(1:12, nrow = 3, ncol = 4)) # créer la matrice

230. length(x)                       # 'x' est un vecteur...

231. dim(x)                          # ... avec un attribut 'dim'...

232. class(x)                        # ... et classe implicite "matrix"

234. ## Une manière moins naturelle mais équivalente --- et parfois

235. ## plus pratique --- de créer une matrice consiste à ajouter

236. ## un attribut 'dim' à un vecteur.

237. x <- 1:12                    # vecteur simple

238. dim(x) <- c(3, 4)            # ajout d'un attribut 'dim'

239. x; class(x)                     # 'x' est une matrice!

241. ## Les matrices sont remplies par colonne par défaut. Utiliser

242. ## l'option 'byrow' pour remplir par ligne.

243. matrix(1:12, nrow = 3, byrow = TRUE)

245. ## Indicer la matrice ou le vecteur sous-jacent est

246. ## équivalent. Utiliser l'approche la plus simple selon le

247. ## contexte.

248. x[1, 3]                         # l'élément en position (1, 3)...

249. x[7]                            # ... est le 7e élément du vecteur

250. x[1, ]                          # première ligne

251. x[, 2]                          # deuxième colonne

252. nrow(x)                         # nombre de lignes

253. dim(x)[1]                       # idem

254. ncol(x)                         # nombre de colonnes

255. dim(x)[2]                       # idem

257. ## Fusion de matrices et vecteurs.

258. x <- matrix(1:12, 3, 4)      # 'x' est une matrice 3 x 4

259. y <- matrix(1:8, 2, 4)       # 'y' est une matrice 2 x 4

260. z <- matrix(1:6, 3, 2)       # 'z' est une matrice 3 x 2

261. rbind(x, 1:4)                   # ajout d'une ligne à 'x'

262. rbind(x, y)                     # fusion verticale de 'x' et 'y'

263. cbind(x, 1:3)                   # ajout d'une colonne à 'x'

264. cbind(x, z)                     # concaténation de 'x' et 'z'

265. rbind(x, z)                     # dimensions incompatibles

266. cbind(x, y)                     # dimensions incompatibles

268. ## Les vecteurs ligne et colonne sont rarement nécessaires. On

269. ## peut les créer avec les fonctions 'rbind' et 'cbind',

270. ## respectivement.

271. rbind(1:3)                      # un vecteur ligne

272. cbind(1:3)                      # un vecteur colonne

274. ## Un tableau (array) est un vecteur avec un attribut 'dim' de

275. ## longueur supérieure à 2 et une classe implicite "array".

276. ## Quant au reste, la manipulation des tableaux est en tous

277. ## points identique à celle des matrices. Ne pas oublier:

278. ## les tableaux sont remplis de la première dimension à la

279. ## dernière!

280. x <- array(1:60, 3:5)        # tableau 3 x 4 x 5

281. length(x)                       # 'x' est un vecteur...

282. dim(x)                          # ... avec un attribut 'dim'...

283. class(x)                        # ... une classe implicite "array"

284. x[1, 3, 2]                      # l'élément en position (1, 3, 2)...

285. x[19]                           # ... est l'élément 19 du vecteur

287. ## Le tableau ci-dessus est un prisme rectangulaire 3 unités

288. ## de haut, 4 de large et 5 de profond. Indicer ce prisme avec

289. ## un seul indice équivaut à en extraire des « tranches », alors

290. ## qu'utiliser deux indices équivaut à en tirer des « carottes »

291. ## (au sens géologique du terme). Il est laissé en exercice de

292. ## généraliser à plus de dimensions...

293. x                               # les cinq matrices

294. x[, , 1]                        # tranches de haut en bas

295. x[, 1, ]                        # tranches d'avant à l'arrière

296. x[1, , ]                        # tranches de gauche à droite

297. x[, 1, 1]                       # carotte de haut en bas

298. x[1, 1, ]                       # carotte d'avant à l'arrière

299. x[1, , 1]                       # carotte de gauche à droite

301. ###

302. ### LISTES

303. ###

305. ## La liste est l'objet le plus général en R. C'est un objet

306. ## récursif qui peut contenir des objets de n'importe quel

307. ## mode et longueur.

308. (x <- list(joueur = c("V", "C", "C", "M", "A"),

309. score = c(10, 12, 11, 8, 15),

310. expert = c(FALSE, TRUE, FALSE, TRUE, TRUE),

311. niveau = 2))

312. is.vector(x)                    # vecteur...

313. length(x)                       # ... de quatre éléments...

314. mode(x)                         # ... de mode "list"

315. is.recursive(x)                 # objet récursif

317. ## Comme tout autre vecteur, une liste peut être concaténée

318. ## avec un autre vecteur avec la fonction 'c'.

319. y <- list(TRUE, 1:5)         # liste de deux éléments

320. c(x, y)                         # liste de six éléments

322. ## Pour initialiser une liste d'une longueur déterminée, mais

323. ## dont chaque élément est vide, utiliser la fonction

324. ## 'vector'.

325. vector("list", 5)               # liste de NULL

327. ## Pour extraire un élément d'une liste, il faut utiliser les

328. ## doubles crochets [[ ]]. Les simples crochets [ ]

329. ## fonctionnent aussi, mais retournent une sous-liste -- ce

330. ## qui est rarement ce que l'on souhaite.

331. x[[1]]                          # premier élément de la liste...

332. mode(x[[1]])                    # ... un vecteur

333. x[1]                            # aussi le premier élément...

334. mode(x[1])                      # ... mais une sous-liste...

335. length(x[1])                    # ... d'un seul élément

336. x[[2]][1]                       # 1er élément du 2e élément

337. x[[c(2, 1)]]                    # idem, par indiçage récursif

339. ## Les éléments d'une liste étant généralement nommés (c'est

340. ## une bonne habitude à prendre! ), il est souvent plus simple

341. ## et sûr d'extraire les éléments d'une liste par leur

342. ## étiquette.

343. x$joueur                        # équivalent à a[[1]]

344. x$score[1]                      # équivalent à a[[c(2, 1)]]

345. x[["expert"]]                   # aussi valide, mais peu usité

346. x$level <- 1                 # aussi pour l'affectation

348. ## Une liste peut contenir n'importe quoi...

349. x[[5]] <- matrix(1, 2, 2)    # ... une matrice...

350. x[[6]] <- list(20:25, TRUE)  # ... une autre liste...

351. x[[7]] <- seq                # ... même le code d'une fonction!

352. x                               # eh ben!

353. x[[c(6, 1, 3)]]                 # de quel élément s'agit-il ?

355. ## Pour supprimer un élément d'une liste, lui assigner la

356. ## valeur 'NULL'.

357. x[[7]] <- NULL; length(x)    # suppression du 7e élément

359. ## Il est parfois utile de convertir une liste en un simple

360. ## vecteur. Les éléments de la liste sont alors « déroulés », y

361. ## compris la matrice en position 5 (qui n'est rien d'autre

362. ## qu'un vecteur, on s'en souviendra).

363. unlist(x)                       # remarquer la conversion

364. unlist(x, recursive = FALSE)    # ne pas appliquer aux sous-listes

365. unlist(x, use.names = FALSE)    # éliminer les étiquettes

367. ###

368. ### DATA FRAMES

369. ###

371. ## Un data frame est une liste dont les éléments sont tous de

372. ## même longueur. Il comporte un attribut 'dim', ce qui fait

373. ## qu'il est représenté comme une matrice. Cependant, les

374. ## colonnes peuvent être de modes différents.

375. (DF <- data.frame(Noms = c("Pierre", "Jean", "Jacques"),

376. Age = c(42, 34, 19),

377. Fumeur = c(TRUE, TRUE, FALSE)))

378. mode(DF)                        # un data frame est une liste...

379. class(DF)                       # ... de classe 'data.frame'

380. dim(DF)                         # dimensions implicites

381. names(DF)                       # titres des colonnes

382. row.names(DF)                   # titres des lignes (implicites)

383. DF[1, ]                         # première ligne

384. DF[, 1]                         # première colonne

385. DF$Noms                         # idem, mais plus simple

387. ## Lorsque l'on doit travailler longtemps avec les différentes

388. ## colonnes d'un data frame, il est pratique de pouvoir y

389. ## accéder directement sans devoir toujours indicer. La

390. ## fonction 'attach' permet de rendre les colonnes

391. ## individuelles visibles dans l'espace de travail. Une fois

392. ## le travail terminé, 'detach' masque les colonnes.

393. exists("Noms")                  # variable n'existe pas

394. attach(DF)                      # rendre les colonnes visibles

395. exists("Noms")                  # variable existe

396. Noms                            # colonne accessible

397. detach(DF)                      # masquer les colonnes

398. exists("Noms")                  # variable n'existe plus

400. ###

401. ### INDIÇAGE

402. ###

404. ## Les opérations suivantes illustrent les différentes

405. ## techniques d'indiçage d'un vecteur pour l'extraction et

406. ## l'affectation, c'est-à-dire que l'on utilise à la fois la

407. ## fonction '[' et la fonction '[<-'. Les mêmes techniques

408. ## existent aussi pour les matrices, tableaux et listes.

409. ##

410. ## On crée d'abord un vecteur quelconque formé de vingt

411. ## nombres aléatoires entre 1 et 100 avec répétitions

412. ## possibles.

413. (x <- sample(1:100, 20, replace = TRUE))

415. ## On ajoute des étiquettes aux éléments du vecteur à partir

416. ## de la variable interne 'letters'.

417. names(x) <- letters[1:20]

419. ## On génère ensuite cinq nombres aléatoires entre 1 et 20

420. ## (sans répétitions).

421. (y <- sample(1:20, 5))

423. ## On remplace maintenant les éléments de 'x' correspondant

424. ## aux positions dans le vecteur 'y' par des données

425. ## manquantes.

426. x[y] <- NA

427. x

429. ## Les cinq méthodes d'indiçage de base.

430. x[1:10]                         # avec des entiers positifs

431. "["(x, 1:10)                    # idem, avec la fonction '['

432. x[-(1:3)]                       # avec des entiers négatifs

433. x[x < 10]                    # avec un vecteur booléen

434. x[c("a", "k", "t")]             # par étiquettes

435. x[]                             # aucun indice...

436. x[numeric(0)]                   # ... différent d'indice vide

438. ## Il arrive souvent de vouloir indicer spécifiquement les

439. ## données manquantes d'un vecteur (pour les éliminer ou les

440. ## remplacer par une autre valeur, par exemple). Pour ce

441. ## faire, on utilise la fonction 'is.na' et l'indiçage par un

442. ## vecteur booléen. (Note : l'opérateur '! ' ci-dessous est la

443. ## négation logique.)

444. is.na(x)                        # positions des données manquantes

445. x[!is.na(x)]                    # suppression des données manquantes

446. x[is.na(x)] <- 0; x          # remplace les NA par des 0

447. "[<-"(x, is.na(x), 0)        # idem, mais très peu usité

449. ## On laisse tomber les étiquettes de l'objet.

450. names(x) <- NULL

452. ## Quelques cas spéciaux d'indiçage.

453. length(x)                       # un rappel

454. x[1:25]                         # allonge le vecteur avec des NA

455. x[25] <- 10; x               # remplit les trous avec des NA

456. x[0]                            # n'extrait rien

457. x[0] <- 1; x                 # n'affecte rien

458. x[c(0, 1, 2)]                   # le 0 est ignoré

459. x[c(1, NA, 5)]                  # indices NA retourne NA

460. x[2.6]                          # fractions tronquées vers 0

462. ## On laisse tomber les 5 derniers éléments et on convertit le

463. ## vecteur en une matrice 4 x 5.

464. x <- x[1:20]                 # ou x[-(21:25)]

465. dim(x) <- c(4, 5); x         # ajouter un attribut 'dim'

467. ## Dans l'indiçage des matrices et tableaux, l'indice de

468. ## chaque dimension obéit aux mêmes règles que ci-dessus. On

469. ## peut aussi indicer une matrice (ou un tableau) avec une

470. ## matrice. Si les exemples ci-dessous ne permettent pas d'en

471. ## comprendre le fonctionnement, consulter la rubrique d'aide

472. ## de la fonction '[' (ou de 'Extract').

473. x[1, 2]                         # élément en position (1, 2)

474. x[1, -2]                        # 1re rangée sans 2e colonne

475. x[c(1, 3), ]                    # 1re et 3e rangées

476. x[-1, ]                         # supprimer 1re rangée

477. x[, -2]                         # supprimer 2e colonne

478. x[x[, 1] > 10, ]                # lignes avec 1er élément > 10

479. x[rbind(c(1, 1), c(2, 2))]      # éléments x[1, 1] et x[2, 2]

480. x[cbind(1:4, 1:4)]              # éléments x[i, i] (diagonale)

481. diag(x)                         # idem et plus explicite

2-10 - Exercices

2.1 (solution)

1. Écrire une expression R pour créer la liste suivante :

1. > x

2. [[1]]

3. [1] 1 2 3 4 5

4. $data

5. [,1] [,2] [,3]

6. [1,]    1    3    5

7. [2,]    2    4    6

8. [[3]]

9. 2.10. Exercices 47

10. [1] 0 0 0

11. $test

12. [1] FALSE FALSE FALSE FALSE

2. Extraire les étiquettes de la liste.

3. Trouver le mode et la longueur du quatrième élément de la liste.

4. Extraire les dimensions du second élément de la liste.

5. Extraire les deuxième et troisième éléments du second élément de la liste.

6. Remplacer le troisième élément de la liste par le vecteur 3:8.

2.2 (solution) Soit x un vecteur contenant les valeurs d'un échantillon :

1. > x

2. [1] 1 18  2  1 5 2 6 1 12 3 13 8 20 1 5 7

3. [17] 7  4 14 10

Écrire une expression R permettant d'extraire les éléments suivants :

1. le deuxième élément de l'échantillon ;

2. les cinq premiers éléments de l'échantillon ;

3. les éléments strictement supérieurs à 14 ;

4. tous les éléments sauf les éléments en positions 6, 10 et 12.

2.3 (solution) Soit x une matrice 10 × 7 obtenue aléatoirement avec

> x <- matrix(sample(1:100, 70), 7, 10)

Écrire des expressions R permettant d'obtenir les éléments de la matrice demandés ci-dessous.

1. L'élément (4, 3).

2. Le contenu de la sixième ligne.

3. Les première et quatrième colonnes (simultanément).

4. Les lignes dont le premier élément est supérieur à 50.

3 - Opérateurs et fonctions

Objectifs du chapitre

- Tirer profit de l'arithmétique vectorielle caractéristique du langage R dans les calculs.

- Utiliser les opérateurs R les plus courants, notamment pour le traitement des vecteurs, le calcul de sommaires et la manipulation des matrices et tableaux.

- Faire l'appel d'une fonction dans R ; concevoir comment les arguments sont passés à la fonction et le traitement des valeurs par défaut.

- Utiliser la fonction if pour l'exécution conditionnelle de commandes R.

- Distinguer la construction if() … else de la fonction ifelse.

- Faire des boucles en R.

- Choisir entre les opérateurs for, while et repeat lors de la construction d'une boucle R.

Ce chapitre présente les principaux opérateurs arithmétiques, fonctions mathématiques et structures de contrôle disponibles dans R. La liste est évidemment loin d'être exhaustive, surtout étant donné l'évolution rapide du langage. Un des meilleurs endroits pour découvrir de nouvelles fonctions demeure la section See Also des rubriques d'aide, qui offre des hyperliens vers des fonctions apparentées au sujet de la rubrique.

Énoncé du problème

On s'intéresse à la somme des résultats du lancer de deux dés. Le premier dé compte 8 faces et le deuxième, 6. On souhaite calculer la moyenne des résultats de la somme supérieurs à 7.

3-1 - Opérations arithmétiques

L'unité de base en R est le vecteur.

- Les opérations sur les vecteurs sont effectuées élément par élément :

1. > c(1, 2, 3) + c(4, 5, 6)

2. [1] 5 7 9

3. > 1:3 * 4:6

4. [1] 4 10 18

- Si les vecteurs impliqués dans une expression arithmétique ne sont pas de la même longueur, les plus courts sont recyclés de façon à correspondre au plus long vecteur. Cette règle est particulièrement apparente avec les vecteurs de longueur 1 :

1. > 1:10 + 2

2. [1] 3 4 5 6 7 8 9 10 11 12

3. > 1:10 + rep(2, 10)

4. [1] 3 4 5 6 7 8 9 10 11 12

- Si la longueur du plus long vecteur est un multiple de celle du ou des autres vecteurs, ces derniers sont recyclés un nombre entier de fois :

1. > 1:10 + 1:5 + c(2, 4) # vecteurs recyclés 2 et 5 fois

2. [1] 4 8 8 12 12 11 11 15 15 19

3. > 1:10 + rep(1:5, 2) + rep(c(2, 4), 5) # équivalent

4. [1] 4 8 8 12 12 11 11 15 15 19

- Sinon, le plus court vecteur est recyclé un nombre fractionnaire de fois, mais comme ce résultat est rarement souhaité et provient généralement d'une erreur de programmation, un avertissement est affiché :

1. > 1:10 + c(2, 4, 6)

2. [1] 3 6 9 6 9 12 9 12 15 12

3. Message d′avis :

4. In 1:10 + c(2, 4, 6) :

5. la taille d′un objet plus long n′est pas un multiple de la taille d′un objet plus court

Tab. 3.1 - Principaux opérateurs du langage R, en ordre décroissant de priorité

Opérateur

Fonction

$

extraction d'une liste

^

puissance

-

changement de signe

:

génération de suites

%*%%%%/%

produit matriciel, modulo, division entière

*/

multiplication, division

+-

addition, soustraction

<<===>=>!=

plus petit, plus petit ou égal, égal, plus grand ou égal, plus grand, différent de

!

négation logique

&&&

« et » logique

|, ||

« ou » logique

->->>

assignation

<-<<-

assignation

Grâce à la propriété d'opération élément par élément de R, il devrait être possible — en fait, il est souhaitable — de résoudre le problème sans avoir recours à des boucles.

3-2 - Opérateurs

Le tableau 3.1 présente les opérateurs mathématiques et logiques les plus fréquemment employés, en ordre décroissant de priorité des opérations. Le tableau contient également les opérateurs d'assignation et d'extraction présentés au chapitre précédent ; il est utile de connaître leur niveau de priorité dans les expressions R.

Les opérateurs de puissance (^) et d'assignation à gauche (<-, <<-) sont évalués de droite à gauche ; tous les autres de gauche à droite. Ainsi, 2^2^3 est 2^8, et non 4^3, alors que 1 - 1 - 1 vaut -1, et non 1.

3-3 - Appels de fonctions

Les opérateurs du tableau 3.1 constituent des raccourcis utiles pour accéder aux fonctions les plus courantes de R. Pour toutes les autres, il faut appeler la fonction directement. Cette section passe en revue les règles d'appels d'une fonction et la façon de spécifier les arguments, qu'il s'agisse d'une fonction interne de R ou d'une fonction personnelle (voir le chapitre 5).

- Il n'y a pas de limite pratique quant au nombre d'arguments que peut avoir une fonction.

- Les arguments d'une fonction peuvent être spécifiés selon l'ordre établi dans la définition de la fonction. Cependant, il est beaucoup plus prudent et fortement recommandé de spécifier les arguments par leur nom, avec une construction de la forme nom = valeur, surtout après les deux ou trois premiers arguments.

- L'ordre des arguments est important ; il est donc nécessaire de les nommer s'ils ne sont pas appelés dans l'ordre.

- Certains arguments ont une valeur par défaut qui sera utilisée si l'argument n'est pas spécifié dans l'appel de la fonction.

Par exemple, la définition de la fonction matrix est la suivante :

matrix(data = NA, nrow = 1, ncol = 1, byrow = FALSE,
dimnames = NULL)

- La fonction compte cinq arguments : data, nrow, ncol, byrow et dimnames.

- Ici, chaque argument a une valeur par défaut (ce n'est pas toujours le cas). Ainsi, un appel à matrix sans argument résulte en une matrice 1 × 1 remplie par colonne (sans importance, ici) de l'objet NA et dont les dimensions sont dépourvues d'étiquettes :

> matrix()
[,1]
[1,] NA

- Appel plus élaboré utilisant tous les arguments. Le premier argument est rarement nommé :

1. > matrix(1:6, nrow = 2, ncol = 3, byrow = TRUE,

2. + dimnames = list(c("Gauche", "Droit"),

3. + c("Rouge", "Vert", "Bleu")))

4. Rouge Vert Bleu

5. Gauche     1    2    3

6. Droit     4    5    6

3-4 - Quelques fonctions utiles

Le langage R compte un très grand nombre (des milliers !) de fonctions internes. Cette section en présente quelques-unes seulement, les fonctions de base les plus souvent utilisées pour programmer en R et pour manipuler des données.

Pour chaque fonction présentée dans les sections suivantes, on fournit un ou deux exemples d'utilisation. Ces exemples sont souvent loin de couvrir toutes les utilisations possibles d'une fonction. La section 3.7 fournit des exemples additionnels, mais il est recommandé de consulter les diverses rubriques d'aide pour connaître toutes les options des fonctions.

3-4-1 - Manipulation de vecteurs

seq

génération de suites de nombres

> seq(1, 9, by = 2)
[1] 1 3 5 7 9

seq_len

version plus rapide de seq pour générer la suite des nombres de 1 à la valeur de l'argument

> seq_len(10)
[1] 1 2 3 4 5 6 7 8 9 10

rep

répétition de valeurs ou de vecteurs

> rep(2, 10)
[1] 2 2 2 2 2 2 2 2 2 2

sort

tri en ordre croissant ou décroissant

> sort(c(4, -1, 2, 6))
[1] -1 2 4 6

rank

rang des éléments d'un vecteur dans l'ordre croissant ou décroissant

> rank(c(4, -1, 2, 6))
[1] 3 1 2 4

order

ordre d'extraction des éléments d'un vecteur pour les placer en ordre croissant ou décroissant

> order(c(4, -1, 2, 6))
[1] 2 3 1 4

rev

renverser un vecteur

> rev(1:10)
[1] 10 9 8 7 6 5 4 3 2 1

head

extraction des n premiers éléments d'un vecteur (n > 0) ou suppression des n derniers (n < 0)

> head(1:10, 3); head(1:10, -3)
[1] 1 2 3
[1] 1 2 3 4 5 6 7

tail

extraction des n derniers éléments d'un vecteur (n > 0) ou suppression des n premiers (n < 0)

> tail(1:10, 3); tail(1:10, -3)
[1]  8  9 10
[1]  4  5  6  7  8  9 10

unique

extraction des éléments différents d'un vecteur

> unique(c(2, 4, 2, 5, 9, 5, 0))
[1] 2 4 5 9 0

3-4-2 - Recherche d'éléments dans un vecteur

Les fonctions de cette sous-section sont toutes illustrées avec le vecteur

> x
[1] 4 -1 2 -3 6

which

positions des valeurs TRUE dans un vecteur booléen

> which(x < 0)
[1] 2 4

which.min

position du minimum dans un vecteur

> which.min(x)
[1] 4

which.max

position du maximum dans un vecteur

> which.max(x)
[1] 5

match

position de la première occurrence d'un élément dans un vecteur

> match(2, x)
[1] 3

%in%

appartenance d'une ou plusieurs valeurs à un vecteur

> -1:2 %in% x
[1] TRUE FALSE FALSE TRUE

3-4-3 - Arrondi

Les fonctions de cette sous-section sont toutes illustrées avec le vecteur

> x
[1] -3.6800000 -0.6666667 3.1415927 0.3333333
[5]  2.5200000

round

arrondi à un nombre défini de décimales (par défaut 0)

> round(x)
[1] -4 -1 3 0 3
> round(x, 3)
[1] -3.680 -0.667 3.142 0.333 2.520

floor

plus grand entier inférieur ou égal à l'argument

> floor(x)
[1] -4 -1 3 0 2

ceiling

plus petit entier supérieur ou égal à l'argument

> ceiling(x)
[1] -3 0 4 1 3

trunc

troncature vers zéro ; différent de floor pour les nombres négatifs

> trunc(x)
[1] -3 0 3 0 2

3-4-4 - Sommaires et statistiques descriptives

Les fonctions de cette sous-section sont toutes illustrées avec le vecteur

> x
[1] 14 17 7 9 3 4 25 21 24 11

sum, prod

somme et produit des éléments d'un vecteur

> sum(x); prod(x)
[1] 135
[1] 24938020800

diff

différences entre les éléments d'un vecteur (opérateur mathématique ∇)

> diff(x)
[1] 3 -10 2 -6 1 21 -4 3 -13

mean

moyenne arithmétique (et moyenne tronquée avec l'argument trim)

> mean(x)
[1] 13.5

var, sd

variance et écart type (versions sans biais)

> var(x)
[1] 64.5

min, max

minimum et maximum d'un vecteur

> min(x); max(x)
[1] 3
[1] 25

range

vecteur contenant le minimum et le maximum d'un vecteur

> range(x)
[1] 3 25

median

médiane empirique

> median(x)
[1] 12.5

quantile

quantiles empiriques

> quantile(x)
0% 25%  50%  75% 100%
3.0 7.5 12.5 20.0 25.0

summary

statistiques descriptives d'un échantillon

> summary(x)
Min. 1st Qu. Median Mean 3rd Qu. Max.
3.0     7.5   12.5 13.5    20.0 25.0

La fonction mean permettra de calculer la moyenne des résultats supérieurs à 7.

3-4-5 - Sommaires cumulatifs et comparaisons élément par élément

Les fonctions de cette sous-section sont toutes illustrées avec le vecteur

> x
[1] 14 17 7 9 3

cumsum, cumprod

somme et produit cumulatif d'un vecteur

> cumsum(x); cumprod(x)
[1] 14  31   38    47    50
[1] 14 238 1666 14994 44982

cummin, cummax

minimum et maximum cumulatif

> cummin(x); cummax(x)
[1] 14 14  7  7  3
[1] 14 17 17 17 17

pmin, pmax

minimum et maximum élément par élément (en parallèle) entre deux vecteurs ou plus

> pmin(x, 12)
[1] 12 12 7 9 3
> pmax(x, c(16, 23, 4, 12, 3))
[1] 16 23 7 12 3

3-4-6 - Opérations sur les matrices

Les fonctions de cette sous-section sont toutes illustrées avec la matrice

> x
[,1] [,2]
[1,]    2    4
[2,]    1    3

nrow, ncol

nombre de lignes et de colonnes d'une matrice

> nrow(x); ncol(x)
[1] 2
[1] 2

rowSums, colSums

sommes par ligne et par colonne, respectivement, des éléments d'une matrice ; voir aussi la fonction apply à la section 6.2

> rowSums(x)
[1] 6 4

rowMeans, colMeans

moyennes par ligne et par colonne, respectivement, des éléments d'une matrice ; voir aussi la fonction apply à la section 6.2

> colMeans(x)
[1] 1.5 3.5

t

transposée

> t(x)
[,1] [,2]
[1,]    2    1
[2,]    4    3

det

déterminant

> det(x)
[1] 2

solve

1) avec un seul argument (une matrice carrée) : inverse d'une matrice ; 2) avec deux arguments (une matrice carrée et un vecteur) : solution du système d'équations linéaires Ax = b

> solve(x)
[,1] [,2]
[1,]  1.5   -2
[2,] -0.5 1
> solve(x, c(1, 2))
[1] -2.5 1.5

diag

1) avec une matrice en argument : diagonale de la matrice ; 2) avec un vecteur en argument : matrice diagonale formée avec le vecteur ; 3) avec un scalaire p en argument : matrice identité p × p

> diag(x)
[1] 2 3

3-4-7 - Produit extérieur

La fonction outer calcule le produit extérieur entre deux vecteurs. Ce n'est pas la fonction la plus intuitive à utiliser, mais elle s'avère extrêmement utile pour faire plusieurs opérations en une seule expression tout en évitant les boucles. La syntaxe de outer est :

outer(X, Y, FUN)

Le résultat est l'application la fonction FUN ("*" par défaut) entre chacun des éléments de X et chacun des éléments de Y, autrement dit

FUN(X[i], Y[j])

pour toutes les valeurs des indices i et j.

- La dimension du résultat est par conséquent c(dim(X), dim(Y)).

- Par exemple, le résultat du produit extérieur entre deux vecteurs est une matrice contenant tous les produits entre les éléments des deux vecteurs :

> outer(c(1, 2, 5), c(2, 3, 6))
[,1] [,2] [,3]
[1,]    2    3    6
[2,]    4    6   12
[3,]   10   15   30

- Lorsque FUN est un opérateur arithmétique du tableau 3.1, on place le symbole entre guillemets : "*", "+", "<=", etc.

- L'opérateur %o% est un raccourci de outer(X, Y, "*").

Nous avons trouvé ici la fonction clé pour résoudre notre problème sans faire de boucles. La fonction outer permet de calculer facilement tous les résultats possibles de la somme de deux dés :

1. > x <- 1:8 # résultats possibles du premier dé

2. > y <- 1:6 # résultats possibles du second dé

3. > outer(x, y, "+") # sommes de toutes les combinaisons

4. [,1] [,2] [,3] [,4] [,5] [,6]

5. [1,]    2    3    4    5    6    7

6. [2,]    3    4    5    6    7    8

7. [3,]    4    5    6    7    8    9

8. [4,]    5    6    7    8    9   10

9. [5,]    6    7    8    9   10   11

10. [6,]    7    8    9   10   11   12

11. [7,]    8    9   10   11   12   13

12. [8,]    9   10   11   12   13   14

3-5 - Structures de contrôle

Les structures de contrôle sont des commandes qui permettent de déterminer le flux d'exécution d'un programme : choix entre des blocs de code, répétition de commandes ou sortie forcée.

On se contente, ici, de mentionner les structures de contrôle disponibles en R. La section 3.7 fournit des exemples d'utilisation.

3-5-1 - Exécution conditionnelle

if (condition) branche.vrai else branche.faux

Si condition est vraie, branche.vrai est exécutée, sinon ce sera branche.faux. Dans le cas où l'une ou l'autre de branche.vrai ou branche.faux comporte plus d'une expression, regrouper celles-ci dans des accolades { }.

ifelse(condition, expression.vrai, expression.faux)

Fonction vectorielle qui retourne un vecteur de la même longueur que condition formé ainsi : pour chaque élément TRUE de condition on choisit l'élément correspondant de expression.vrai et pour chaque élément FALSE on choisit l'élément correspondant de expression.faux. L'utilisation n'est pas très intuitive, alors examiner attentivement les exemples de la rubrique d'aide.

switch(test, cas.1 = action.1, cas.2 = action.2, ...)

Structure utilisée plutôt rarement. Consulter la rubrique d'aide au besoin.

3-5-2 - Boucles

Les boucles sont et doivent être utilisées avec parcimonie en R, car elles sont généralement inefficaces. Dans la majeure partie des cas, il est possible de vectoriser les calculs pour éviter les boucles explicites, ou encore de s'en remettre aux fonctions outer, apply, lapplysapply et mapply (section 6.2) pour réaliser les boucles de manière plus efficace.

for (variable in suite) expression

Exécuter expression successivement pour chaque valeur de variable contenue dans suite. Encore ici, on groupera les expressions dans des accolades { }. À noter que suite n'a pas à être composée de nombres consécutifs, ni même de nombres, en fait.

while (condition) expression

Exécuter expression tant que condition est vraie. Si condition est fausse lors de l'entrée dans la boucle, celle-ci n'est pas exécutée. Une boucle while n'est par conséquent pas nécessairement toujours exécutée.

repeat expression

Répéter expression. Cette dernière devra comporter un test d'arrêt qui utilisera la commande break. Une boucle repeat est toujours exécutée au moins une fois.

break

Sortie immédiate d'une boucle for, while ou repeat.

next

Passage immédiat à la prochaine itération d'une boucle for, while ou repeat.

3-6 - Fonctions additionnelles

La bibliothèque des fonctions internes de R est divisée en ensembles de fonctions et de jeux de données apparentés nommés packages (terme que l'équipe de traduction française de R a choisi de conserver tel quel). On démarrage, R charge automatiquement quelques packages de la bibliothèque, ceux contenant les fonctions les plus fréquemment utilisées. On peut voir la liste des packages déjà en mémoire avec la fonction search et le contenu de toute la bibliothèque avec la fonction library (résultat non montré ici) :

1. > search()

2. [1]       ".GlobalEnv"     "package:stats"

3. [3] "package:graphics" "package:grDevices"

4. [5]    "package:utils"  "package:datasets"

5. [7]  "package:methods"         "Autoloads"

6. [9]     "package:base"

Une des grandes forces de R est la facilité avec laquelle on peut ajouter des fonctionnalités au système par le biais de packages externes. Dès les débuts de R, les développeurs et utilisateurs ont mis sur pied le dépôt central de packages Comprehensive R Archive Network (CRAN ; http://cran.r-project.org). Ce site compte aujourd'hui des milliers d'extensions et le nombre ne cesse de croître.

Le système R rend simple le téléchargement et l'installation de nouveaux packages avec la fonction install.packages. L'annexe D explique plus en détail comment gérer sa bibliothèque personnelle et installer des packages externes.

Solution du problème

Nous savons déjà que l'expression

> m <- outer(x, y, "+")

calcule tous les résultats de la somme du lancer des deux dés. Il faut maintenant extraire de cette matrice les résultats supérieurs à 7 et en faire la moyenne.
Par la propriété d'opération élément par élément, la comparaison

> m > 7
[,1]  [,2]  [,3]  [,4]  [,5]  [,6]
[1,] FALSE FALSE FALSE FALSE FALSE FALSE
[2,] FALSE FALSE FALSE FALSE FALSE  TRUE
[3,] FALSE FALSE FALSE FALSE  TRUE  TRUE
[4,] FALSE FALSE FALSE  TRUE  TRUE  TRUE
[5,] FALSE FALSE  TRUE  TRUE  TRUE  TRUE
[6,] FALSE  TRUE  TRUE  TRUE  TRUE  TRUE
[7,]  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE
[8,]  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE

retourne une matrice booléenne de la même dimension que m. Pour extraire les résultats de m supérieurs à 7, il suffit d'indicer m avec la matrice booléenne ci-dessus :

> m[m > 7]
[1]  8  9  8  9 10 8  9 10 11  8  9 10 11 12 8 9
[17] 10 11 12 13  8 9 10 11 12 13 14

Le résultat recherché est donc

> mean(m[m > 7])
[1] 10.07407

Le tout peut s'écrire en une seule expression — quelque peu alambiquée — ainsi :

> mean((m <- outer(x, y, "+"))[m > 7])
[1] 10.07407

3-7 - Exemples

1. ###

2. ### OPÉRATIONS ARITHMÉTIQUES

3. ###

4. ## L'arithmétique vectorielle caractéristique du langage R

5. ## rend très simple et intuitif de faire des opérations

6. ## mathématiques courantes. Là où plusieurs langages de

7. ## programmation exigent des boucles, R fait le calcul

8. ## directement. En effet, les règles de l'arithmétique en R

9. ## sont globalement les mêmes qu'en algèbre vectorielle et

10. ## matricielle.

11. 5 * c(2, 3, 8, 10)          # multiplication par une constante

12. c(2, 6, 8) + c(1, 4, 9)     # addition de deux vecteurs

13. c(0, 3, -1, 4)^2            # élévation à une puissance

15. ## Dans les règles de l'arithmétique vectorielle, les

16. ## longueurs des vecteurs doivent toujours concorder. R permet

17. ## plus de flexibilité en recyclant les vecteurs les plus

18. ## courts dans une opération. Il n'y a donc à peu près jamais

19. ## d'erreurs de longueur en R! C'est une arme à deux

20. ## tranchants: le recyclage des vecteurs facilite le codage,

21. ## mais peut aussi résulter en des réponses complètement

22. ## erronées sans que le système ne détecte d'erreur.

23. 8 + 1:10                    # 8 est recyclé 10 fois

24. c(2, 5) * 1:10              # c(2, 5) est recyclé 5 fois

25. c(-2, 3, -1, 4)^1:4         # quatre puissances différentes

27. ## On se rappelle que les matrices (et les tableaux) sont des

28. ## vecteurs. Les règles ci-dessus s'appliquent donc aussi aux

29. ## matrices, ce qui résulte en des opérateurs qui ne sont pas

30. ## définis en algèbre linéaire usuelle.

31. (x <- matrix(1:4, 2))       # matrice 2 x 2

32. (y <- matrix(3:6, 2))       # autre matrice 2 x 2

33. 5 * x                       # multiplication par une constante

34. x + y                       # addition matricielle

35. x * y                       # produit *élément par élément*

36. x %*% y                     # produit matriciel

37. x / y                       # division *élément par élément*

38. x * c(2, 3)                 # produit par colonne

40. ###

41. ### OPÉRATEURS

42. ###

44. ## Seuls les opérateurs %%, %/% et logiques sont illustrés

45. ## ici. Premièrement, l'opérateur modulo retourne le reste

46. ## d'une division.

47. 5 %% 2                      # 5/2 = 2 reste 1

48. 5 %% 1:5                    # remarquer la périodicité

49. 10 %% 1:15                  # x %% y = x si x < y

51. ## Le modulo est pratique dans les boucles, par exemple pour

52. ## afficher un résultat à toutes les n itérations seulement.

53. for (i in 1:50)

54. {

55. ## Affiche la valeur du compteur toutes les 5 itérations.

56. if (0 == i %% 5)

57. print(i)

58. }

60. ##La division entière retourne la partie entière de la

61. ## division d'un nombre par un autre.

62. 5 %/% 1:5

63. 10 %/% 1:15

65. ## Le ET logique est vrai seulement lorsque les deux

66. ## expressions sont vraies.

67. c(TRUE, TRUE, FALSE) & c(TRUE, FALSE, FALSE)

69. ## Le OU logique est faux seulement lorsque les deux

70. ## expressions sont fausses.

71. c(TRUE, TRUE, FALSE) | c(TRUE, FALSE, FALSE)

73. ## La négation logique transforme les vrais en faux et vice

74. ## versa.

75. ! c(TRUE, FALSE, FALSE, TRUE)

77. ## On peut utiliser les opérateurs logiques &, | et !

78. ## directement avec des nombres. Dans ce cas, le nombre zéro

79. ## est traité comme FALSE et tous les autres nombres comme

80. ## TRUE.

81. 0:5 & 5:0

82. 0:5 | 5:0

83. !0:5

85. ## Ainsi, dans une expression conditionnelle, inutile de

86. ## vérifier si, par exemple, un nombre est égal à zéro. On

87. ## peut utiliser le nombre directement et sauver des

88. ## opérations de comparaison qui peuvent devenir coûteuses en

89. ## temps de calcul.

90. x <- 1                      # valeur quelconque

91. if (x != 0) x + 1           # TRUE pour tout x != 0

92. if (x) x + 1                # tout à fait équivalent!

94. ## L'exemple de boucle ci-dessus peut donc être légèrement

95. ## modifié.

96. for (i in 1:50)

97. {

98. ## Affiche la valeur du compteur toutes les 5 itérations.

99. if (!i %% 5)

100. print (i)

101. }

103. ## Dans les calculs numériques, TRUE vaut 1 et FALSE vaut 0.

104. a <- c("Impair", "Pair")

105. x <- c(2, 3, 6, 8, 9, 11, 12)

106. x %% 2

107. (!x %% 2) + 1

108. a[(!x %% 2) + 1]

110. ## Un mot en terminant sur l'opérateur '=='. C'est l'opérateur

111. ## à utiliser pour vérifier si deux valeurs sont égales, et

112. ## non '='. C'est là une erreur commune --- et qui peut être

113. ## difficile à détecter --- lorsque l'on programme en R.

114. 5 = 2                       # erreur de syntaxe

115. 5 == 2                      # comparaison

117. ###

118. ### APPELS DE FONCTIONS

119. ###

121. ## Les invocations de la fonction 'matrix' ci-dessous sont

122. ## toutes équivalentes. On remarquera, entre autres, comment

123. ## les arguments sont spécifiés (par nom ou par position).

124. matrix(1:12, 3, 4)

125. matrix(1:12, ncol = 4, nrow = 3)

126. matrix(nrow = 3, ncol = 4, data = 1:12)

127. matrix(nrow = 3, ncol = 4, byrow = FALSE, 1:12)

128. matrix(nrow = 3, ncol = 4, 1:12, FALSE)

130. ###

131. ### QUELQUES FONCTIONS UTILES

132. ###

134. ## MANIPULATION DE VECTEURS

135. x <- c(50, 30, 10, 20, 60, 30, 20, 40) # vecteur non ordonné

137. ## Séquences de nombres.

138. seq(from = 1, to = 10)          # équivalent à 1:10

139. seq_len(10)                     # plus rapide que 'seq'

140. seq(-10, 10, length = 50)       # incrément automatique

141. seq(-2, by = 0.5, along = x)    # même longueur que 'x'

142. seq_along(x)                    # plus rapide que 'seq'

144. ## Répétition de nombres ou de vecteurs complets.

145. rep(1, 10)                  # utilisation de base

146. rep(x, 2)                   # répéter un vecteur

147. rep(x, times = 2, each = 4) # combinaison des arguments

148. rep(x, times = 1:8)         # nombre de répétitions différent

149. # pour chaque élément de 'x'

151. ## Classement en ordre croissant ou décroissant.

152. sort(x)                     # classement en ordre croissant

153. sort(x, decr = TRUE)        # classement en ordre décroissant

154. sort(c("abc", "B", "Aunt", "Jemima")) # chaînes de caractères

155. sort(c(TRUE, FALSE))        # FALSE vient avant TRUE

157. ## La fonction 'order' retourne la position, dans le vecteur

158. ## donné en argument, du premier élément selon l'ordre

159. ## croissant, puis du deuxième, etc. Autrement dit, on obtient

160. ## l'ordre dans lequel il faut extraire les données du vecteur

161. ## pour les obtenir en ordre croissant.

162. order(x)                    # regarder dans le blanc des yeux

163. x[order(x)]                 # équivalent à 'sort(x)'

165. ## Rang des éléments d'un vecteur dans l'ordre croissant.

166. rank(x)                     # rang des éléments de 'x'

168. ## Renverser l'ordre d'un vecteur.

169. rev(x)

171. ## Extraction ou suppression en tête ou en queue de vecteur.

172. head(x, 3)                  # trois premiers éléments

173. head(x, -2)                 # tous sauf les deux derniers

174. tail(x, 3)                  # trois derniers éléments

175. tail(x, -2)                 # tous sauf les deux premiers

177. ## Expressions équivalentes sans 'head' et 'tail'

178. x[1:3]                      # trois premiers éléments

179. x[1:(length(x) - 2)]        # tous sauf les deux derniers

180. x[(length(x)-2):length(x)]  # trois derniers éléments

181. rev(rev(x)[1:3])            # avec petits vecteurs seulement

182. x[c(-1, -2)]                # tous sauf les deux premiers

184. ## Seulement les éléments différents d'un vecteur.

185. unique(x)

187. ## RECHERCHE D'ÉLÉMENTS DANS UN VECTEUR

188. which(x >= 30)              # positions des éléments >= 30

189. which.min(x)                # position du minimum

190. which.max(x)                # position du maximum

191. match(20, x)                # position du premier 20 dans 'x'

192. match(c(20, 30), x)         # aussi pour plusieurs valeurs

193. 60 %in% x                   # 60 appartient à 'x'

194. 70 %in% x                   # 70 n'appartient pas à 'x'

196. ## ARRONDI

197. (x <- c(-21.2, -pi, -1.5, -0.2, 0, 0.2, 1.7823, 315))

198. round(x)                    # arrondi à l'entier

199. round(x, 2)                 # arrondi à la seconde décimale

200. round(x, -1)                # arrondi aux dizaines

201. ceiling(x)                  # plus petit entier supérieur

202. floor(x)                    # plus grand entier inférieur

203. trunc(x)                    # troncature des décimales

205. ## SOMMAIRES ET STATISTIQUES DESCRIPTIVES

206. sum(x)                      # somme des éléments

207. prod(x)                     # produit des éléments

208. diff(x)                     # x[2] - x[1], x[3] - x[2], etc.

209. mean(x)                     # moyenne des éléments

210. mean(x, trim = 0.125)       # moyenne sans minimum et maximum

211. var(x)                      # variance (sans biais)

212. (length(x) - 1)/length(x) * var(x) # variance biaisée

213. sd(x)                       # écart type

214. max(x)                      # maximum

215. min(x)                      # minimum

216. range(x)                    # c(min(x), max(x))

217. diff(range(x))              # étendue de 'x'

218. median(x)                   # médiane (50e quantile) empirique

219. quantile(x)                 # quantiles empiriques

220. quantile(x, 1:10/10)        # on peut spécifier les quantiles

221. summary(x)                  # plusieurs des résultats ci-dessus

223. ## SOMMAIRES CUMULATIFS ET COMPARAISONS ÉLÉMENT PAR ÉLÉMENT

224. (x <- sample(1:20, 6))

225. (y <- sample(1:20, 6))

226. cumsum(x)                   # somme cumulative de 'x'

227. cumprod(y)                  # produit cumulatif de 'y'

228. rev(cumprod(rev(y)))        # produit cumulatif renversé

229. cummin(x)                   # minimum cumulatif

230. cummax(y)                   # maximum cumulatif

231. pmin(x, y)                  # minimum élément par élément

232. pmax(x, y)                  # maximum élément par élément

234. ## OPÉRATIONS SUR LES MATRICES

235. (A <- sample(1:10, 16, replace = TRUE)) # avec remise

236. dim(A) <- c(4, 4)           # conversion en une matrice 4 x 4

237. b <- c(10, 5, 3, 1)         # un vecteur quelconque

238. A                           # la matrice 'A'

239. t(A)                        # sa transposée

240. solve(A)                    # son inverse

241. solve(A, b)                 # la solution de Ax = b

242. A %*% solve(A, b)           # vérification de la réponse

243. diag(A)                     # extraction de la diagonale de 'A'

244. diag(b)                     # matrice diagonale formée avec 'b'

245. diag(4)                     # matrice identité 4 x 4

246. (A <- cbind(A, b))          # matrice 4 x 5

247. nrow(A)                     # nombre de lignes de 'A'

248. ncol(A)                     # nombre de colonnes de 'A'

249. rowSums(A)                  # sommes par ligne

250. colSums(A)                  # sommes par colonne

251. apply(A, 1, sum)            # équivalent à 'rowSums(A)'

252. apply(A, 2, sum)            # équivalent à 'colSums(A)'

253. apply(A, 1, prod)           # produit par ligne avec 'apply'

255. ## PRODUIT EXTÉRIEUR

256. x <- c(1, 2, 4, 7, 10, 12)

257. y <- c(2, 3, 6, 7, 9, 11)

258. outer(x, y)                 # produit extérieur

259. x %o% y                     # équivalent plus court

260. outer(x, y, "+")            # « somme extérieure »

261. outer(x, y, "<=")           # toutes les comparaisons possibles

262. outer(x, y, pmax)           # idem

264. ###

265. ### STRUCTURES DE CONTRÔLE

266. ###

268. ## Pour illustrer les structures de contrôle, on a recours à

269. ## un petit exemple tout à fait artificiel : un vecteur est

270. ## rempli des nombres de 1 à 100, à l'exception des multiples

271. ## de 10. Ces derniers sont affichés à l'écran.

272. ##

273. ## À noter qu'il est possible --- et plus efficace --- de

274. ## créer le vecteur sans avoir recours à des boucles.

275. (1:100)[-((1:10) * 10)]              # sans boucle !

276. rep(1:9, 10) + rep(0:9*10, each = 9) # une autre façon !

278. ## Bon, l'exemple proprement dit...

279. x <- numeric(0)             # initialisation du contenant 'x'

280. j <- 0                      # compteur pour la boucle

281. for (i in 1:100)

282. {

283. if (i %% 10)            # si i n'est pas un multiple de 10

284. x[j <- j + 1] <- i  #     stocker sa valeur dans 'x'

285. else                    # sinon

286. print(i)            #     afficher la valeur à l'écran

287. }

288. x                           # vérification

290. ## Même chose que ci-dessus, mais sans le compteur 'j' et les

291. ## valeurs manquantes aux positions 10, 20, ... , 100 sont

292. ## éliminées à la sortie de la boucle.

293. x <- numeric(0)

294. for (i in 1:100)

295. {

296. if (i %% 10)

297. x[i] <- i

298. else

299. print(i)

300. }

301. x <- x[!is.na(x)]

302. x

304. ## On peut refaire l'exemple avec une boucle 'while', mais

305. ## cette structure n'est pas naturelle ici puisque l'on sait

306. ## d'avance qu'il faudra faire la boucle exactement 100

307. ## fois. Le 'while' est plutôt utilisé lorsque le nombre de

308. ## répétitions est inconnu. De plus, une boucle 'while' n'est

309. ## pas nécessairement exécutée puisque le critère d'arrêt est

310. ## évalué dès l'entrée dans la boucle.

311. x <- numeric(0)

312. j <- 0

313. i <- 1                      # pour entrer dans la boucle [*]

314. while (i <= 100)

315. {

316. if (i %% 10)

317. x[j <- j + 1] <- i

318. else

319. print(i)

320. i <- i + 1              # incrémenter le compteur!

321. }

322. x

324. ## La remarque faite au sujet de la boucle 'while' s'applique

325. ## aussi à la boucle 'repeat'. Par contre, le critère d'arrêt

326. ## de la boucle 'repeat' étant évalué à la toute fin, la

327. ## boucle est exécutée au moins une fois. S'il faut faire la

328. ## manœuvre marquée [*] ci-dessus pour s'assurer qu'une

329. ## boucle 'while' est exécutée au moins une fois... c'est

330. ## qu'il faut utiliser 'repeat'.

331. x <- numeric(0)

332. j <- 0

333. i <- 1

334. repeat

335. {

336. if (i %% 10)

337. x[j <- j + 1] <- i

338. else

339. print(i)

340. if (100 < (i <- i + 1)) # incrément et critère d'arrêt

341. break

342. }

343. x

345. ###

346. ### FONCTIONS ADDITIONNELLES

347. ###

349. ## La fonction 'search' retourne la liste des environnements

350. ## dans lesquels R va chercher un objet (en particulier une

351. ## fonction). '.GlobalEnv' est l'environnement de travail.

352. search()

354. ## Liste de tous les packages installés sur votre système.

355. library()

357. ## Chargement du package 'MASS', qui contient plusieurs

358. ## fonctions statistiques très utiles.

359. library("MASS")

3-8 - Exercices

3.1 (solution) À l'aide des fonctions rep, seq et c seulement, générer les séquences suivantes.

1. 0 6 0 6 0 6

2. 1 4 7 10

3. 1 2 3 1 2 3 1 2 3 1 2 3

4. 1 2 2 3 3 3

5. 1 1 1 2 2 3

6. 1 5.5 10

7. 1 1 1 1 2 2 2 2 3 3 3 3

3.2 (solution) Générer les suites de nombres suivantes à l'aide des fonctions : et rep seulement, donc sans utiliser la fonction seq.

1. 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2

2. 1 3 5 7 9 11 13 15 17 19

3. -2 -1 0 1 2 -2 -1 0 1 2

4. -2 -2 -1 -1 0 0 1 1 2 2

5. 10 20 30 40 50 60 70 80 90 100

3.3 (solution) À l'aide de la commande apply, écrire des expressions R qui remplaceraient les fonctions suivantes.

1. rowSums

2. colSums

3. rowMeans

4. colMeans

3.4 (solution) Sans utiliser les fonctions factorial, lfactorial, gamma ou lgamma, générer la séquence 1! , 2! … , 10!.

3.5 (solution) Trouver une relation entre x, y, x %% y (modulo) et x %/% y (division entière), où y != 0.

3.6 (solution) (Simuler un échantillon ) avec la fonction sample.

Écrire une expression R permettant d'obtenir ou de calculer chacun des résultats demandés ci-dessous.

1. Les cinq premiers éléments de l'échantillon.

2. La valeur maximale de l'échantillon.

3. La moyenne des cinq premiers éléments de l'échantillon.

4. La moyenne des cinq derniers éléments de l'échantillon.

3.7 (solution)

1. Trouver une formule pour calculer la position, dans le vecteur sous-jacent, de l'élément d'une matrice remplie par colonne.

2. Répéter la partie a) pour l'élément d'un tableau .

3.8 (solution) Simuler une matrice mat 10×7, puis écrire des expressions R permettant d'effectuer les tâches demandées ci-dessous.

1. Calculer la somme des éléments de chacune des lignes de la matrice.

2. Calculer la moyenne des éléments de chacune des colonnes de la matrice.

3. Calculer la valeur maximale de la sous-matrice formée par les trois premières lignes et les trois premières colonnes.

4. Extraire toutes les lignes de la matrice dont la moyenne des éléments est supérieure à 7.

3.9 (solution) On vous donne la liste et la date des 31 meilleurs temps enregistrés au 100 mètres homme entre 1964 et 2005 :

1. > temps <- c(10.06, 10.03, 10.02,  9.95, 10.04, 10.07,

2. +            10.08, 10.05,  9.98, 10.09, 10.01, 10.00,

3. +             9.97,  9.93,  9.96,  9.99,  9.92,  9.94,

4. +             9.90,  9.86,  9.88,  9.87,  9.85,  9.91,

5. +             9.84,  9.89,  9.79,  9.80,  9.82,  9.78,

6. +             9.77)

7. > names(temps) <- c("1964-10-15", "1968-06-20",

8. +     "1968-10-13", "1968-10-14", "1968-10-14",

9. +     "1968-10-14", "1968-10-14", "1975-08-20",

10. +     "1977-08-11", "1978-07-30", "1979-09-04",

11. +     "1981-05-16", "1983-05-14", "1983-07-03",

12. +     "1984-05-05", "1984-05-06", "1988-09-24",

13. +     "1989-06-16", "1991-06-14", "1991-08-25",

14. +     "1991-08-25", "1993-08-15", "1994-07-06",

15. +     "1994-08-23", "1996-07-27", "1996-07-27",

16. +     "1999-06-16", "1999-08-22", "2001-08-05",

17. +     "2002-09-14", "2005-06-14")

Extraire de ce vecteur les records du monde seulement, c'est-à-dire la première fois que chaque temps a été réalisé.

4 - Exemples résolus

Objectifs du chapitre

- Mettre en pratique les connaissances acquises dans les chapitres précédents.

- Tirer profit de l'arithmétique vectorielle de R pour effectuer des calculs complexes sans boucles.

- Utiliser l'initialisation de vecteurs et leur indiçage de manière à réduire le temps de calcul.

Ce chapitre propose de faire le point sur les concepts étudiés jusqu'à maintenant par le biais de quelques exemples résolus. On y met particulièrement en évidence les avantages de l'approche vectorielle du langage R.

Les exemples font appel à quelques connaissances de base en mathématiques financières et en théorie des probabilités.

Énoncé du problème

On modélise la distance atteinte par le lancer du poids d'un adulte à l'aide d'une variable aléatoire . On rappelle que la fonction de répartition de cette loi est :

En compétition, les juges inscrivent un résultat seulement si le lancer dépasse la distance La probabilité qu'un lancer dépasse la distance sachant que ce lancer a dépassé la distance est de 0 ,75. Calculer la constante pour laquelle la relation précédente est vraie, autrement dit la valeur de tel que

4-1 - Calcul de valeurs actuelles

La valeur actuelle d'une série de paiements à la fin des années est

où est le taux d'intérêt effectif annuellement durant l'année . Lorsque le taux d'intérêt est constant au cours des années, cette formule se simplifie en

Un prêt est remboursé par une série de cinq paiements, le premier étant dû dans un an. On doit trouver le montant du prêt pour chacune des hypothèses ci-dessous.

1. Paiement annuel de 1000, taux d'intérêt de 6 % effectif annuellement. Avec un paiement annuel et un taux d'intérêt constants, on utilise la formule (4.3) avec  :

> 1000 * sum((1 + 0.06)^(-(1:5)))
[1] 4212.364

Remarquer comme l'expression R se lit exactement comme la formule mathématique. De plus, le terme constant 1000 est sorti de la somme pour réduire le nombre de multiplications de cinq à une seule.

2. Paiements annuels de 500, 800, 900, 750 et 1000, taux d'intérêt de 6 % effectif annuellement.

Les paiements annuels sont différents, mais le taux d'intérêt est toujours le même. La formule (4.3) s'applique donc directement :

> sum(c(500, 800, 900, 750, 1000) * (1 + 0.06)^(-(1:5)))
[1] 3280.681

3. Paiements annuels de 500, 800, 900, 750 et 1000, taux d'intérêt de 5 %, 6 %, 5.5 %, 6.5 % et 7 % effectifs annuellement.

Avec différents paiements annuels et des taux d'intérêt différents, il faut employer la formule (4.2). On obtient le résultat voulu sans aucune boucle en réalisant le produit cumulatif des taux d'intérêt avec la fonction cumprod :

> sum(c(500, 800, 900, 750, 1000) /
+     cumprod(1 + c(0.05, 0.06, 0.055, 0.065, 0.07)))
[1] 3308.521

4-2 - Fonctions de masse de probabilité

On souhaite calculer toutes ou la majeure partie des probabilités de deux lois de probabilité, puis vérifier que la somme des probabilités est bien égale à 1.

Cet exemple est quelque peu artificiel dans la mesure où il existe dans R des fonctions internes pour calculer les principales caractéristiques des lois de probabilité les plus usuelles. Nous utiliserons d'ailleurs ces fonctions pour vérifier nos calculs.

1. Calculer toutes les masses de probabilité de la distribution binomiale pour des valeurs des paramètres et quelconques. La fonction de masse de probabilité de la binomiale est

Soit n = 10 et p = 0.8. Les coefficients binomiaux sont calculés avec la fonction choose :

1. > n <- 10

2. > p <- 0.8

3. > x <- 0:n

4. > choose(n, x) * p^x * (1 - p)^rev(x)

5. [1]  0.0000001024 0.0000040960 0.0000737280

6. [4]  0.0007864320 0.0055050240 0.0264241152

7. [7]  0.0880803840 0.2013265920 0.3019898880

8. [10] 0.2684354560 0.1073741824

On vérifie les réponses obtenues avec la fonction interne dbinom :

1. > dbinom(x, n, prob = 0.8)

2. [1]  0.0000001024 0.0000040960 0.0000737280

3. [4]  0.0007864320 0.0055050240 0.0264241152

4. [7]  0.0880803840 0.2013265920 0.3019898880

5. [10] 0.2684354560 0.1073741824

On vérifie enfin que les probabilités somment à 1 :

> sum(choose(n, x) * p^x * (1 - p)^rev(x))
[1] 1

2. Calculer la majeure partie des masses de probabilité de la distribution de Poisson, dont la fonction de masse de probabilité est

où .

La loi de Poisson ayant un support infini, on calcule les probabilités en seulement avec . On calcule les factorielles avec la fonction factorial. On notera au passage que factorial(x) == gamma(x + 1), où la fonction R gamma calcule les valeurs de la fonction mathématique du même nom

avec . Pour entier, on a donc .

1. > lambda <- 5

2. > x <- 0:10

3. > exp(-lambda) * (lambda^x / factorial(x))

4. [1]  0.006737947 0.033689735 0.084224337

5. [4]  0.140373896 0.175467370 0.175467370

6. [7]  0.146222808 0.104444863 0.065278039

7. [10] 0.036265577 0.018132789

Vérification avec la fonction interne dpois :

1. > dpois(x, lambda)

2. [1]  0.006737947 0.033689735 0.084224337

3. [4]  0.140373896 0.175467370 0.175467370

4. [7]  0.146222808 0.104444863 0.065278039

5. [10] 0.036265577 0.018132789

Pour vérifier que les probabilités somment à 1, il faudra d'abord tronquer le support infini de la Poisson à une « grande » valeur. Ici, 170 est suffisamment éloigné de la moyenne de la distribution, 5. Remarquer que le produit par est placé à l'extérieur de la somme pour ainsi faire un seul produit plutôt que 171.

1. > x <- 0:170

2. > exp(-lambda) * sum((lambda^x / factorial(x)))

3. [1] 1

4-3 - Fonction de répartition de la loi gamma

La loi gamma est fréquemment employée pour modéliser des événements ne pouvant prendre que des valeurs positives et pour lesquels les petites valeurs sont plus fréquentes que les grandes. Par exemple, on utilise parfois la loi gamma en sciences actuarielles pour la modélisation des montants de sinistres. Nous utiliserons la paramétrisation où la fonction de densité de probabilité est

où est la fonction gamma définie dans l'exemple précédent.

Il n'existe pas de formule explicite de la fonction de répartition de la loi gamma. Néanmoins, la valeur de la fonction de répartition d'une loi gamma de paramètre entier et peut être obtenue à partir de la formule

1. Évaluer .

Cet exercice est simple puisqu'il s'agit de calculer une seule valeur de la fonction de répartition avec un paramètre fixe. Par une application directe de (4.5), on a :

1. > alpha <- 5

2. > x <- 4

3. > 1 - exp(-x) * sum(x^(0:(alpha - 1))/gamma(1:alpha))

4. [1] 0.3711631

Vérification avec la fonction interne pgamma :

> pgamma(x, alpha)
[1] 0.3711631

On peut éviter de générer essentiellement la même suite de nombres à deux reprises en ayant recours à une variable intermédiaire. Au risque de rendre le code un peu moins lisible (mais plus compact !), l'affectation et le calcul final peuvent même se faire dans une seule expression.

> 1 - exp(-x) * sum(x^(-1 + (j <- 1:alpha))/gamma(j))
[1] 0.3711631

2. Évaluer pour en une seule expression.

Cet exercice est beaucoup plus compliqué qu'il n'y paraît au premier abord. Ici, la valeur de demeure fixe, mais on doit calculer, en une seule expression, la valeur de la fonction de répartition en plusieurs points. Or, cela exige de faire d'un coup le calcul pour plusieurs valeurs de et plusieurs valeurs de . C'est un travail pour la fonction outer :

1. > x <- 2:10

2. > 1 - exp(-x) * colSums(t(outer(x, 0:(alpha-1), "^")) /

3. + gamma(1:alpha))

4. [1] 0.05265302 0.18473676 0.37116306 0.55950671

5. [5] 0.71494350 0.82700839 0.90036760 0.94503636

6. [9] 0.97074731

Vérification avec la fonction interne pgamma :

1. > pgamma(x, alpha)

2. [1] 0.05265302 0.18473676 0.37116306 0.55950671

3. [5] 0.71494350 0.82700839 0.90036760 0.94503636

4. [9] 0.97074731

Il est laissé en exercice de déterminer pourquoi la transposée est nécessaire dans l'expression ci-dessus. Exécuter l'expression étape par étape, de l'intérieur vers l'extérieur, pour mieux comprendre comment on arrive à faire le calcul en (4.5).

Simplifions déjà l'équation à résoudre. On a

On cherche donc à résoudre pour l'équation

Celle-ci n'admet pas de solution explicite.

4-4 - Algorithme du point fixe

Trouver la racine d'une fonction — c'est-à-dire le point où — est un problème classique en mathématiques. Très souvent, il est possible de reformuler le problème de façon à plutôt chercher le point où . La solution d'un tel problème est appelée point fixe.

L'algorithme du calcul numérique du point fixe d'une fonction est très simple :

1. Choisir une valeur de départ  ;

2. Calculer pour  ;

3. Répéter l'étape 2 jusqu'à ce que ou .

On doit trouver, à l'aide de la méthode du point fixe, la valeur de telle que

c'est-à-dire le taux de rendement d'une série de 10 versements de 1 pour laquelle on a payé un montant de 8,21.

Puisque, d'une part, nous ignorons combien de fois la procédure itérative devra être répétée et que, d'autre part, il faut exécuter la procédure au moins une fois, le choix logique pour la structure de contrôle à utiliser dans cette procédure itérative est repeat. De plus, il faut comparer deux valeurs successives du taux d'intérêt, nous devrons donc avoir recours à deux variables.

On a :

1. > i <- 0.05

2. > repeat

3. + {

4. + it <- i

5. + i <- (1 - (1 + it)^(-10))/8.21

6. + if (abs(i - it)/it < 1E-10)

7. + break

8. + }

9. > i

10. [1] 0.03756777

Vérification :

> (1 - (1 + i)^(-10))/i
[1] 8.21

Nous verrons au chapitre 5 comment créer une fonction à partir de ce code.

L'algorithme du point fixe peut nous permettre de résoudre numériquement l'équation

Pour utiliser l'algorithme, il faut réécrire l'équation sous une forme où l'inconnue se trouve de part et d'autre de l'égalité. En isolant, pour des raisons numériques, le se trouvant au dénominateur de l'équation ci-dessus, on obtient :

4-5 - Suite de Fibonacci

La suite de Fibonacci est une suite de nombres entiers très connue. Les deux premiers termes de la suite sont 0 et 1 et tous les autres sont la somme des deux termes précédents. Mathématiquement, les valeurs de la suite de Fibonacci sont données par la fonction

Le quotient de deux termes successifs converge vers , le nombre d'or.

On veut calculer les premiers termes de la suite de Fibonacci. Ce problème étant intrinsèquement récursif, nous devons utiliser une boucle.

Voici une première solution pour  :

1. > n <- 10

2. > x <- c(0, 1)

3. > for (i in 3:n) x[i] <- x[i - 1] + x[i - 2]

4. > x

5. [1] 0 1 1 2 3 5 8 13 21 34

La procédure ci-dessus a un gros défaut : la taille de l'objet x est constamment augmentée pour stocker une nouvelle valeur de la suite. Tentons une analogie alimentaire pour cette manière de procéder. Pour ranger des biscuits frais sortis du four, on prend un premier biscuit et on le range dans un plat ne pouvant contenir qu'un seul biscuit. Arrivé au second biscuit, on constate que le contenant n'est pas assez grand, alors on sort un plat pouvant contenir deux biscuits, on change le premier biscuit de plat et on y range aussi le second biscuit. Arrivé au troisième biscuit, le petit manège recommence, et ainsi de suite jusqu'à ce que le plateau de biscuits soit épuisé. C'est ce que nous nommerons, non sans un sourire en coin, le Syndrôme de la plaque à biscuits™.

Le manège décrit ci-dessus se reproduit à l'identique dans la mémoire de l'ordinateur, l'odeur des bons biscuits chauds en moins. En effet, l'ordinateur doit constamment allouer de la nouvelle mémoire et déplacer les termes déjà sauvegardés au fur et à mesure que le vecteur x grandit. On aura compris qu'une telle façon de faire est à éviter absolument lorsque c'est possible — et ça l'est la plupart du temps.

Quand on sait quelle sera la longueur d'un objet, comme c'est le cas dans cet exemple, il vaut mieux créer un contenant vide de la bonne longueur et le remplir par la suite. Cela nous donne une autre façon de calculer la suite de Fibonacci :

1. > n <- 10

2. > x <- numeric(n) # création du contenant

3. > x[2] <- 1 # x[1] vaut déjà 0

4. > for (i in 3:n) x[i] <- x[i - 1] + x[i - 2]

5. > x

6. [1] 0 1 1 2 3 5 8 13 21 34

Dans le code informatique du chapitre 5, nous composerons des fonctions avec ces deux solutions et nous comparerons les temps de calcul.

Solution du problème

On résout le problème par la méthode du point fixe en procédant de la même façon qu'à la section 4.4. On prend toutefois soin de modifier la deuxième ligne à l'intérieur de la boucle repeat :

1. > p <- 0.75

2. > alpha <- 20

3. > lambda <- 75

4. > d <- 2

5. > repeat

6. + {

7. +      dt <- d

8. +      d <- log(lambda * (1 - p^(1/alpha)) + d) - log(p)/alpha

9. +      if (abs(d - dt)/dt < 1E-10)

10. +          break

11. + }

12. > d

13. [1] 0.4012012

Vérifions que l'équation de départ (4.1) est vraie :

> ((lambda + d) / (lambda + exp(d)))^alpha
[1] 0.75

Le package actuar fournit une fonction ppareto pour calculer la fonction de répartition de la loi de Pareto. (Consulter l'annexe D pour installer le package à partir du site CRAN.) La fonction permet de vérifier directement l'équation de départ :

1. > library(actuar)

2. > ppareto(exp(d), alpha, lambda, lower.tail = FALSE) /

3. +     ppareto(d, alpha, lambda, lower.tail = FALSE)

4. [1] 0.75

4-6 - Exercices

Dans chacun des exercices ci-dessous, écrire une expression R pour faire le calcul demandé. Parce qu'elles ne sont pas nécessaires, il est interdit d'utiliser des boucles.

4.1 (solution) Calculer la valeur actuelle d'une série de paiements fournie dans un vecteur P en utilisant les taux d'intérêt annuels d'un vecteur i.

4.2 (solution) Étant donné un vecteur d'observations et un vecteur de poids correspondants , calculer la moyenne pondérée des observations,

où . Tester l'expression avec les vecteurs de données

1. x = (7 , 13 , 3 , 8 , 12 , 12 , 20 , 11)

et

1. w = (0 ,15 , 0 ,04 , 0 ,05 , 0 ,06 , 0 ,17 , 0 ,16 , 0 ,11 , 0 ,09).

4.3 (solution) Soit un vecteur d'observations . Calculer la moyenne harmonique de ce vecteur, définie comme

Tester l'expression avec les valeurs de l'exercice 4.2.

4.4 (solution) Calculer la fonction de répartition en d'une loi de Poisson avec paramètre , qui est donnée par

où .

4.5 (solution) a) Calculer l'espérance d'une variable aléatoire dont le support est et les probabilités correspondantes sont , dans l'ordre.

b) Calculer la variance de la variable aléatoire définie en a).

4.6 (solution) Calculer le taux d'intérêt nominal composé quatre fois par année, , équivalent à un taux de = 6 % effectif annuellement.

4.7 (solution) La valeur actuelle d'une série de paiements de fin d'année à un taux d'intérêt effectif annuellement est

où . Calculer en une seule expression, toujours sans boucle, un tableau des valeurs actuelles de séries de paiements à chacun des taux d'intérêt effectifs annuellement .

4.8 (solution) Calculer la valeur actuelle d'une annuité croissante de 1 $ payable annuellement en début d'année pendant dix ans si le taux d'actualisation est de 6 %. Cette valeur actuelle est donnée par

toujours avec .

4.9 (solution) Calculer la valeur actuelle de la suite de paiements 1, 2, 2, 3, 3, 3, 4, 4, 4, 4 si les paiements sont effectués en fin d'année et que le taux d'actualisation est de 7 %.

4.10 (solution) Calculer la valeur actuelle de la suite de paiements de l'exercice 4.9 en supposant que le taux d'intérêt d'actualisation alterne successivement entre 5 % et 8 % chaque année, c'est-à-dire que le taux d'intérêt est de 5 %, 8 %, 5 %, 8 %, etc.

5 - Fonctions définies par l'usager

Objectifs du chapitre

- Définir une fonction R, ses divers arguments et, le cas échéant, les valeurs par défaut de ceux-ci.

- Déboguer une fonction R.

- Adopter un style de codage correspondant à la pratique reconnue en R.

La possibilité pour l'usager de définir facilement et rapidement de nouvelles fonctions — et donc des extensions au langage — est une des grandes forces de R. Les fonctions personnelles définies dans l'espace de travail ou dans un package sont traitées par le système exactement comme les fonctions internes.

Ce court chapitre passe en revue la syntaxe et les règles pour créer des fonctions dans R. On discute également brièvement de débogage et de style de codage.

Énoncé du problème

Stock Ticker est un jeu canadien datant de 1937 dans lequel on brasse une série de dés pour simuler les mouvements boursiers et les dividendes de six titres financiers.
Pour les fins de cet exercice, nous ne nous intéresserons qu'aux mouvements boursiers, et ce, pour un seul titre — nous ignorerons donc le versement des dividendes.
La valeur du titre est établie à partir de la valeur au tour de jeu précédent et du résultat du lancer de deux dés.

1. Un premier dé à trois faces détermine la direction du mouvement boursier : hausse, baisse ou nul.

2. Un second dé à trois faces détermine l'amplitude de la hausse ou de la baisse, le cas échéant. Les valeurs possibles sont 5, 10 et 20.

La valeur de départ du titre est 100.
Par exemple, si les résultats des dés au premier tour de jeu sont « baisse » et 20, la valeur du titre après ce tour sera 80. Si, au tour suivant, les résultats des dés sont « nul » et « 5 », la valeur du titre demeurera de 80. On fournit les résultats de 10 lancers des dés sous forme de data frame :

1. > x

2. direction amplitude

3. 1    baisse        20

4. 2    baisse        10

5. 3       nul        20

6. 4       nul        20

7. 5       nul        20

8. 6    hausse        20

9. 7       nul        20

10. 8    hausse         5

11. 9    baisse        10

12. 10    hausse        20

Écrire une fonction valeurs() servant à calculer les valeurs successives du titre.
La fonction prend en arguments x, un data frame contenant les résultats des lancers des dés, ainsi que start, la valeur de départ du titre (100 par défaut). La fonction retourne un vecteur contenant les valeurs successives du titre.
Par exemple, avec les données ci-dessus, on obtiendrait :

> valeurs(x)            # start = 100 par défaut
[1]  80  70  70  70  70  90  90  95  85 105

Le calcul de la valeur du titre étant récursif, l'utilisation d'une boucle est ici inévitable. Comme le nombre de répétitions est connu d'avance (le nombre de valeurs à calculer correspond au nombre de lancers de dés), une boucle for() serait le choix approprié.

5-1 - Définition d'une fonction

On définit une nouvelle fonction avec la syntaxe suivante :

fun <- function(arguments) expression

où

- fun est le nom de la fonction (les règles pour les noms de fonctions étant les mêmes que celles présentées à la section 2.2 pour tout autre objet) ;

- arguments est la liste des arguments, séparés par des virgules ;

- expression constitue le corps de la fonction, soit une expression ou un groupe d'expressions réunies par des accolades.

5-2 - Retourner des résultats

La plupart des fonctions sont écrites dans le but de retourner un résultat. Or, les règles d'interprétation d'un groupe d'expressions présentées à la section 2.1 s'appliquent ici au corps de la fonction.

- Une fonction retourne tout simplement le résultat de la dernière expression du corps de la fonction.

- On évitera donc que la dernière expression soit une affectation, car la fonction ne retournera alors rien et on ne pourra utiliser une construction de la forme x <- f() pour affecter le résultat de la fonction à une variable.

- Si on doit retourner un résultat sans être à la dernière ligne de la fonction (à l'intérieur d'un bloc conditionnel, par exemple), on utilise la fonction return. L'utilisation de return à la toute fin d'une fonction est tout à fait inutile et considérée comme du mauvais style en R.

- Lorsqu'une fonction doit retourner plusieurs résultats, il est en général préférable d'avoir recours à une liste nommée.

5-3 - Variables locales et globales

Comme la majorité des langages de programmation, R comporte des concepts de variable locale et de variable globale.

- Toute variable définie dans une fonction est locale à cette fonction, c'est-à-dire qu'elle :

- n'apparaît pas dans l'espace de travail ;

- n'écrase pas une variable du même nom dans l'espace de travail.

- Il est possible de définir une variable dans l'espace de travail depuis une fonction avec l'opérateur d'affectation <<-. Il est très rare — et généralement non recommandé — de devoir recourir à de telles variables globales.

- On peut définir une fonction à l'intérieur d'une autre fonction. Cette fonction sera locale à la fonction dans laquelle elle est définie.

Le lecteur intéressé à en savoir plus pourra consulter les sections de la documentation de R portant sur la portée lexicale (lexical scoping). C'est un sujet important et intéressant, mais malheureusement trop avancé pour ce document d'introduction à la programmation en R.

5-4 - Exemple de fonction

Le code développé pour l'exemple de point fixe de la section 4.4 peut être intégré dans une fonction ; voir la figure 5.1.

- Le nom de la fonction est fp.

- La fonction compte quatre arguments : k, n, start et TOL.

- Les deux derniers arguments ont respectivement des valeurs par défaut de 0.05 et 10−10.

- La fonction retourne la valeur de la variable i puisque l'on évalue celle-ci à la dernière ligne (ou expression) de la fonction.

La définition de la fonction devra inclure start = 100 comme deuxième argument afin de spécifier que :

1. start est un argument de la fonction ;

2. La valeur par défaut de l'argument est 100.

Fig. 5.1 - Exemple de fonction de point fixe

1. fp <- function(k, n, start = 0.05, TOL = 1E-10)

2. {

3. ## Fonction pour trouver par la méthode du point

4. ## fixe le taux d'intérêt pour lequel une série de

5. ## 'n' paiements vaut 'k'.

6. ##

7. ## ARGUMENTS

8. ##

9. ##     k: la valeur présente des paiements

10. ##     n: le nombre de paiements

11. ## start: point de départ des itérations

12. ##   TOL: niveau de précision souhaité

13. ##

14. ## RETOURNE

15. ##

16. ## Le taux d'intérêt

18. i <- start

19. repeat

20. {

21. it <- i

22. i <- (1 - (1 + it)^(-n))/k

23. if (abs(i - it)/it < TOL)

24. break

25. }

26. i

27. }

5-5 - Fonctions anonymes

Il est parfois utile de définir une fonction sans lui attribuer un nom — d'où la notion de fonction anonyme. Il s'agira en général de fonctions courtes utilisées dans une autre fonction. Par exemple, pour calculer la valeur de xy^2 pour toutes les combinaisons de x et y stockées dans des vecteurs du même nom, on pourrait utiliser la fonction outer ainsi :

1. > x <- 1:3; y <- 4:6

2. > f <- function(x, y) x * y^2

3. > outer(x, y, f)

4. [,1] [,2] [,3]

5. [1,]   16   25   36

6. [2,]   32   50   72

7. [3,]   48   75  108

Cependant, si la fonction f ne sert à rien ultérieurement, on peut se contenter de passer l'objet fonction à outer sans jamais lui attribuer un nom :

1. > outer(x, y, function(x, y) x * y^2)

2. [,1] [,2] [,3]

3. [1,]   16   25   36

4. [2,]   32   50   72

5. [3,]   48   75  108

On a alors utilisé dans outer une fonction anonyme.

5-6 - Débogage de fonctions

Il est assez rare d'arriver à écrire un bout de code sans bogue du premier coup. Par conséquent, qui dit programmation dit séances de débogage.

Les techniques de débogage les plus simples et naïves sont parfois les plus efficaces et certainement les plus faciles à apprendre. Loin d'un traité sur le débogage de code R, nous offrons seulement ici quelques trucs que nous utilisons régulièrement.

- Les erreurs de syntaxe sont les plus fréquentes (en particulier l'oubli de virgules). Lors de la définition d'une fonction, une vérification de la syntaxe est effectuée par l'interprète R. Attention, cependant : une erreur peut prendre sa source plusieurs lignes avant celle que l'interprète pointe comme causant problème.

- Les messages d'erreur de l'interprète ne sont pas toujours d'un grand secours… tant que l'on n'a pas appris à les reconnaître. Un exemple de message d'erreur fréquemment rencontré :

1. valeur manquante là où TRUE / FALSE est requis

Cette erreur provient généralement d'une commande if dont l'argument vaut NA plutôt que TRUE ou FALSE. La raison : des valeurs manquantes se sont faufilées dans les calculs à notre insu jusqu'à l'instruction if, faisant en sorte que l'argument de if vaut NA alors qu'il ne peut être que booléen.

- Lorsqu'une fonction ne retourne pas le résultat attendu, placer des commandes print à l'intérieur de la fonction, de façon à pouvoir suivre les valeurs prises par les différentes variables.

Par exemple, la modification suivante à la boucle de la fonction fp permet d'afficher les valeurs successives de la variable i et de détecter, par exemple, une procédure itérative divergente :

1. repeat

2. {

3. it <- i

4. i <- (1 - (1 + it)^(-n))/k

5. print(i)

6. if (abs((i - it)/it < TOL))

7. break

8. }

- Quand ce qui précède ne fonctionne pas, ne reste souvent qu'à exécuter manuellement la fonction. Pour ce faire, définir dans l'espace de travail tous les arguments de la fonction, puis exécuter le corps de la fonction ligne par ligne. La vérification du résultat de chaque ligne permet généralement de retrouver la ou les expressions qui causent problème.

5-7 - Styles de codage

Si tous conviennent que l'adoption d'un style propre et uniforme favorise le développement et la lecture de code, il existe plusieurs chapelles dans le monde des programmeurs quant à la « bonne façon » de présenter et, surtout, d'indenter le code informatique.

Par exemple, Emacs reconnaît et supporte les styles de codage suivants, entre autres :

C++/Stroustrup

1. for (i in 1:10)

2. {

3. expression

4. }

K&R (1TBS)

1. for (i in 1:10){

2. expression

3. }

Whitesmith

1. for (i in 1:10)

2. {

3. expression

4. }

GNU

1. for (i in 1:10)

2. {

3. expression

4. }

- Pour des raisons générales de lisibilité et de popularité, le style C++, avec les accolades sur leurs propres lignes et une indentation de quatre (4) espaces est considéré comme standard pour la programmation en R.

- Consulter la documentation de votre éditeur de texte pour savoir s'il est possible de configurer le niveau d'indentation. La plupart des bons éditeurs pour programmeurs le permettent.

- Surtout, éviter de ne pas du tout indenter le code.

Solution du problème

Il existe bien évidemment une multitude de solutions valides. Celle que nous proposons à la figure 5.2 repose sur deux idées principales :

1. Le vecteur créé pour accueillir les résultats contient la valeur de départ en première position afin d'éviter de traiter la première boucle comme une exception ; cette valeur de départ est supprimée du vecteur au moment de retourner les résultats ;

2. Les étiquettes de mouvement du titre sont rapidement converties en valeurs numériques qui permettent de calculer les valeurs successives du titre.

Ensemble, ces deux stratégies permettent d'en arriver à une fonction compacte et efficace.
On remarquera également que la création d'un contenant pour les résultats permet d'éviter le Syndrome de la plaque à biscuits™.

Fig. 5.2 - Fonction valeurs solution du problème du chapitre

1. valeurs <- function(x, start = 100)

2. {

3. ## Création d'un vecteur pour les résultats. La valeur

4. ## de départ est placée au début du vecteur pour faire

5. ## la boucle. Elle sera supprimée à la fin.

6. res <- c(start, numeric(nrow(x)))

8. ## Conversion des étiquettes ("hausse", "nul", "baisse")

9. ## de la première colonne des données en

10. ## valeurs numériques (1, 0, -1).

11. d <- (x[, 1] == "hausse") - (x[, 1] == "baisse")

13. ## Calcul des valeurs successives du titre.

14. for(i in seq(length(res) - 1))

15. res[i + 1] <- res[i] + d[i] * x[i, 2]

17. ## Résultats sans la valeur de départ

18. res[-1]

19. }

5-8 - Exemples

1. ### POINT FIXE

3. ## Comme premier exemple de fonction, on réalise une mise en

4. ## œuvre de l'algorithme du point fixe pour trouver le taux

5. ## d'intérêt tel que a_angle{n} = k pour 'n' et 'k' donnés.

6. ## Cette mise en œuvre est peu générale puisqu'il faudrait

7. ## modifier la fonction chaque fois que l'on change la

8. ## fonction f(x) dont on cherche le point fixe.

9. fp1 <- function(k, n, start = 0.05, TOL = 1E-10)

10. {

11. i <- start

12. repeat

13. {

14. it <- i

15. i <- (1 - (1 + it)^(-n))/k

16. if (abs(i - it)/it < TOL)

17. break

18. }

19. i

20. }

22. fp1(7.2, 10)            # valeur de départ par défaut

23. fp1(7.2, 10, 0.06)      # valeur de départ spécifiée

24. i                       # les variables n'existent pas...

25. start                   # ... dans l'espace de travail

27. ## Généralisation de la fonction 'fp1': la fonction f(x) dont

28. ## on cherche le point fixe (c'est-à-dire la valeur de 'x'

29. ## tel que f(x) = x) est passée en argument. On peut faire

30. ## ça ? Bien sûr, puisqu'une fonction est un objet comme un

31. ## autre en R. On ajoute également à la fonction un argument

32. ## 'echo' qui, lorsque TRUE, fera en sorte d'afficher à

33. ## l'écran les valeurs successives de 'x'.

34. ##

35. ## Ci-dessous, il est implicite que le premier argument, FUN,

36. ## est une fonction.

37. fp2 <- function(FUN, start, echo = FALSE, TOL = 1E-10)

38. {

39. x <- start

40. repeat

41. {

42. xt <- x

44. if (echo)       # inutile de faire 'if (echo == TRUE)'

45. print(xt)

47. x <- FUN(xt)    # appel de la fonction

48. if (abs(x - xt)/xt < TOL)

49. break

50. }

51. x

52. }

54. f <- function(i) (1 - (1+i)^(-10))/7.2 # définition de f(x)

55. fp2(f, 0.05)              # solution

56. fp2(f, 0.05, echo = TRUE) # avec résultats intermédiaires

57. fp2(function(x) 3^(-x), start = 0.5) # avec fonction anonyme

59. ## Amélioration mineure à la fonction 'fp2': puisque la

60. ## valeur de 'echo' ne change pas pendant l'exécution de la

61. ## fonction, on peut éviter de refaire le test à chaque

62. ## itération de la boucle. Une solution élégante consiste à

63. ## utiliser un outil avancé du langage R : les expressions.

64. ##

65. ## L'objet créé par la fonction 'expression' est une

66. ## expression non encore évaluée (comme si on n'avait pas

67. ## appuyé sur Entrée à la fin de la ligne). On peut ensuite

68. ## évaluer l'expression (appuyer sur Entrée) avec 'exec'.

69. fp3 <- function(FUN, start, echo = FALSE, TOL = 1E-10)

70. {

71. x <- start

73. ## Choisir l'expression à exécuter plus loin

74. if (echo)

75. expr <- expression(print(xt <- x))

76. else

77. expr <- expression(xt <- x)

79. repeat

80. {

81. eval(expr) # évaluer l'expression

82. x <- FUN(xt) # appel de la fonction

83. if (abs(x - xt)/xt < TOL)

84. break

85. }

86. x

87. }

89. fp3(f, 0.05, echo = TRUE) # avec résultats intermédiaires

90. fp3(function(x) 3^(-x), start = 0.5) # avec une fonction anonyme

92. ### SUITE DE FIBONACCI

94. ## On a présenté au chapitre 4 deux manières différentes de

95. ## pour calculer les 'n' premières valeurs de la suite de

96. ## Fibonacci. On crée d'abord des fonctions à partir de ce

97. ## code. Avantage d'avoir des fonctions : elles sont valides

98. ## pour tout 'n' > 2.

99. ##

100. ## D'abord la version inefficace parce qu'elle souffre du

101. ## Syndrome de la plaque à biscuits décrit au chapitre 4.

102. fib1 <- function(n)

103. {

104. res <- c(0, 1)

105. for (i in 3:n)

106. res[i] <- res[i - 1] + res[i - 2]

107. res

108. }

109. fib1(10)

110. fib1(20)

112. ## Puis la version qui devrait s'avérer plus efficace parce

113. ## que l'on initialise d'entrée de jeu un contenant de la

114. ## bonne longueur qu'on remplit par la suite.

115. fib2 <- function(n)

116. {

117. res <- numeric(n)       # contenant créé

118. res[2] <- 1             # res[1] vaut déjà 0

119. for (i in 3:n)

120. res[i] <- res[i - 1] + res[i - 2]

121. res

122. }

123. fib2(5)

124. fib2(20)

126. ## A-t-on vraiment gagné en efficacité? Comparons le temps

127. ## requis pour générer une longue suite de Fibonacci avec les

128. ## deux fonctions.

129. system.time(fib1(10000)) # version inefficace

130. system.time(fib2(10000)) # version efficace, ~5x plus rapide

132. ## Variation sur un même thème : une fonction pour calculer non

133. ## pas les 'n' premières valeurs de la suite de Fibonacci,

134. ## mais uniquement la 'n'ième valeur.

135. ##

136. ## Mais il y a un mais : la fonction 'fib3' est truffée

137. ## d'erreurs (de syntaxe, d'algorithmique, de conception). À

138. ## vous de trouver les bogues. (Afin de préserver cet

139. ## exemple, copier le code erroné plus bas ou dans un autre

140. ## fichier avant d'y faire les corrections.)

141. fib3 <- function(nb)

142. {

143. x <- 0

144. x1 _ 0

145. x2 <- 1

146. while (n > 0)

147. {

148. x <- x1 + x2

149. x2 <- x1

150. x1 <- x

151. n <- n - 1

152. }

153. }

154. fib3(1)         # devrait donner 0

155. fib3(2)         # devrait donner 1

156. fib3(5)         # devrait donner 3

157. fib3(10)        # devrait donner 34

158. fib3(20)        # devrait donner 4181

5-9 - Exercices

5.1 (solution) La fonction var calcule l'estimateur sans biais de la variance d'une population à partir de l'échantillon donné en argument. Écrire une fonction variance qui calculera l'estimateur biaisé ou sans biais selon que l'argument biased sera TRUE ou FALSE, respectivement. Le comportement par défaut de variance devrait être le même que celui de var. L'estimateur sans biais de la variance à partir d'un échantillon est

alors que l'estimateur biaisé est

où .

5.2 (solution) Écrire une fonction matrix2 qui, contrairement à la fonction matrix, remplira par défaut la matrice par ligne. La fonction ne doit pas utiliser matrix. Les arguments de la fonction matrix2 seront les mêmes que ceux de matrix, sauf que l'argument byrow sera remplacé par bycol.

5.3 (solution) Écrire une fonction phi servant à calculer la fonction de densité de probabilité d'une loi normale centrée réduite, soit

La fonction devrait prendre en argument un vecteur de valeurs de x.Comparer les résultats avec ceux de la fonction dnorm.

5.4 (solution) Écrire une fonction Phi servant à calculer la fonction de répartition d'une loi normale centrée réduite, soit

Supposer, pour le moment, que . L'évaluation numérique de l'intégrale ci-dessus peut se faire avec l'identité

Utiliser la fonction phi de l'exercice 5.3 et tronquer la somme infinie à une grande valeur, 50 par exemple. La fonction ne doit pas utiliser de boucles, mais peut ne prendre qu'une seule valeur de x à la fois. Comparer les résultats avec ceux de la fonction pnorm.

5.5 (solution) Modifier la fonction Phi de l'exercice 5.4 afin qu'elle admette des valeurs négatives. Lorsque . La solution simple consiste à utiliser une structure de contrôle if ... else, mais les curieux chercheront à s'en passer.

5.6 (solution) Généraliser maintenant la fonction de l'exercice 5.5 pour qu'elle prenne en argument un vecteur de valeurs de . Ne pas utiliser de boucle. Comparer les résultats avec ceux de la fonction pnorm.

5.7 (solution) Sans utiliser l'opérateur %*%, écrire une fonctionprod.mat qui effectuera le produit matriciel de deux matrices seulement si les dimensions de celles-ci le permettent. Cette fonction aura deux arguments (mat1 et mat2) et devra tout d'abord vérifier si le produit matriciel est possible. Si celui-ci est impossible, la fonction retourne un message d'erreur.

1. Utiliser une structure de contrôle if ... else et deux boucles.

2. Utiliser une structure de contrôle if ... else et une seule boucle.

Dans chaque cas, comparer le résultat avec l'opérateur %*%.

5.8 (solution) Vous devez calculer la note finale d'un groupe d'étudiants à partir de deux informations : 1) une matrice contenant la note sur 100 des étudiants à chacune des évaluations, et 2) un vecteur contenant la pondération des évaluations. Un de vos collègues a composé la fonction notes.finales ci-dessous afin de faire le calcul de la note finale pour chacun de ses étudiants. Votre collègue vous mentionne toutefois que sa fonction est plutôt lente et inefficace pour de grands groupes d'étudiants. Modifiez la fonction afin d'en réduire le nombre d'opérations et faire en sorte qu'elle n'utilise aucune boucle.

1. notes.finales <- function(notes, p)

2. {

3. netud <- nrow(notes)

4. neval <- ncol(notes)

5. final <- (1:netud) * 0

6. for(i in 1:netud)

7. {

8. for(j in 1:neval)

9. {

10. final[i] <- final[i] + notes[i, j] * p[j]

11. }

12. }

13. final

14. }

5.9 (solution) Trouver les erreurs qui empêchent la définition de la fonction ci-dessous.

1. AnnuiteFinPeriode <- function(n, i)

2. {{

3. v <- 1/1 + i)

4. ValPresChaquePmt <- v^(1:n)

5. sum(ValPresChaquepmt)

6. }

5.10 (solution) La fonction ci-dessous calcule la valeur des paramètres d'une loi normale, gamma ou Pareto à partir de la moyenne et de la variance, qui sont connues par l'utilisateur.

1. param <- function(moyenne, variance, loi)

2. {

3. loi <- tolower(loi) {

4. if (loi == "normale")

5. param1 <- moyenne

6. param2 <- sqrt(variance)

7. return(list(mean = param1, sd = param2))

8. }

9. if (loi == "gamma") {

10. param2 <- moyenne/variance

11. param1 <- moyenne * param2

12. return(list(shape = param1, scale = param2))

13. }

14. if (loi == "pareto") {

15. cte <- variance/moyenne^2

16. param1 <- 2 * cte/(cte-1)

17. param2 <- moyenne * (param1 - 1)

18. return(list(alpha = param1, lambda = param2))

19. }

20. stop("La loi doit être une de \"normale\", \"gamma\" ou \"pareto\"")

21. }

L'utilisation de la fonction pour diverses lois donne les résultats suivants :

1. > param(2, 4, "normale")

2. $mean

3. [1] 2

4. $sd

5. [1] 2

6. > param(50, 7500, "gamma")

7. Erreur dans param(50, 7500, "gamma") : Objet "param1" introuvable

8. > param(50, 7500, "pareto")

9. Erreur dans param(50, 7500, "pareto") : Objet "param1" introuvable

1. Expliquer pour quelle raison la fonction se comporte ainsi.

2. Appliquer les correctifs nécessaires à la fonction pour que celle-ci puisse calculer les bonnes valeurs. (Les erreurs ne se trouvent pas dans les mathématiques de la fonction.) Astuce : tirer profit du moteur d'indentation de votre éditeur de texte pour programmeur.

6 - Concepts avancés

Objectifs du chapitre

- Passer des valeurs à une fonction via l'argument ....

- Effectuer des sommaires sur des tableaux à l'aide de la fonction apply.

- Réduire des listes avec les fonctions lapply, sapply et mapply ; comparer l'effet de ces fonctions.

- Concevoir comment la classe d'un objet peut modifier le traitement qu'en feront les fonctions génériques.

Ce chapitre traite de divers concepts et fonctions un peu plus avancés du langage R, dont les fonctions de la famille apply auxquelles nous avons fait référence à quelques reprises dans les chapitres précédents. Ce sont des fonctions d'une grande importance en R.

Énoncé du problème

Soit la fonction

On essaie d'en faire une mise en œuvre en R :

1. > g <- function(x)

2. + {

3. +     if (x <= 5)

4. +         2 * x

5. +     else

6. +         0

7. + }

On évalue . On devrait obtenir (cas ).

> g(1)
[1] 2

On évalue . On devrait obtenir 0 (cas ).

> g(10)
[1] 0

On désire maintenant évaluer la fonction à tous les entiers de 1 à 10, inclusivement. On devrait obtenir (2, 4, 6, 8, 10, 0, 0, 0, 0, 0).

1. > g(1:10)

2. [1] 2 4 6 8 10 12 14 16 18 20

3. Warning message:

4. In if (x <= 5) 2 * x else 0 :

5. la condition a une longueur > 1 et seul le premier élément est utilisé

Que se passe-t-il ? Les résultats sont erronés et le message d'avertissement nous indique que seulement la première valeur a été utilisée pour évaluer la condition. La fonction g évalue donc correctement son argument seulement lorsqu'il s'agit d'un vecteur de longueur 1. En d'autres termes, la fonction n'est pas vectorielle.

6-1 - Argument '...'

La mention ... apparaît dans la définition de plusieurs fonctions en R. Il ne faut pas voir là de la paresse de la part des rédacteurs des rubriques d'aide, mais bel et bien un argument formel dont ... est le nom.

- Cet argument signifie qu'une fonction peut accepter un ou plusieurs arguments autres que ceux faisant partie de sa définition.

- Le contenu de l'argument ... n'est ni pris en compte, ni modifié par la fonction. Il est généralement simplement passé tel quel à une autre fonction qui, elle, saura traiter les arguments qui lui sont ainsi passés.

- Pour des exemples, voir les définitions des fonctions apply, lapply et sapply, ci-dessous.

6-2 - Fonction apply

La fonction apply sert à appliquer une fonction quelconque sur une partie d'une matrice ou, plus généralement, d'un tableau. La syntaxe de la fonction est la suivante :

apply(X, MARGIN, FUN, ...)

où

- X est une matrice ou un tableau ;

- MARGIN est un vecteur d'entiers contenant la ou les dimensions de la matrice ou du tableau sur lesquelles la fonction doit s'appliquer ;

- FUN est la fonction à appliquer ;

- ... est un ensemble d'arguments supplémentaires, séparés par des virgules, à passer à la fonction FUN.

Lorsque X est une matrice, apply sert principalement à calculer des sommaires par ligne (dimension 1) ou par colonne (dimension 2) autres que la somme ou la moyenne (puisque les fonctions rowSums, colSums, rowMeans et colMeans existent pour ce faire).

- Utiliser la fonction apply plutôt que des boucles puisque celle-ci est plus efficace.

- Considérer les exemples suivants :

1. > (x <- matrix(sample(1:100, 20, rep = TRUE), 5, 4))

2. [,1] [,2] [,3] [,4]

3. [1,]   27   90   21   50

4. [2,]   38   95   18   72

5. [3,]   58   67   69  100

6. [4,]   91   63   39   39

7. [5,]   21    7   77   78

8. > apply(x, 1, var) # variance par ligne

9. [1]  978.000 1181.583  335.000  612.000 1376.917

10. > apply(x, 2, min) # minimum par colonne

11. [1]  21   7  18  39

12. > apply(x, 1, mean, trim = 0.2) # moy. tronquée par ligne

13. [1] 47.00 55.75 73.50 58.00 45.75

Puisqu'il n'existe pas de fonctions internes pour effectuer des sommaires sur des tableaux, il faut toujours utiliser la fonction apply.

- Si X est un tableau de plus de deux dimensions, alors l'argument passé à FUN peut être une matrice ou un tableau.

- Lorsque X est un tableau à trois dimensions et que MARGIN est de longueur 1, cela équivaut à appliquer la fonction FUN sur des « tranches » (des matrices) de X. Si MARGIN est de longueur 2, on applique FUN sur des « carottes » (des vecteurs) tirées de X.

- Truc mnémotechnique : la ou les dimensions absentes de MARGIN sont celles qui disparaissent après le passage de apply.

- Considérer les exemples suivants :

1. > (x <- array(sample(1:10, 80, rep = TRUE), c(3, 3, 4)))

2. , , 1

3. [,1] [,2] [,3]

4. [1,]   10    2    1

5. [2,]    3    3    4

6. [3,]    7    4    9

7. , , 2

8. [,1] [,2] [,3]

9. [1,]    4    5    7

10. [2,]    5    2    8

11. [3,]    6    9    2

12. , , 3

13. [,1] [,2] [,3]

14. [1,]    8    7    6

15. [2,]    5    8    8

16. [3,]    9    6    1

17. , , 4

18. [,1] [,2] [,3]

19. [1,]    5    5    3

20. [2,]    8    9    1

21. [3,]    7    5    1

22. > apply(x, 3, det) # déterminants matrices 3 x 3

23. [1] 103 149 -103 -54

24. > apply(x, 1, sum) # sommes tranches horizontales

25. [1] 63 64 66

26. > apply(x, c(1, 2), sum) # sommes carottes horizontales

27. [,1] [,2] [,3]

28. [1,]   27   19   17

29. [2,]   21   22   21

30. [3,]   29   24   13

31. > apply(x, c(1, 3), sum) # sommes carottes transversales

32. [,1] [,2] [,3] [,4]

33. [1,]   13   16   21   13

34. [2,]   10   15   21   18

35. [3,]   20   17   16   13

36. > apply(x, c(2, 3), sum) # sommes carottes verticales

37. [,1] [,2] [,3] [,4]

38. [1,]   20   15   22   20

39. [2,]    9   16   21   19

40. [3,]   14   17   15    5

6-3 - Fonctions lapply et sapply

Les fonctions lapply et sapply sont similaires à la fonction apply en ce qu'elles permettent d'appliquer une fonction aux éléments d'une structure — le vecteur ou la liste en l'occurrence. Leur syntaxe est similaire :

lapply(X, FUN, …)
sapply(X, FUN, ...)

- La fonction lapply applique une fonction FUN à tous les éléments d'un vecteur ou d'une liste X et retourne le résultat sous forme de liste. Le résultat est donc :

list(FUN(X[[1]], ...),
FUN(X[[2]], ...),
FUN(X[[3]], ...),
...)

- Les éléments de X sont passés comme à la fonction FUN sans être nommés. Les règles habituelles d'évaluation d'un appel de fonction s'appliquent. Par conséquent, les éléments de X seront considérés comme les premiers arguments de FUN à moins que des arguments nommés dans ... aient préséance.

- Par exemple, on crée une liste formée de quatre vecteurs aléatoires de taille 5, 6, 7 et 8 :

1. > (x <- lapply(5:8, sample, x = 1:10))

2. [[1]]

3. [1] 7 4 3 10 9

4. [[2]]

5. [1] 3 2 4 7 8 5

6. [[3]]

7. [1] 8 4 9 2 1 10 6

8. [[4]]

9. [1] 5 6 8 4 3 1 7 2

Le premier argument de la fonction sample est x. Dans l'expression ci-dessus, cet argument est passé à la fonction via l'argument ... de lapply. Par conséquent, les valeurs successives de 5:8 servent comme deuxième argument à la fonction sample, soit la taille de l'échantillon.

- On peut ensuite calculer la moyenne de chacun des vecteurs obtenus ci-dessus, toujours sans faire de boucle :

1. > lapply(x, mean)

2. [[1]]

3. [1] 6.6

4. [[2]]

5. [1] 4.833333

6. [[3]]

7. [1] 5.714286

8. [[4]]

9. [1] 4.5

La fonction sapply est similaire à lapply, sauf que le résultat est retourné sous forme de vecteur, si possible. Le résultat est donc simplifié par rapport à celui de lapply, d'où le nom de la fonction.

- Dans l'exemple ci-dessus, il est souvent plus utile d'obtenir les résultats sous la forme d'un vecteur :

> sapply(x, mean)
[1] 6.600000 4.833333 5.714286 4.500000

- Si le résultat de chaque application de la fonction est un vecteur et que les vecteurs sont tous de la même longueur, alors sapply retourne une matrice, remplie comme toujours par colonne :

1. > (x <- lapply(rep(5, 3), sample, x = 1:10))

2. [[1]]

3. [1] 6 2 9 5 4

4. [[2]]

5. [1] 1 10 6 7 4

6. [[3]]

7. [1] 6 5 8 4 9

8. > sapply(x, sort)

9. [,1] [,2] [,3]

10. [1,]    2    1    4

11. [2,]    4    4    5

12. [3,]    5    6    6

13. [4,]    6    7    8

14. [5,]    9   10    9

Dans un grand nombre de cas, il est possible de remplacer les boucles for par l'utilisation de lapply ou sapply. On ne saurait donc trop insister sur l'importance de ces fonctions.

La fonction sapply permet de vectoriser une fonction R qui n'est pas vectorielle. On peut procéder de deux façons : en utilisant la fonction non vectorielle dans une application de sapply() ou en vectorisant la fonction à y ajoutant un appel à sapply.

6-4 - Fonction mapply

La fonction mapply est une version multidimensionnelle de sapply. Sa syntaxe est, essentiellement,

mapply(FUN, ...)

- Le résultat de mapply est l'application de la fonction FUN aux premiers éléments de tous les arguments contenus dans ..., puis à tous les seconds éléments, et ainsi de suite.

- Ainsi, si v et w sont des vecteurs, mapply(FUN, v, w) retourne sous forme de liste, de vecteur ou de matrice, selon le cas, FUN(v[1], w[1]), FUN(v[2], w[2]), etc.

- Par exemple :

1. > mapply(rep, 1:4, 4:1)

2. [[1]]

3. [1] 1 1 1 1

4. [[2]]

5. [1] 2 2 2

6. [[3]]

7. [1] 3 3

8. [[4]]

9. [1] 4

- Les éléments de '...' sont recyclés au besoin.

1. > mapply(seq, 1:6, 6:8)

2. [[1]]

3. [1] 1 2 3 4 5 6

4. [[2]]

5. [1] 2 3 4 5 6 7

6. [[3]]

7. [1] 3 4 5 6 7 8

8. [[4]]

9. [1] 4 5 6

10. [[5]]

11. [1] 5 6 7

12. [[6]]

13. [1] 6 7 8

6-5 - Fonction replicate

La fonction replicate est une fonction enveloppante de sapply simplifiant la syntaxe pour l'exécution répétée d'une expression.

- Son usage est particulièrement indiqué pour les simulations. Ainsi, on peut construire une fonction fun qui fait tous les calculs d'une simulation, puis obtenir les résultats pour, disons, 10 000 simulations avec

> replicate(10000, fun(...))

- L'annexe C présente en détail différentes stratégies — dont l'utilisation de la fonction replicate — pour la réalisation d'études de simulation en R.

6-6 - Classes et fonctions génériques

Dans le langage R, tous les objets ont une classe. La classe est parfois implicite ou dérivée du mode de l'objet (consulter la rubrique d'aide de class pour de plus amples détails).

- Certaines fonctions, dites fonctions génériques, se comportent différemment selon la classe de l'objet donné en argument. Les fonctions génériques les plus fréquemment employées sont print, plot et summary.

- Une fonction générique possède une méthode correspondant à chaque classe qu'elle reconnaît et, généralement, une méthode default pour les autres objets. La liste des méthodes existant pour une fonction générique s'obtient avec la fonction methods :

1. > methods(plot)

2. [1] plot.acf*           plot.data.frame*

3. [3] plot.decomposed.ts* plot.default

4. [5] plot.dendrogram*    plot.density*

5. [7] plot.ecdf           plot.factor*

6. [9] plot.formula*       plot.function

7. [11] plot.hclust*        plot.histogram*

8. [13] plot.HoltWinters*   plot.isoreg*

9. [15] plot.lm*            plot.medpolish*

10. [17] plot.mlm*           plot.ppr*

11. [19] plot.prcomp*        plot.princomp*

12. [21] plot.profile.nls*   plot.raster*

13. [23] plot.spec*          plot.stepfun

14. [25] plot.stl*           plot.table*

15. [27] plot.ts             plot.tskernel*

16. [29] plot.TukeyHSD*

17. see ' ?methods' for accessing help and source code

- À chaque méthode methode d'une fonction générique fun correspond une fonction fun.methode. C'est donc la rubrique d'aide de cette dernière fonction qu'il faut consulter au besoin, et non celle de la fonction générique, qui contient en général peu d'informations.

- Il est intéressant de savoir que lorsque l'on tape le nom d'un objet à la ligne de commande pour voir son contenu, c'est la fonction générique print qui est appelée. On peut donc complètement modifier la représentation à l'écran du contenu d'un objet en créant une nouvelle classe et une nouvelle méthode pour la fonction print.

Solution du problème

Débutons par expliquer ce qui cloche avec notre fonction g de départ :

1. > g

2. function(x)

3. {

4. if (x <= 5)

5. 2 * x

6. else

7. 0

8. }

Il s'avère que la structure de contrôle if n'est pas une fonction vectorielle. Elle n'accepte qu'une seule valeur. De la rubrique d'aide :

cond: A length-one logical vector that is not 'NA'.
Conditions of length greater than one are accepted
with a warning, but only the first element is used.

Par conséquent, lorsque l'argument de la fonction g est 1:10, l'exécution de la fonction est la suivante :

1. L'argument x vaut 1:10 ;

2. La condition if (x <= 5) n'est évaluée que pour la première valeur de x ;

3. Puisque x[1] <= 5, c'est l'expression vectorielle 2 * x qui est évaluée ;

4. Le résultat de la fonction est 2 * x.

Comme indiqué précédemment, la fonction sapply peut rendre vectorielles une fonction ou une expression qui ne le sont pas. Ici, la solution la plus simple à court terme serait :

> sapply(1:10, g)
[1] 2 4 6 8 10 0 0 0 0 0

Une solution plus complète consisterait à transformer la fonction g pour la rendre vectorielle avec la fonction sapply :

1. > g2 <- function(x)

2. + {

3. +     h <- function(x)

4. +     {

5. +         if (x <= 5)

6. +             2 * x

7. +         else

8. +             0

9. +     }

10. +     sapply(x, h)

11. + }

12. > g2(1:10)

13. [1] 2 4 6 8 10 0 0 0 0 0

Attention, toutefois, cette solution n'est pas une panacée. Par exemple, dans ce cas bien précis où la fonction est une fonction en branches, une autre solution fait appel à la fonction ifelse :

1. > g3 <- function(x)

2. +     ifelse (x <= 5, 2 * x, 0)

3. > g3(1:10)

4. [1] 2 4 6 8 10 0 0 0 0 0

Pour le problème sous étude, on peut faire encore beaucoup mieux en tirant directement profit de l'approche vectorielle de R. Par une judicieuse utilisation des valeurs booléennes, il est possible d'éliminer complètement la condition if :

1. > g4 <- function(x)

2. +     2 * x * (x <= 5)

3. > g4(1:10)

4. [1] 2 4 6 8 10 0 0 0 0 0

On préférera, dans l'ordre, les fonctions g4, g3 (la fonction ifelse est lente, mais néanmoins plus rapide que sapply), g2 et l'approche combinant sapply et la fonction g de départ.

6-7 - Exemples

1. ###

2. ### FONCTION 'apply'

3. ###

5. ## Création d'une matrice et d'un tableau à trois dimensions

6. ## pour les exemples.

7. m <- matrix(sample(1:100, 20), nrow = 4, ncol = 5)

8. a <- array(sample(1:100, 60), dim = 3:5)

10. ## Les fonctions 'rowSums', 'colSums', 'rowMeans' et

11. ## 'colMeans' sont des raccourcis pour des utilisations

12. ## fréquentes de 'apply'.

13. rowSums(m)                      # somme par ligne

14. apply(m, 1, sum)                # idem, mais moins lisible

15. colMeans(m)                     # somme par colonne

16. apply(m, 2, mean)               # idem, mais moins lisible

18. ## Puisqu'il n'existe pas de fonctions comme 'rowMax' ou

19. ## 'colProds', il faut utiliser 'apply'.

20. apply(m, 1, max)                # maximum par ligne

21. apply(m, 2, prod)               # produit par colonne

23. ## L'argument '... ' de 'apply' permet de passer des arguments

24. ## à la fonction FUN.

25. m[sample(1:20, 5)] <- NA        # ajout de données manquantes

26. apply(m, 1, var, na.rm = TRUE)  # variance par ligne sans NA

28. ## Lorsque 'apply' est utilisée sur un tableau, son résultat

29. ## est de dimension dim(X)[MARGIN], d'où le truc

30. ## mnémotechnique donné dans le texte du chapitre.

31. apply(a, c(2, 3), sum) # le résultat est une matrice

32. apply(a, 1, prod) # le résultat est un vecteur

34. ## L'utilisation de 'apply' avec les tableaux peut rapidement

35. ## devenir confondante si l'on ne visualise pas les calculs

36. ## qui sont réalisés. On reprend ici les exemples du chapitre

37. ## en montrant comment l'on calculerait le premier élément de

38. ## chaque utilisation de 'apply'. Au besoin, retourner à

39. ## l'indiçage des tableaux au chapitre 2.

40. (x <- array(sample(1:10, 80, rep = TRUE), c(3, 3, 4)))

41. apply(x, 3, det)            # déterminants des quatre matrices 3 x 3

42. det(x[, , 1])               # équivalent pour le premier déterminant

44. apply(x, 1, sum)            # sommes des trois tranches horizontales

45. sum(x[1, , ])               # équivalent pour la première somme

47. apply(x, c(1, 2), sum)      # sommes des neuf carottes horizontales

48. sum(x[1, 1, ])              # équivalent pour la première somme

50. apply(x, c(1, 3), sum)      # sommes des 12 carottes transversales

51. sum(x[1, , 1])              # équivalent pour la première somme

53. apply(x, c(2, 3), sum)      # sommes des 12 carottes verticales

54. sum(x[, 1, 1])              # équivalent pour la première somme

56. ###

57. ### FONCTIONS 'lapply' ET 'sapply'

58. ###

60. ## La fonction 'lapply' applique une fonction à tous les

61. ## éléments d'un vecteur ou d'une liste et retourne une liste,

62. ## peu importe les dimensions des résultats. La fonction

63. ## 'sapply' retourne un vecteur ou une matrice, si possible.

64. ##

65. ## Somme « interne » des éléments d'une liste.

66. (x <- list(1:10, c(-2, 5, 6), matrix(3, 4, 5)))

67. sum(x) # erreur

68. lapply(x, sum)              # sommes internes (liste)

69. sapply(x, sum)              # sommes internes (vecteur)

71. ## Création de la suite 1, 1, 2, 1, 2, 3, 1, 2, 3, 4, ... , 1,

72. ## 2, ... , 9, 10.

73. lapply(1:10, seq)           # le résultat est une liste

74. unlist(lapply(1:10, seq))   # le résultat est un vecteur

76. ## Soit une fonction calculant la moyenne pondérée d'un

77. ## vecteur. Cette fonction prend en argument une liste de deux

78. ## éléments : 'donnees' et 'poids'.

79. fun <- function(x)

80. sum(x$donnees * x$poids)/sum(x$poids)

82. ## On peut maintenant calculer la moyenne pondérée de

83. ## plusieurs ensembles de données réunis dans une liste

84. ## itérée.

85. (x <- list(list(donnees = 1:7,

86. poids = (5:11)/56),

87. list(donnees = sample(1:100, 12),

88. poids = 1:12),

89. list(donnees = c(1, 4, 0, 2, 2),

90. poids = c(12, 3, 17, 6, 2))))

91. sapply(x, fun)              # aucune boucle explicite !

93. ###

94. ### EXEMPLES ADDITIONNELS SUR L'UTILISATION DE L'ARGUMENT

95. ### '... ' AVEC 'lapply' ET 'sapply'

96. ###

98. ## Aux fins des exemples ci-dessous, on crée d'abord une liste

99. ## formée de nombres aléatoires. Cette expression fait usage

100. ## de l'argument '... ' de 'lapply'. Pouvez-vous la décoder?

101. ## Nous y reviendrons plus loin, ce qui compte pour le moment

102. ## c'est simplement de l'exécuter.

103. x <- lapply(c(8, 12, 10, 9), sample, x = 1:10, replace = TRUE)

105. ## Soit maintenant une fonction qui calcule la moyenne

106. ## arithmétique des données d'un vecteur 'x' supérieures à une

107. ## valeur 'y'. On remarquera que cette fonction n'est pas

108. ## vectorielle pour 'y', c'est-à-dire qu'elle n'est valide que

109. ## lorsque 'y' est un vecteur de longueur 1.

110. fun <- function(x, y) mean(x[x > y])

112. ## Pour effectuer ce calcul sur chaque élément de la liste

113. ## 'x', nous pouvons utiliser 'sapply' plutôt que 'lapply',

114. ## car chaque résultat est de longueur 1. Cependant, il faut

115. ## passer la valeur de 'y' à la fonction 'fun'. C'est là

116. ## qu'entre en jeu l'argument '... ' de 'sapply'.

117. sapply(x, fun, 7) # moyennes des données > 7

119. ## Les fonctions 'lapply' et 'sapply' passent tour à tour les

120. ## éléments de leur premier argument comme *premier* argument

121. ## à la fonction, sans le nommer explicitement. L'expression

122. ## ci-dessus est donc équivalente à

123. ##

124. ##   c(fun(x[[1]], 7), ... , fun(x[[4]], 7)

125. ##

126. ## Que se passe-t-il si l'on souhaite passer les valeurs à un

127. ## argument de la fonction autre que le premier ? Par exemple,

128. ## supposons que l'ordre des arguments de la fonction 'fun'

129. ## ci-dessus est inversé.

130. fun <- function(y, x) mean(x[x > y])

132. ## Les règles d'appariement des arguments des fonctions en R

133. ## font en sorte que lorsque les arguments sont nommés dans

134. ## l'appel de fonction, leur ordre n'a pas d'importance. Par

135. ## conséquent, un appel de la forme

136. ##

137. ##   fun(x, y = 7)

138. ##

139. ## est tout à fait équivalent à fun(7, x). Pour effectuer les calculs

140. ##

141. ##   c(fun(x[[1]], y = 7), ... , fun(x[[4]], y = 7)

142. ##

143. ## avec la liste définie plus haut, il s'agit de nommer

144. ## l'argument 'y' dans '... ' de 'sapply'.

145. sapply(x, y = 7)

147. ## Décodons maintenant l'expression

148. ##

149. ## lapply(c(8, 12, 10, 9), sample, x = 1:10, replace = TRUE)

150. ##

151. ## qui a servi à créer la liste. La définition de la fonction

152. ## 'sample' est la suivante :

153. ##

154. ##   sample(x, size, replace = FALSE, prob = NULL)

155. ##

156. ## L'appel à 'lapply' est équivalent à

157. ##

158. ##   list(sample(8, x = 1:10, replace = TRUE),

159. ##      ... ,

160. ##      sample(9, x = 1:10, replace = TRUE))

161. ##

162. ## Toujours selon les règles d'appariement des arguments, on

163. ## voit que les valeurs 8, 12, 10, 9 seront attribuées à

164. ## l'argument 'size', soit la taille de l'échantillon.

165. ## L'expression crée donc une liste comprenant quatre

166. ## échantillons aléatoires de tailles différentes des nombres

167. ## de 1 à 10 pigés avec remise.

168. ##

169. ## Une expression équivalente, quoique moins élégante, aurait

170. ## recours à une fonction anonyme pour replacer les arguments

171. ## de 'sample' dans l'ordre voulu.

172. lapply(c(8, 12, 10, 9),

173. function(x) sample(1:10, x, replace = TRUE))

175. ## La fonction 'sapply' est aussi très utile pour vectoriser

176. ## une fonction qui n'est pas vectorielle. Supposons que l'on

177. ## veut généraliser la fonction 'fun' pour qu'elle accepte un

178. ## vecteur de seuils 'y'.

179. fun <- function(x, y)

180. sapply(y, function(y) mean(x[x > y]))

182. ## Utilisation sur la liste 'x' avec trois seuils.

183. sapply(x, fun, y = c(3, 5, 7))

185. ###

186. ### FONCTION 'mapply'

187. ###

189. ## Création de quatre échantillons aléatoires de taille 12.

190. x <- lapply(rep(12, 4), sample, x = 1:100)

192. ## Moyennes tronquées à 0, 10, 20 et 30%, respectivement, de

193. ## ces quatre échantillons aléatoires.

194. mapply(mean, x, 0:3/10)

196. ###

197. ### FONCTION 'replicate'

198. ###

200. ## La fonction 'replicate' va répéter un certain nombre de

201. ## fois une expression quelconque. Le principal avantage de

202. ## 'replicate' sur 'sapply' est qu'on n'a pas à se soucier des

203. ## arguments à passer à une fonction.

204. ##

205. ## Par exemple, on veut simuler dix échantillons aléatoires

206. ## indépendants de longueur 12. On peut utiliser 'sapply',

207. ## mais la syntaxe n'est ni élégante, ni facile à lire

208. ## (l'argument 'i' ne sert à rien).

209. sapply(rep(1, 10), function(i) sample(1:100, 12))

211. ## En utilisant 'replicate', on sait tout de suite de quoi il

212. ## s'agit. À noter que les échantillons se trouvent dans les

213. ## colonnes de la matrice résultante.

214. replicate(10, sample(1:100, 12))

216. ## Vérification que la moyenne arithmétique (bar{X}) est un

217. ## estimateur sans biais de la moyenne de la loi normale. On

218. ## doit calculer la moyenne de plusieurs échantillons

219. ## aléatoires, puis la moyenne de toutes ces moyennes.

220. ##

221. ## On définit d'abord une fonction pour faire une simulation.

222. ## Remarquer que dans la fonction ci-dessous, 'mean' est tour

223. ## à tour le nom d'un argument (qui pourrait aussi bien être

224. ## «toto») et la fonction pour calculer une moyenne.

225. fun <- function(n, mean, sd)

226. mean(rnorm(n, mean = mean, sd = sd))

228. ## Avec 'replicate', on fait un grand nombre de simulations.

229. x <- replicate(10000, fun(100, 0, 1)) # 10000 simulations

230. hist(x)                     # distribution de bar{X}

231. mean(x)                     # moyenne de bar{X}

233. ###

234. ### CLASSES ET FONCTIONS GÉNÉRIQUES

235. ###

237. ## Pour illustrer les classes et fonctions génériques, on

238. ## reprend la fonction de point fixe 'fp3' des exemples du

239. ## chapitre 5 en y faisant deux modifications :

240. ##

241. ##   1. Ajout d'un compteur pour le nombre d'itérations ;

242. ##   2. La fonction retourne une liste de classe 'fp'

243. ##      contenant diverses informations relatives à la

244. ##      procédure de point fixe.

245. ##

246. ## Ainsi, la fonction 'fp4' retourne un objet qui peut ensuite

247. ## être manipulé par des méthodes de fonctions génériques.

248. ## C'est l'approche de programmation objet favorisée dans le

249. ## langage R.

250. fp4 <- function(FUN, start, echo = FALSE, TOL = 1E-10)

251. {

252. x <- start              # valeur de départ

253. i <- 0                  # compteur des itérations

255. if (echo)

256. expr <- expression(print(xt <- x))

257. else

258. expr <- expression(xt <- x)

260. repeat

261. {

262. eval(expr)

264. x <- FUN(xt)        # nouvelle valeur

265. i <- i + 1          # incrémenter le compteur

267. if (abs(x - xt)/xt < TOL)

268. break

269. }

271. structure(list(fixed.point = x, # point fixe

272. nb.iter = i,     # nombre d'itérations

273. fun = FUN,       # fonction f(x)

274. x0 = start,      # valeur de départ

275. TOL = TOL),      # précision relative

276. class = "fp")

277. }

279. ## On crée maintenant des méthodes pour la classe 'fp' pour

280. ## les fonctions génériques les plus courantes, soit 'print',

281. ## 'summary' et 'plot'.

282. ##

283. ## La méthode de 'print' sera utilisée pour afficher seulement

284. ## la valeur du point fixe. C'est en quelque sorte

285. ## l'utilisation la plus simple de la fonction 'fp4'.

286. ##

287. ## La méthode de 'summary' fournira un peu plus d'informations

288. ## sur la procédure de point fixe.

289. ##

290. ## Enfin, la méthode de 'plot' fera un graphique de la

291. ## fonction f(x) et son intersection avec la droite y = x.

292. print.fp <- function(x)

293. print(x$fixed.point)

294. summary.fp <- function(x)

295. {

296. if (class(x) != "fp")

297. stop("object is not of class 'fp'")

298. cat("Function:\n ")

299. print(x$fun)

300. cat("\n")

301. cat("Fixed point:\n ", x$fixed.point, fill = TRUE)

302. cat("\n")

303. cat("Number of iterations:\n ", x$nb.iter, fill = TRUE)

304. cat("\n")

305. cat("Precision:\n ", x$TOL, fill = TRUE)

306. }

308. plot.fp <- function(x, ...)

309. {

310. ## Valeur du point fixe

311. fp <- x$fixed.point

313. ## Il faut déterminer un intervalle pour lequel tracer la

314. ## fonction. Celui-ci est déterminé de façon arbitraire

315. ## comme un multiple de la distance entre la valeur de

316. ## départ et la valeur du point fixe.

317. r <- abs(x$x0 - fp)

319. ## Fonction à tracer

320. FUN <- x$fun

322. ## Fonction y = x. 'FUN2' est nécessaire parce que 'curve'

323. ## n'admet pas de fonctions anonymes en argument.

324. FUN2 <- function(x) x

326. ## Graphique de la fonction 'FUN'

327. curve(FUN, from = fp - 3 * r, to = fp + 3 * r,

328. xlab = "x", ylab = "f(x)", lwd = 2)

330. ## Ajout de la droite 'FUN2' au graphique

331. curve(FUN2, add = TRUE, lwd = 1)

333. ## Ajout d'un point sur le point fixe

334. points(fp, FUN(fp), ...)

335. }

337. ## Exemples d'utilisation

338. x <- fp4(function(x) 3^(-x), start = 0.5)

339. x                           # affichage de 'print.fp'

340. summary(x)                  # plus d'information

341. plot(x)                     # graphique de base

342. plot(x, pch = 21,           # graphique plus élaboré...

343. bg = "orange",         # ... consulter la rubrique

344. cex = 2, lwd = 2)      # ... d'aide de 'par'

346. ###

347. ### OPÉRATEURS EN TANT QUE FONCTIONS

348. ###

350. ## Les opérateurs représentés par des caractères spéciaux sont

351. ## des fonctions comme les autres. On peut donc les appeler

352. ## comme toute autre fonction. (En fait, l'interprète R fait

353. ## cette traduction à l'interne. )

354. x <- sample(1:100, 12)      # un vecteur

355. x + 2                       # appel usuel

356. "+"(x, 2)                   # équivalent

357. x[c(3, 5)]                  # extraction usuelle

358. "["(x, c(3, 5))             # équivalent

359. x[1] <- 0; x                # assignation usuelle

360. "[<-"(x, 2, 0)              # équivalent (à x[2] <- 0)

362. ## D'une part, cela explique pourquoi il faut placer les

363. ## opérateurs entre guillemets (" ") lorsqu'on les utilise

364. ## dans les fonctions comme 'outer', 'lapply', etc.

365. outer(x, x, +)              # erreur de syntaxe

366. outer(x, x, "+")            # correct

368. ## D'autre part, cela permet d'utiliser les opérateurs

369. ## d'extraction "[" et "[[" dans de telles fonctions. Par

370. ## exemple, voici comment extraire le deuxième élément de

371. ## chaque élément d'une liste.

372. (x <- list(1:4, 8:2, 6:12, -2:2)) # liste quelconque

373. x[[1]][2]                   # 2e élément du 1er élément

374. x[[2]][2]                   # 2e élément du 2e élément

375. x[[3]][2]                   # 2e élément du 3e élément

376. x[[4]][2]                   # 2e élément du 4e élément

377. lapply(x, "[", 2)           # même chose en une ligne

378. sapply(x, "[", 2)           # résultat sous forme de vecteur

380. ###

381. ### COMMENT JOUER DES TOURS AVEC R

382. ###

384. ## Redéfinir un opérateur dans l'espace de travail de

385. ## quelqu'un...

386. "+" <- function(x, y) x * y # redéfinition de "+"

387. 5 + 2 # ouch!

388. ls()                        # traîtrise dévoilée...

389. rm("+")                     # ... puis éliminée

390. 5 + 2                       # c'est mieux

392. ## Faire croire qu'une fonction fait autre chose que ce

393. ## qu'elle fait en réalité. Si l'attribut "source" d'une

394. ## fonction existe, c'est son contenu qui est affiché lorsque

395. ## l'on examine une fonction.

396. f <- function(x, y) x + y   # vraie fonction

397. attr(f, "source") <- "function(x, y) x * y" # ce qui est affiché

398. f                           # une fonction pour faire le produit?

399. f(2, 3)                     # non!

400. str(f)                      # structure de l'objet

401. attr(f, "source") <- NULL   # attribut "source" effacé

402. f                           # c'est mieux

404. ## Redéfinir la méthode de 'print' pour une classe d'objet...

405. ## Ici, l'affichage d'un objet de classe "lm" cause la

406. ## fermeture de R !

407. print.lm <- function(x) q("ask")

408. x <- rnorm(10)              # échantillon aléatoire

409. y <- x + 2 + rnorm(10)      # modèle de régression linéaire

410. lm(y ~ x)                   # répondre "c"!

6-8 - Exercices

6.1 (solution) À l'exercice 4.2, on a calculé la moyenne pondérée d'un vecteur d'observations

où Si l'on a plutôt une matrice d'observations on peut définir les moyennes pondérées

et

De même, on peut définir des moyennes pondérées calculées à partir d'un tableau de données de dimensions dont la notation suit la même logique que ci-dessus. Écrire des expressions R pour calculer, sans boucle, les moyennes pondérées suivantes.

1. en supposant une matrice de données .

2. en supposant une matrice de données .

3. en supposant une matrice de données .

4. en supposant un tableau de données .

5. en supposant un tableau de données .

6. en supposant un tableau de données .

7. en supposant un tableau de données .

6.2 (solution) Générer les suites de nombres suivantes à l'aide d'expressions R. (Évidemment, il faut trouver un moyen de générer les suites sans simplement concaténer les différentes sous-suites.)

1. 0, 0, 1, 0, 1, 2,…, 0, 1, 2, 3,…, 10.

2. 10, 9, 8,…, 2, 1, 10, 9, 8,… 3, 2,…, 10, 9, 10.

3. 10, 9, 8,…, 2, 1, 9, 8,…, 2, 1,…, 2, 1, 1.

6.3 (solution) La fonction de densité de probabilité et la fonction de répartition de la loi de Pareto de paramètres et sont, respectivement,

et

La fonction suivante simule un échantillon aléatoire de taille issu d'une distribution de Pareto de paramètres et  :

rpareto <- function(n, shape, scale)
scale * (runif(n)^(-1/shape) - 1)

1. Écrire une expression R utilisant la fonction rpareto ci-dessus qui permet de simuler cinq échantillons aléatoires de tailles 100, 150, 200, 250 et 300 d'une loi de Pareto avec et . Les échantillons aléatoires devraient être stockés dans une liste.

2. On vous donne l'exemple suivant d'utilisation de la fonction paste :

> paste("a", 1:5, sep = "")
[1] "a1" "a2" "a3" "a4" "a5"

Nommer les éléments de la liste créée en a) sample1…, sample5.

3. Calculer la moyenne de chacun des échantillons aléatoires obtenus en a). Retourner le résultat dans un vecteur.

4. Évaluer la fonction de répartition de la loi de Pareto(2, 5000) en chacune des valeurs de chacun des échantillons aléatoires obtenus en a). Retourner les valeurs de la fonction de répartition en ordre croissant.

5. Faire l'histogramme des données du cinquième échantillon aléatoire avec la fonction hist.

6. Ajouter 1000 à toutes les valeurs de tous les échantillons simulés en a), ceci afin d'obtenir des observations d'une distribution de Pareto translatée.

6.4 (solution) Une base de données contenant toutes les informations sur les assurés est stockée dans une liste de la façon suivante :

1. > x[[1]]

2. $num.police

3. [1] 1001

5. $franchise

6. [1] 500

8. $nb.acc

9. [1] 0 1 1 0

11. $montants

12. [1] 3186.864 3758.389

14. > x[[2]]

15. $num.police

16. [1] 1002

18. $franchise

19. [1] 250

21. $nb.acc

22. [1] 4 0 0 4 1 1 0

24. $montants

25. [1]  16728.7354  1414.7264 1825.7495  282.5609

26. [5]   1684.6686 14869.1731 7668.4196 2501.7257

27. [9] 108979.3725  2775.3161

Ainsi, x[[i]] contient les informations relatives à l'assuré i. Sans utiliser de boucles, écrire des expressions ou des fonctions R qui permettront de calculer les quantités suivantes.

1. La franchise moyenne dans le portefeuille.

2. Le nombre annuel moyen de réclamations par assuré.

3. Le nombre total de réclamations dans le portefeuille.

4. Le montant moyen par accident dans le portefeuille.

5. Le nombre d'assurés n'ayant eu aucune réclamation.

6. Le nombre d'assurés ayant eu une seule réclamation dans leur première année.

7. La variance du nombre total de sinistres.

8. La variance du nombre de sinistres pour chaque assuré.

9. La probabilité empirique qu'une réclamation soit inférieure à x (un scalaire) dans le portefeuille.

10. La probabilité empirique qu'une réclamation soit inférieure à x (un vecteur) dans le portefeuille.

7 - Fonctions d'optimisation

Objectifs du chapitre

- Connaître et savoir utiliser les différentes fonctions de calcul de racines de R.

- Connaître et savoir utiliser les différentes fonctions d'optimisation de R.

- Savoir reformuler un problème d'optimisation en base logarithmique pour éviter les difficultés numériques.

Les méthodes de bissection, du point fixe, de Newton-Raphson et consorts permettent de résoudre des équations à une variable de la forme ou Il existe également des versions de ces méthodes pour les systèmes à plusieurs variables de la forme

De tels systèmes d'équations surviennent plus souvent qu'autrement lors de l'optimisation d'une fonction. Par exemple, en recherchant le maximum ou le minimum d'une fonction , on souhaitera résoudre le système d'équations

En inférence statistique, les fonctions d'optimisation sont fréquemment employées pour calculer numériquement des estimateurs du maximum de vraisemblance.

La grande majorité des suites logicielles de calcul comportent des outils d'optimisation de fonctions. Ce chapitre passe en revue les fonctions disponibles dans R. Comme pour les chapitres précédents, des exemples d'utilisation de chacune des fonctions se trouvent dans le code informatique de la section 7.4.

7-1 - Fonctions d'optimisation et de calcul de racines

Le système R compte un certain nombre de fonctions pour trouver numériquement le minimum ou le maximum d'une fonction ainsi que pour calculer la racine d'une fonction dans un intervalle ou toutes les racines d'un polynôme. Ces fonctions diffèrent par certaines de leurs fonctionnalités et leur interface, mais aussi par les algorithmes utilisés. Consulter les rubriques d'aide pour les détails.

7-1-1 - Fonction uniroot

La fonction uniroot recherche la racine d'une fonction dans un intervalle. C'est donc la fonction de base pour trouver la solution (unique) de l'équation dans un intervalle déterminé.

7-1-2 - Fonction optimize

La fonction optimize recherche le minimum local (par défaut) ou le maximum local d'une fonction dans un intervalle donné.

7-1-3 - Fonction nlm

La fonction nlm minimise une fonction non linéaire sur un nombre arbitraire de paramètres.

7-1-4 - Fonction nlminb

La fonction nlminb est similaire à nlm, sauf qu'elle permet de spécifier des bornes inférieure ou supérieure pour les paramètres. Attention, toutefois : les arguments de la fonction ne sont ni les mêmes, ni dans le même ordre que ceux de nlm.

7-1-5 - Fonction optim

La fonction optim est l'outil d'optimisation tout usage de R. À ce titre, la fonction est souvent utilisée par d'autres fonctions. Elle permet de choisir parmi plusieurs algorithmes d'optimisation différents et, selon l'algorithme choisi, de fixer des seuils minimum ou maximum aux paramètres à optimiser.

7-1-6 - polyroot

En terminant, un mot sur polyroot(), qui n'est pas à proprement parler une fonction d'optimisation, mais qui pourrait être utilisée dans ce contexte. La fonction polyroot calcule toutes les racines (complexes) du polynôme . Le premier argument est le vecteur des coefficients dans cet ordre.

7-2 - Astuce Ripley

Brian Ripley — un important développeur de R — a publié le truc suivant dans les forums de discussion de R. Puisqu'il est très utile, nous nous permettons de le disséminer.

Une application statistique fréquente de l'optimisation est la maximisation numérique d'une fonction de vraisemblance ou, plus communément, la minimisation de la log-vraisemblance négative

Les fonctions d'optimisation sont d'ailleurs illustrées dans ce contexte dans le code informatique de la section 7.4.

Plusieurs lois de probabilité ont des paramètres strictement positifs. Or, en pratique, il n'est pas rare que les fonctions d'optimisation s'égarent dans les valeurs négatives des paramètres. La fonction de densité n'étant pas définie, la log-vraisemblance vaut alors NaN et cela peut faire complètement dérailler la procédure d'optimisation ou, à tout le moins, susciter des doutes sur la validité de la réponse.

Afin de pallier ce problème, l'Astuce Ripley™ propose d'estimer non pas les paramètres de la loi eux-mêmes, mais plutôt leur logarithme. Si l'on définit alors on peut écrire la fonction de log-vraisemblance ci-dessus sous la forme

Dès lors, (qui peut représenter un ou plusieurs paramètres) demeure valide sur tout l'axe des réels, ce qui permet d'éviter bien des soucis de nature numérique lors de la minimisation de .

Évidemment, le résultat de l'optimisation est l'estimateur du maximum de vraisemblance de . Il faudra donc veiller à faire la transformation inverse ̃ pour retrouver l'estimateur de .

L'utilisation de l'astuce est illustrée à la section 7.4.

7-3 - Pour en savoir plus

Les packages disponibles sur CRAN fournissent plusieurs autres outils d'optimisation pour R. Pour un bon résumé des options disponibles, consulter la CRAN Task View consacrée à l'optimisation : http://cran.r-project.org/web/views/Optimization.html.

7-4 - Exemples

1. ###

2. ### FONCTION 'uniroot'

3. ###

5. ## La fonction 'uniroot' recherche la racine d'une fonction

6. ## 'f' dans un intervalle spécifié soit comme une paire de

7. ## valeurs dans un argument 'interval', soit via des arguments

8. ## 'lower' et 'upper'.

9. ##

10. ## On calcule la solution de l'équation x - 2^(-x) = 0 dans

11. ## l'intervalle [0, 1].

12. f <- function(x) x - 2^(-x)      # fonction

13. uniroot(f, c(0, 1))              # appel simple

14. uniroot(f, lower = 0, upper = 1) # équivalent

16. ## On peut aussi utiliser 'uniroot' avec une fonction anonyme.

17. uniroot(function(x) x - 2^(-x), lower = 0, upper = 1)

19. ###

20. ### FONCTION 'optimize'

21. ###

23. ## On cherche le maximum local de la densité d'une loi bêta

24. ## dans l'intervalle (0, 1), son domaine de définition. (Ce

25. ## problème est facile à résoudre explicitement.)

26. ##

27. ## Les arguments de 'optimize' sont essentiellement les mêmes

28. ## que ceux de 'uniroot'. Ici, on utilise aussi l'argument

29. ## '... ' pour passer les paramètres de la loi bêta à 'dbeta'.

30. ##

31. ## Par défaut, la fonction recherche un minimum. Il faut donc

32. ## lui indiquer de rechercher plutôt un maximum.

33. optimize(dbeta, interval = c(0, 1), maximum = TRUE,

34. shape1 = 3, shape2 = 2)

36. ## On pourrait aussi avoir recours à une fonction auxiliaire.

37. ## Moins élégant et moins flexible.

38. f <- function(x) dbeta(x, 3, 2)

39. optimize(f, lower = 0, upper = 1, maximum = TRUE)

41. ###

42. ### FONCTION 'nlm'

43. ###

45. ## Pour la suite, nous allons donner des exemples

46. ## d'utilisation des fonctions d'optimisation dans un contexte

47. ## d'estimation des paramètres d'une loi gamma par la méthode

48. ## du maximum de vraisemblance.

49. ##

50. ## On commence par se donner un échantillon aléatoire de la

51. ## loi. Évidemment, pour ce faire nous devons connaître les

52. ## paramètres de la loi. C'est un exemple fictif.

53. set.seed(1) # toujours le même échantillon

54. x <- rgamma(10, 5, 2)

56. ## Les estimateurs du maximum de vraisemblance des paramètres

57. ## 'shape' et 'rate' de la loi gamma sont les valeurs qui

58. ## maximisent la fonction de vraisemblance

59. ##

60. ##     prod(dgamma(x, shape, rate))

61. ##

62. ## ou, de manière équivalente, qui minimisent la fonction de

63. ## log-vraisemblance négative

64. ##

65. ##     -sum(log(dgamma(x, shape, rate))).

66. ##

67. ## On remarquera au passage que les fonctions de calcul de

68. ## densités de lois de probabilité dans R ont un argument

69. ## 'log' qui, lorsque TRUE, retourne la valeur du logarithme

70. ## (naturel) de la densité de manière plus précise qu'en

71. ## prenant le logarithme après coup. Ainsi, pour faire le

72. ## calcul ci-dessus, on optera plutôt, pour l'expression

73. ##

74. ##     -sum(dgamma(x, shape, rate, log = TRUE))

75. ##

76. ## La fonction 'nlm' suppose que la fonction à optimiser

77. ## passée en premier argument a elle-même comme premier

78. ## argument le vecteur 'p' des paramètres à optimiser. Le

79. ## second argument de 'nlm' est un vecteur de valeurs de

80. ## départ, une pour chaque paramètre.

81. ##

82. ## Ainsi, pour trouver les estimateurs du maximum de

83. ## vraisemblance avec la fonction 'nlm' pour l'échantillon

84. ## ci-dessus, on doit d'abord définir une fonction auxiliaire

85. ## conforme aux attentes de 'nlm' pour calculer la fonction de

86. ## log-vraisemblance (à un signe près).

87. f <- function(p, x) -sum(dgamma(x, p[1], p[2], log = TRUE))

89. ## L'appel de 'nlm' est ensuite tout simple. Remarquer comment

90. ## on passe notre échantillon aléatoire (contenu dans l'objet

91. ## 'x') comme second argument à 'f' via l'argument '... ' de

92. ## 'nlm'. Le fait que l'argument de 'f' et l'objet contenant

93. ## les valeurs portent le même nom est sans importance. R sait

94. ## faire la différence entre l'un et l'autre.

95. nlm(f, c(1, 1), x = x)

97. ## === ASTUCE RIPLEY ===

98. ## L'optimisation ci-dessus a généré des avertissements ? C'est

99. ## parce que la fonction d'optimisation s'est égarée dans les

100. ## valeurs négatives, alors que les paramètres d'une gamma

101. ## sont strictement positifs. Cela arrive souvent en pratique

102. ## et cela peut faire complètement dérailler la procédure

103. ## d'optimisation (c'est-à-dire : pas de convergence).

104. ##

105. ## L'Astuce Ripley consiste à pallier ce problème en

106. ## estimant plutôt les logarithmes des paramètres. Pour ce

107. ## faire, il s'agit de réécrire la log-vraisemblance comme une

108. ## fonction du logarithme des paramètres, mais de la calculer

109. ## avec les véritables paramètres.

110. f2 <- function(logp, x)

111. {

112. p <- exp(logp) # retour aux paramètres originaux

113. -sum(dgamma(x, p[1], p[2], log = TRUE))

114. }

115. nlm(f2, c(0, 0), x = x)

117. ## Les valeurs obtenues ci-dessus sont toutefois les

118. ## estimateurs des logarithmes des paramètres de la loi gamma.

119. ## On retrouve les estimateurs des paramètres en prenant

120. ## l'exponentielle des réponses.

121. exp(nlm(f2, c(0, 0), x = x)$estimate)

122. ## ====================

124. ###

125. ### FONCTION 'nlminb'

126. ###

128. ## L'utilisation de la fonction 'nlminb' peut s'avérer

129. ## intéressante dans notre contexte puisque l'on sait que les

130. ## paramètres d'une loi gamma sont strictement positifs.

131. nlminb(c(1, 1), f, x = x, lower = 0, upper = Inf)

133. ###

134. ### FONCTION 'optim'

135. ###

137. ## La fonction 'optim' est très puissante, mais requiert aussi

138. ## une bonne dose de prudence. Ses principaux arguments sont :

139. ##

140. ## par : un vecteur contenant les valeurs initiales des

141. ##      paramètres ;

142. ## fn : la fonction à minimiser. Le premier argument de fn

143. ##     doit être le vecteur des paramètres.

144. ##

145. ## Comme pour les autres fonctions étudiées ci-dessus, on peut

146. ## passer des arguments à 'fn' (les données, par exemple) par

147. ## le biais de l'argument '... ' de 'optim'.

148. optim(c(1, 1), f, x = x)

150. ## L'estimation par le maximum de

151. ## vraisemblance\index{vraisemblance} est de beaucoup

152. ## simplifiée par l'utilisation de la fonction

153. ## \fonction{fitdistr} du package

154. ## \texttt{MASS}\index{package! MASS@\texttt{MASS}}.

156. ###

157. ### FONCTION 'polyroot'

158. ###

160. ## Racines du polynôme x^3 + 4 x^2 - 10. Les réponses sont

161. ## données sous forme de nombre complexe. Utiliser les

162. ## fonctions 'Re' et 'Im' pour extraire les parties réelles et

163. ## imaginaires des nombres, respectivement.

164. polyroot(c(-10, 0, 4, 1))     # racines

165. Re(polyroot(c(-10, 0, 4, 1))) # parties réelles

166. Im(polyroot(c(-10, 0, 4, 1))) # parties imaginaires

7-5 - Exercices

7.1 (solution) Trouver la solution des équations suivantes à l'aide des fonctions R appropriées.

1. pour

2. pour

3. pour

4. pour

5. pour

7.2 (solution) En théorie de la crédibilité, l'estimateur d'un paramètre est donné sous forme de point fixe

où

et , et sont des données. Calculer la valeur de si et que les valeurs de et sont telles qu'elles apparaissent dans le tableau ci-dessous.

1

2

3

4

5

2061

1511

1806

1353

1600

100 155

19 895

13 735

4 152

36 110

7.3 (solution) Les fonctions de densité de probabilité et de répartition de la distribution de Pareto sont données à l'exercice 6.3. Calculer les estimateurs du maximum de vraisemblance des paramètres de la Pareto à partir d'un échantillon aléatoire obtenu par simulation avec la commande

> x <- lambda * (runif(100)^(-1/alpha) - 1)

pour des valeurs de alpha et lambda choisies.

8 - Générateurs de nombres aléatoires

Objectifs du chapitre

- Générer des nombres aléatoires uniformes avec la fonction runif.

- Générer des nombres aléatoires non uniformes provenant de lois de probabilité discrètes et continues.

- Générer des nombres aléatoires provenant d'une distribution discrète quelconque.

- Tirer profit de la nature vectorielle des fonctions de simulation de R.

Avant d'utiliser pour quelque tâche de simulation moindrement importante un générateur de nombres aléatoires inclus dans un logiciel, il importe de s'assurer de la qualité de cet outil. On trouvera en général relativement facilement de l'information dans Internet.

On présente ici, sans entrer dans les détails, les générateurs de nombres uniformes utilisés dans R ainsi que la liste des différentes fonctions de simulation de variables aléatoires.

8-1 - Générateurs de nombres aléatoires

On obtient des nombres uniformes sur un intervalle quelconque avec la fonction runif dans R. La fonction set.seed permet de spécifier la valeur de l'amorce du générateur aléatoire, ce qui est utile si on veut répéter une simulation absolument à l'identique.

R offre la possibilité de choisir entre plusieurs générateurs de nombres aléatoires différents, ou encore de spécifier son propre générateur. Par défaut, R utilise le générateur Marsenne-Twister, considéré comme le plus avancé en ce moment. La période de ce générateur est 219 937 − 1 (rien de moins !) et la distribution des nombres est uniforme dans 623 dimensions consécutives sur toute la période.

Pour de plus amples détails et les dernières informations sur les générateurs disponibles et la procédure de réglage de l'amorce, consulter les rubriques d'aide des fonctions .Random.seed et set.seed.

8-2 - Fonctions de simulation de variables aléatoires non uniformes

Un large éventail de fonctions donnent directement accès aux caractéristiques de plusieurs lois de probabilité dans R. Pour chaque racine loi, il existe quatre fonctions différentes :

1. dloi calcule la fonction de densité de probabilité (loi continue) ou la fonction de masse de probabilité (loi discrète) ;

2. ploi calcule la fonction de répartition ;

3. qloi calcule la fonction de quantile ;

4. rloi simule des observations de cette loi.

Les différentes lois de probabilité disponibles dans le système R de base, leur racine et le nom de leurs paramètres sont rassemblés au tableau 8.1. Des packages fournissent des fonctions pour d'autres lois dont, entre autres, actuar (Dutang et collab., 2008) et SuppDists (Wheeler, 2013).

Tab. 8.1 - Lois de probabilité pour lesquelles il existe des fonctions dans le système R de base

Loi de probabilité

Racine dans R

Noms des paramètres

Bêta

beta

shape1, shape2

Binomiale

binom

size, prob

Binomiale négative

nbinom

size, prob ou mu

Cauchy

cauchy

location, scale

Exponentielle

exp

rate

F (Fisher)

f

df1, df2

Gamma

gamma

shape, rate ou scale

Géométrique

geom

prob

Hypergéométrique

hyper

m, n, k

Khi carré

chisq

df

Logistique

logis

location, scale

Log-normale

lnorm

meanlog, sdlog

Normale

norm

mean, sd

Poisson

pois

lambda

t (Student)

t

df

Uniforme

unif

min, max

Weibull

weibull

shape, scale

Wilcoxon

wilcox

m, n

Toutes les fonctions du tableau 8.1 sont vectorielles, c'est-à-dire qu'elles acceptent en argument un vecteur de points où la fonction (de densité, de répartition ou de quantile) doit être évaluée et même un vecteur de paramètres. Par exemple,

> dpois(c(3, 0, 8), lambda = c(1, 4, 10))
[1] 0.06131324 0.01831564 0.11259903

retourne la probabilité que des lois de Poisson de paramètre 1, 4 et 10 prennent les valeurs 3, 0 et 8, dans l'ordre.

Le premier argument de toutes les fonctions de simulation est la quantité de nombres aléatoires désirée. Ainsi,

> rpois(3, lambda = c(1, 4, 10))
[1] 2 6 10

retourne trois nombres aléatoires issus de distributions de Poisson de paramètre 1, 4 et 10, respectivement. Évidemment, passer un vecteur comme premier argument n'a pas tellement de sens, mais, si c'est fait, R retournera une quantité de nombres aléatoires égale à la longueur du vecteur (sans égard aux valeurs contenues dans le vecteur).

La fonction sample permet de simuler des nombres d'une distribution discrète quelconque. Sa syntaxe est

sample(x, size, replace = FALSE, prob = NULL)

où x est un vecteur des valeurs possibles de l'échantillon à simuler (le support de la distribution), size est la quantité de nombres à simuler et prob est un vecteur de probabilités associées à chaque valeur de x (1/length(x) par défaut). Enfin, si replace est TRUE, l'échantillonnage se fait avec remise.

8-3 - Exemples

1. ###

2. ### GÉNÉRATEURS DE NOMBRES ALÉATOIRES

3. ###

5. ## La fonction de base pour simuler des nombres uniformes est

6. ## 'runif'.

7. runif(10) # sur (0, 1) par défaut

8. runif(10, 2, 5) # sur un autre intervalle

9. 2 + 3 * runif(10) # équivalent, moins lisible

11. ## R est livré avec plusieurs générateurs de nombres

12. ## aléatoires. On peut en changer avec la fonction 'RNGkind'.

13. RNGkind("Wichmann-Hill")    # générateur de Excel

14. runif(10)                   # rien de particulier à voir

15. RNGkind("default")          # retour au générateur par défaut

17. ## La fonction 'set.seed' est très utile pour spécifier

18. ## l'amorce d'un générateur. Si deux simulations sont

19. ## effectuées avec la même amorce, on obtiendra exactement les

20. ## mêmes nombres aléatoires et, donc, les mêmes résultats.

21. ## Très utile pour répéter une simulation à l'identique.

22. set.seed(1)                 # valeur sans importance

23. runif(5)                    # 5 nombres aléatoires

24. runif(5)                    # 5 autres nombres

25. set.seed(1)                 # réinitialisation de l'amorce

26. runif(5)                    # les mêmes 5 nombres que ci-dessus

28. ###

29. ### FONCTIONS POUR LA SIMULATION DE VARIABLES ALÉATOIRES NON

30. ### UNIFORMES

31. ###

33. ## Plutôt que de devoir utiliser la méthode de l'inverse ou un

34. ## autre algorithme de simulation pour obtenir des nombres

35. ## aléatoires d'une loi de probabilité non uniforme, R fournit

36. ## des fonctions de simulation pour bon nombre de lois. Toutes

37. ## ces fonctions sont vectorielles. Ci-dessous, P == Poisson

38. ## et G == Gamma pour économiser sur la notation.

39. n <- 10 # taille des échantillons

40. rbinom(n, 5, 0.3)           # Binomiale(5, 0,3)

41. rbinom(n, 1, 0.3)           # Bernoulli(0,3)

42. rnorm(n)                    # Normale(0, 1)

43. rnorm(n, 2, 5)              # Normale(2, 25)

44. rpois(n, c(2, 5))           # P(2), P(5), P(2), ... , P(5)

45. rgamma(n, 3, 2:11)          # G(3, 2), G(3, 3), ... , G(3, 11)

46. rgamma(n, 11:2, 2:11)       # G(11, 2), G(10, 3), ... , G(2, 11)

48. ## La fonction 'sample' sert pour simuler d'une distribution

49. ## discrète quelconque. Le premier argument est le support de

50. ## la distribution et le second, la taille de l'échantillon

51. ## désirée. Par défaut, l'échantillonnage se fait avec remise

52. ## et avec des probabilités égales sur tout le support.

53. sample(1:49, 7)             # numéros pour le 7/49

54. sample(1:10, 10)            # mélange des nombres de 1 à 10

56. ## On peut échantillonner avec remise.

57. sample(1:10, 10, replace = TRUE)

59. ## On peut aussi spécifier une distribution de probabilités

60. ## non uniforme.

61. x <- sample(c(0, 2, 5), 1000, replace = TRUE,

62. prob = c(0.2, 0.5, 0.3))

63. table(x)                    # tableau de fréquences

8-4 - Exercices

8.1 (solution) La loi log-normale est obtenue par transformation de la loi normale : si la distribution de la variable aléatoire est une normale de paramètres et , alors la distribution de est une log-normale. Simuler 1000 observations d'une loi log-normale de paramètres et , puis tracer l'histogramme de l'échantillon aléatoire obtenu.

8.2 (solution) Simuler 10 000 observations d'un mélange continu Poisson/gamma où les paramètres de la loi gamma sont et , puis tracer la distribution de fréquence de l'échantillon aléatoire obtenu à l'aide des fonctions plot et table. Superposer à ce graphique la fonction de probabilité d'une binomiale négative de paramètres et .

8.3 (solution) Simuler 10 000 observations d'un mélange discret de deux distributions log-normales, l'une de paramètres ( = 3 ,5 , = 0 ,6) et l'autre de paramètres () . Utiliser un paramètre de mélange . Tracer ensuite l'histogramme de l'échantillon aléatoire obtenu.

A - GNU Emacs et ESS : la base

Emacs est l'Éditeur de texte des éditeurs de texte. À l'origine un éditeur pour les programmeurs (avec des modes spéciaux pour une multitude de langages différents), Emacs est devenu au fil du temps un environnement logiciel en soi dans lequel on peut réaliser une foule de tâches différentes : rédiger des documents LaTeX, interagir avec R, SAS ou un logiciel de base de données, consulter son courrier électronique, gérer son calendrier ou même jouer à Tetris !

Cette annexe passe en revue les quelques commandes essentielles à connaître pour commencer à travailler avec GNU Emacs et le mode ESS. L'ouvrage de Cameron et collab. (2004) constitue une excellente référence pour l'apprentissage plus poussé de l'éditeur.

A-1 - Mise en contexte

Emacs est le logiciel étendard du projet GNU (« GNU is not Unix »), dont le principal commanditaire est la Free Software Foundation (FSF) à l'origine de tout le mouvement du logiciel libre.

- Richard M. Stallman, président de la FSF et grand apôtre du libre, a écrit la première version de Emacs et il continue à ce jour à contribuer au projet.

- Les origines de Emacs remontent au début des années 1980, une époque où les interfaces graphiques n'existaient pas, le parc informatique était beaucoup plus hétérogène qu'aujourd'hui (les claviers n'étaient pas les mêmes d'une marque d'ordinateur à une autre) et les modes de communication entre les ordinateurs demeuraient rudimentaires.

- L'âge vénérable de Emacs transparaît à plusieurs endroits, notamment dans la terminologie inhabituelle, les raccourcis clavier non conformes aux standards d'aujourd'hui ou la manipulation des fenêtres qui ne se fait pas avec une souris.

Tiré de XKCD.com. La commande M-X butterfly existe vraiment dans Emacs… en référence à cette bande dessinée !

Emacs s'adapte à différentes tâches par l'entremise de modes qui modifient son comportement ou lui ajoutent des fonctionnalités. L'un de ces modes est ESS (Emacs Speaks Statistics).

- ESS permet d'interagir avec des logiciels statistiques (en particulier R, S+ et SAS) directement depuis Emacs.

- Quelques-uns des développeurs de ESS sont aussi des développeurs de R, d'où la grande compatibilité entre les deux logiciels.

- Lorsque ESS est installé, le mode est activé automatiquement en ouvrant dans Emacs un fichier dont le nom se termine par l'extension .R.

A-2 - Installation

GNU Emacs et le mode ESS sont normalement livrés d'office avec toutes les distributions Linux. Pour les environnements Windows et OS X, le plus simple consiste à installer les distributions préparées par le présent auteur. Consulter le site http://vgoulet.act.ulaval.ca/emacs/

A-3 - Description sommaire

Au lancement, Emacs affiche un écran d'information contenant des liens vers différentes ressources. Cet écran disparaît dès que l'on appuie sur une touche. La fenêtre Emacs se divise en quatre zones principales (voir la figure A.1) :

Fig. A.1 - Fenêtre GNU Emacs et ses différentes parties au lancement de l'application sous OS X. Sous Windows et Linux, la barre de menu se trouve à l'intérieur de la fenêtre.

1. Tout au haut de la fenêtre (ou de l'écran sous OS X), on trouve l'habituelle barre de menu dont le contenu change selon le mode dans lequel se trouve Emacs ;

2. L'essentiel de la fenêtre sert à afficher un buffer, soit le contenu d'un fichier ouvert ou l'invite de commande d'un programme externe ;

3. La ligne de mode est le séparateur horizontal contenant diverses informations sur le fichier ouvert et l'état de Emacs ;

4. Le minibuffer est la région au bas de la fenêtre où l'on entre des commandes et reçoit de l'information de Emacs.

Il est possible de séparer la fenêtre Emacs en sous-fenêtres pour afficher plusieurs buffers à la fois. Il y a alors une ligne de mode pour chaque buffer.

A-4 - Emacs-ismes et Unix-ismes

Emacs possède sa propre terminologie qu'il vaut mieux connaître lorsque l'on consulte la documentation. De plus, l'éditeur utilise des conventions du monde Unix qui sont moins usitées sur les plateformes Windows et OS X.

- Dans les définitions de raccourcis clavier :

- C est la touche Contrôle () ;

- M est la touche Meta, qui correspond à la touche Alt de gauche sur un PC ou la touche Option () sur un Mac (toutefois, voir l'encadré) ;

- ESC est la touche Échap () et est équivalente à Meta ;

- SPC est la barre d'espacement ;

- DEL est la touche Retour arrière () — et non la toucheSupprimer ;

- RET est la touche Entrée ().

- Toutes les fonctionnalités de Emacs correspondent à une commande pouvant être tapée dans le minibuffer. M-x démarre l'invite de commande.

- Le caractère ~ représente le dossier vers lequel pointe la variable d'environnement $HOME (Linux, OS X) ou %HOME% (Windows). C'est le dossier par défaut de Emacs.

- La barre oblique (/) est utilisée pour séparer les dossiers dans les chemins d'accès aux fichiers, même sous Windows.

- En général, il est possible d'appuyer sur TAB dans le minibuffer pour compléter les noms de fichiers ou de commandes.

Par défaut sous OS X, la touche Meta est assignée à Option (). Sur les claviers français, cela empêche d'accéder à certains caractères spéciaux tels que [, ], { ou }.

Une solution à cette fâcheuse situation consiste à assigner la touche Meta à Commande (). Cela bloque alors l'accès à certains raccourcis Mac, mais la situation est moins critique ainsi. Pour assigner la touche Meta à Commande () et laisser la touche Option () jouer son rôle usuel, il suffit d'insérer les lignes suivantes dans son fichier de configuration .emacs (voir la section A.7) :

1. ;;; ====================================

2. ;;; Assigner la touche Meta à Commande

3. ;;; et laisser Option être Option

4. ;;; ====================================

5. (setq-default ns-command-modifier 'meta)

6. (setq-default ns-option-modifier 'none)

A-5 - Commandes de base

Emacs comporte une pléthore de commandes, il serait donc futile de tenter d'en faire une liste exhaustive ici. Nous nous contenterons de mentionner les commandes les plus importantes regroupées par tâche.

A-5-1 - Les essentielles

démarrer l'invite de commande.

bouton de panique : annuler, quitter ! Presser plus d'une fois au besoin.

A-5-2 - Manipulation de fichiers

Entre parenthèses, le nom de la commande Emacs correspondante. On peut entrer cette commande dans le minibuffer au lieu d'utiliser le raccourci clavier.

On remarquera qu'il n'existe pas de commande « nouveau fichier » dans Emacs. Pour créer un nouveau fichier, il suffit d'ouvrir un fichier n'existant pas.

ouvrir un fichier (find-file)

sauvegarder (save-buffer)

sauvegarder sous (write-file)

fermer un fichier (kill-buffer)

annuler (pratiquement illimité) ; aussi C-x u (undo)

recherche incrémentale avant (isearch-forward)

Recherche incrémentale arrière (isearch-backward)

rechercher et remplacer (query-replace)

A-5-3 - Déplacements simples du curseur

déplacer d'un caractère vers l'arrière | l'avant (backward-char | forward-char)

aller au début | fin de la ligne (move-beginning-of-line | move-end-of-line)

aller à la ligne précédente | suivante (previous-line | next-line)

aller au début | fin du fichier (beginning-of-buffer | end-of-buffer)

effacer le caractère à gauche | droite du curseur (delete-backward-char | delete-char)

effacer le mot à gauche | droite du curseur (backward-kill-word | kill-word)

supprimer jusqu'à la fin de la ligne (kill-line)

Plusieurs des raccourcis clavier de Emacs composés avec la touche Contrôle () sont valides sous OS X. Par exemple, A et E déplacent le curseur au début et à la fin de la ligne dans les champs texte.

A-5-4 - Sélection de texte, copier, coller, couper

débute la sélection (set-mark-command)

couper la sélection (kill-region)

copier la sélection (kill-ring-save)

coller (yank)

remplacer le dernier texte collé par la sélection précédente (yank-pop)

- Il est possible d'utiliser les raccourcis clavier usuels de Windows (C-c, C-x, C-v) et OS X (C, X, V) en activant le mode CUA dans le menu Options.

- On peut copier-coller directement avec la souris dans Windows en sélectionnant du texte puis en appuyant sur le bouton central (ou la molette) à l'endroit souhaité pour y copier le texte.

A-5-5 - Manipulation de fenêtres

changer de buffer (switch-buffer)

séparer l'écran en deux fenêtres (split-window-vertically)

conserver uniquement la fenêtre courante (delete-other-windows)

fermer la fenêtre courante (delete-window)

aller vers une autre fenêtre lorsqu'il y en a plus d'une (other-window)

A-5-6 - Manipulation de fihiers de script dans le mode ESS

Le mode ESS dispose de fonctions « intelligentes » qui facilitent grandement la manipulation des fichiers de script. Les deux principales commandes à connaître sont les suivantes :

évaluer dans le processus R la ligne sous le curseur ou la région sélectionnée, puis déplacer le curseur à la prochaine expression (ess-eval-region-or-line-and-step) ;

évaluer dans le processus R la région sélectionnée, la fonction ou le paragraphe (tout bloc entre deux lignes blanches) dans lequel se trouve le curseur, puis déplacer le curseur à la prochaine expression

(ess-eval-region-or-function-or-paragraph-and-step).

Les quelques autres fonctions utiles sont :

déplacer le curseur vers le processus R (ess-switch-to-inferior-or-script-buffer) ;

évaluer le code de la fonction courante dans le processus R (ess-eval-function) ;

évaluer le code du fichier courant en entier dans le processus R (ess-load-file) ;

aide sur une commande R (ess-display-help-on-object).

A-5-7 - Interaction avec l'invite de commande R

commande précédente | suivante dans l'historique (previous-matching-history-from-input | next-matching-history-from-input).

déplacer le curseur vers le fichier de script courant (ess-switch-to-inferior-or-script-buffer).

sélectionner le résultat de la dernière commande (mark-paragraph).

effacer le résultat de la dernière commande (comint-delete-output).

aide sur une commande R (ess-display-help-on-object).

terminer le processus R (ess-quit).

A-5-8 - Consultation des rubriques d'aide de R

aller à la section précédente | suivante de la rubrique (ess-skip-to-previous-section | ess-skip-to-next-section).

aller à la section de la liste des arguments (Arguments).

aller à la section des détails sur la fonction (Details).

aller à la section sur la valeur retournée par la fonction (Value).

aller à la section des fonctions apparentée (See Also).

aller à la section des exemples (Examples).

évaluer la ligne sous le curseur ; pratique pour exécuter les exemples (ess-eval-line-and-step).

évaluer la région sélectionnée (ess-eval-region).

ouvrir une nouvelle rubrique d'aide, par défaut pour le mot se trouvant sous le curseur (ess-display-help-on-object).

retourner au processus ESS en laissant la rubrique d'aide visible (ess-switch-to-end-of-ESS).

fermer la rubrique d'aide et retourner au processus ESS (ess-kill-buffer-and-go).

A-6 - Anatomie d'une session de travail (bis)

On reprend ici les étapes d'une session de travail type présentées à la section 1.6, mais en expliquant comment compléter chacune dans Emacs avec le mode ESS.

1. Lancer Emacs et ouvrir un fichier de script avec

1. C-x C-f

ou avec le menu

1. File|Open File…

En spécifiant un nom de fichier qui n'existe pas déjà, on se trouve à créer un nouveau fichier de script. S'assurer de terminer le nom du nouveau fichier par .R pour que Emacs reconnaisse automatiquement qu'il s'agit d'un fichier de script R.

2. Démarrer un processus R à l'intérieur même de Emacs avec

1. M-x R

Emacs demandera alors de spécifier le répertoire de travail (starting data directory). Accepter la valeur par défaut, par exemple

1. ~/

ou indiquer un autre dossier. Un éventuel message de Emacs à l'effet que le fichier .Rhistory n'a pas été trouvé est sans conséquence et peut être ignoré.

3. Composer le code. Lors de cette étape, on se déplacera souvent du fichier de script à la ligne de commande afin d'essayer diverses expressions. On exécutera également des parties seulement du code se trouvant dans le fichier de script. Les commandes les plus utilisées sont alors :

1. C-RET pour exécuter une ligne du fichier de script ou la région sélectionnée ;

2. C-c C-c pour exécuter un paragraphe du fichier de script ;

3. C-c C-z pour se déplacer entre le fichier de script et la ligne de commande R.

4. Sauvegarder le fichier de script :

1. C-x C-s

Les quatrième et cinquième caractères de la ligne de mode changent de ** à --.

5. Sauvegarder si désiré l'espace de travail de R avec save.image(). Cela n'est habituellement pas nécessaire à moins que l'espace de travail ne contienne des objets importants ou longs à recréer.

6. Quitter le processus R avec

1. C-c C-q

Cette commande ESS se chargera de fermer tous les fichiers associés au processus R. On peut ensuite quitter Emacs en fermant l'application de la manière usuelle.

A-7 - Configuration de l'éditeur

Une des grandes forces de Emacs est qu'à peu près chacune de ses facettes est configurable : couleurs, polices de caractères, raccourcis clavier, etc.

- La configuration de Emacs se fait par le biais de commandes réunies dans un fichier de configuration nommé .emacs (le point est important !) que Emacs lit au démarrage.

- Le fichier .emacs doit se trouver dans le dossier ~/, c'est-à-dire dans le dossier de départ de l'utilisateur sous Linux et OS X, et dans le dossier référencé par la variable d'environnement %HOME% sous Windows.

A-8 - Aide et documentation

Emacs possède son propre système d'aide très exhaustif, mais dont la navigation est peu intuitive selon les standards d'aujourd'hui. Consulter le menu Help.

Autrement, on trouvera les manuels de Emacs et de ESS en divers formats dans les sites respectifs des deux projets :

- http://www.gnu.org/software/emacs ;

- http://ess.r-project.org.

Enfin, si le désespoir vous prend au cours d'une séance de codage intensive, vous pouvez toujours consulter le psychothérapeute Emacs. On le trouve, bien entendu, dans le menu Help !

B - RStudio : une introduction

Un environnement de développement intégré (integrated development environment, IDE) est un progiciel de productivité destiné au développement de logiciels ou, plus largement, à la programmation informatique. Il comprend toujours un éditeur de texte adapté au langage de programmation visé, un environnement de compilation ou d'exécution du code et, généralement, des outils de contrôle de versions, de gestion des projets et de navigation dans le code source(5).

Offert au public depuis 2011, RStudio est un IDE convivial conçu spécifiquement pour l'analyse de données et le développement de packages avec R. Il est produit par RStudio Inc. et est offert en version libre ou commerciale, pour une exécution locale (desktop) ou pour une exécution sur un serveur via un navigateur web.

B-1 - Installation

RStudio est disponible à l'identique pour les plateformes Windows, OS X et Linux. Pour une utilisation locale sur son poste de travail, on installera la version libre (Open Source) de RStudio Desktop depuis le site https://www.rstudio.com/products/rstudio/download/.

B-2 - Description sommaire

La fenêtre de RStudio se divise toujours en quatre sous-fenêtres(6) — sauf au lancement, alors que la sous-fenêtre d'édition de code source n'est pas visible ; voir la figure B.1. Dans le sens des aiguilles d'une montre en partant en haut à gauche, on trouve :

Fig. B.1 - Fenêtre de RStudio et trois de ses sous-fenêtres au lancement de l'application sous OS X. Sous Windows et Linux, la fenêtre comporte également une barre de menu.

1. La sous-fenêtre d'édition de code source, avec un onglet par fichier de script ;

2. Le navigateur d'environnement de travail ou d'historique des commandes, selon l'onglet sélectionné ;

3. Le navigateur de fichiers du projet, de packages, de graphiques, etc., selon l'onglet sélectionné ;

4. La console — ou ligne de commande — R.

Au lancement de l'application, la console R occupe toute la gauche de la fenêtre jusqu'à ce qu'un fichier de script soit ouvert.

- Le navigateur d'environnement de travail est particulièrement utile pour voir le contenu, les attributs, le type et la taille de chaque objet sauvegardé dans la session R. Il permet également de visualiser le contenu des objets en cliquant sur leur nom ou sur l'icône de grille à droite de leur nom.

- Il ne peut y avoir qu'un seul processus R (affiché dans la console) actif par fenêtre RStudio. Pour utiliser plusieurs processus R simultanément, il faut démarrer autant de copies de Rstudio.

- La position des sous-fenêtres dans la grille ne peut être modifiée. Par contre, chaque sous-fenêtre peut être redimensionnée.

- On peut modifier la liste des onglets affichés dans les deux navigateurs dans les préférences de l'application ; voir la section B.6.

B-3 - Projets

Il est possible d'utiliser RStudio un peu comme un simple éditeur de texte.

- On ouvre les fichiers de scripts un à un, soit à partir du menu File|Open file…, soit à partir de l'onglet Files du navigateur de fichiers.

- Lorsque nécessaire, on change le répertoire de travail de R à partir du menu Session.

Pour faciliter l'organisation de son travail, l'ouverture des fichiers de script et le lancement d'un processus R dans le bon répertoire de travail, RStudio propose la notion de projet.

- Un projet RStudio est associé à un répertoire de travail de R (section 1.7).

- On crée un nouveau projet à partir du menu Project à l'extrémité droite de la barre d'outils ou à partir du menu File|New Project… On a alors l'option de créer un nouveau dossier sur notre poste de travail ou de créer un projet à partir d'un dossier existant.

- Lors de la création d'un projet, RStudio crée dans le dossier visé un fichier avec une extension .Rproj contenant diverses informations en lien avec le projet. De plus, le projet est immédiatement chargé dans Rstudio.

- L'ouverture d'un projet entraîne : le lancement d'une session R ayant comme répertoire de travail le dossier du projet ; le chargement du fichier .RData (le cas échéant) ; l'ouverture de tous les fichiers de scripts qui étaient ouverts lors de la dernière séance de travail.

- Chaque projet dispose de ses propres réglages. On accède à ceux-ci via la commande Project Options… du menu Project de la barre d'outils.

On trouvera plus d'informations sur les projets dans l'aide en ligne de RStudio.

B-4 - Commandes de base

Comme l'interface de RStudio respecte les standards modernes, nous ne soulignons ici que les commandes particulièrement utiles pour la manipulation des fichiers de script. On accède rapidement à la liste des commandes les plus utiles via le menu Help de l'application.

Les raccourcis clavier sous, d'une part, Windows et Linux et sous, d'autre part, OS X légèrement différents. Nous fournissons les deux jeux ci-dessous, séparés par le symbole •.

insérer le symbole d'assignation ␣<-␣

évaluer dans le processus R la ligne sous le curseur ou la région sélectionnée, puis déplacer le curseur à la prochaine expression.

évaluer le code du fichier courant en entier dans le processus R.

évaluer dans le processus R le code source du début du fichier jusqu'à la ligne sous le curseur.

évaluer dans le processus R le code source de la ligne sous le curseur jusqu'à la fin du fichier.

évaluer le code de la fonction courante dans le processus R

À la console — ou ligne de commande — R, les raccourcis suivants sont particulièrement utiles.

expression précédente | suivante dans l'historique des commandes.

afficher la fenêtre d'historique des commandes.

Le très pratique raccourci clavier servant à insérer le symbole d'assignation ne fonctionne pas avec le clavier canadien français puisque cette combinaison de touches sert déjà à insérer le symbole |.

Au moment d'écrire ces lignes, RStudio ne permet pas de réaffecter la commande d'insertion du symbole d'assignation à une autre combinaison de touches.

Une solution de rechange consiste à utiliser une disposition de clavier anglaise pour travailler dans RStudio.

Pour ce faire, accéder aux préférences système de OS X puis sélectionner l'option Clavier. Dans l'onglet Méthodes de saisie, installer un nouveau clavier anglais. Cocher l'option « Afficher le menu de saisie dans la barre des menus » pour pouvoir rapidement et facilement passer d'un type de clavier à un autre.

B-5 - Anatomie d'une session de travail (ter)

On reprend ici les étapes d'une session de travail type présentées à la section 1.6, mais en expliquant comment compléter chacune dans RStudio.

1. Lancer RStudio et ouvrir soit un nouveau fichier de script avec

1. Ctrl+Shift+N N

ou avec le menu

1. File|New File|R Script… ,

soit un fichier de script existant avec

1. Ctrl+O O

ou

1. File|Open File…

2. C'est une bonne pratique de faire du dossier où se trouve son ou ses fichiers de scripts le répertoire de travail de R. Il suffit de sélectionner le menu

Session|Set Working Directory|To Source File Location.

3. Composer le code. Lors de cette étape de programmation, on se déplacera souvent du fichier de script à la ligne de commande afin d'essayer diverses expressions. On exécutera également des parties seulement du code se trouvant dans le fichier de script. Les commandes les plus utilisées sont alors :

1. Ctrl+Retour pour exécuter une ligne ou une région sélectionnée du fichier de script ;

2. Ctrl+1 1 pour déplacer le curseur vers la sous-fenêtre d'édition de script ;

3. Ctrl+2 2 pour déplacer le curseur vers la console R.

4. Sauvegarder le fichier de script :

1. Ctrl+S S

(S'il s'agit d'un nouveau fichier, s'assurer de terminer son nom par .R.) Le nom du fichier dans l'onglet de la sous-fenêtre passe du rouge au noir.

5. Sauvegarder si désiré l'espace de travail de R avec save.image(). Cela n'est habituellement pas nécessaire à moins que l'espace de travail ne contienne des objets importants ou longs à recréer.

6. Quitter RStudio de la manière usuelle. Par défaut, RStudio devrait demander si l'on souhaite sauvegarder l'espace de travail de R. Nous suggérons de ne pas le faire.

La section B.6 explique comment configurer RStudio afin d'éviter de se faire poser la question à chaque fermeture de l'application.

B-6 - Configuration de l'éditeur

Il est possible de configurer plusieurs facettes de RStudio à partir d'une interface familière.

- On accède aux options de configuration par le menu

1. Tools|Global Options…

sous Windows et Linux et par le menu standard

1. RStudio|Preferences (,)

sous OS X.

- Nous suggérons de régler l'option Save workspace to .RData on exit à Never dans les options de configuration générales ; voir la figure B.2. Avec ce réglage, l'espace de travail de R ne sera pas sauvegardé à la fermeture de Rstudio.

Fig. B.2 - Réglage suggéré de RStudio (encadré en jaune) faisant en sorte que l'espace de travail de R n'est jamais sauvegardé au moment de quitter l'application

B-7 - Aide et documentation

La documentation de RStudio se trouve entièrement en ligne. On y accède par le menu Help. L'onglet Help du navigateur de fichiers (sous-fenêtre en bas à droite) offre également une interface unifiée pour accéder à l'aide de R et à celle de RStudio.

C - Planification d'une simulation en R

Il existe de multiples façons de réaliser la mise en œuvre informatique d'une simulation, mais certaines sont plus efficaces que d'autres. Cette annexe passe en revue diverses façons de faire des simulations avec R à l'aide d'un exemple simple de nature statistique.

C-1 - Contexte

Soit un échantillon aléatoire tiré d'une population distribuée selon une loi uniforme sur l'intervalle . On considère trois estimateurs sans biais du paramètre inconnu :

1. La moyenne arithmétique

2. La médiane empirique

où est la ke statistique d'ordre de l'échantillon aléatoire ;

3. La mi-étendue

À l'aide de la simulation on veut, d'une part, vérifier si les trois estimateurs sont bel et bien sans biais et, d'autre part, déterminer lequel a la plus faible variance.

Pour ce faire, on doit d'abord simuler un grand nombre d'échantillons aléatoires de taille d'une distribution pour une valeur de choisie. Pour chaque échantillon, on calculera ensuite les trois estimateurs ci-dessus, puis la moyenne et la variance, par type d'estimateur, de tous les estimateurs obtenus. Si la moyenne des estimateurs est près de alors on pourra conclure que est sans biais. De même, on déterminera lequel des trois estimateurs a la plus faible variance selon le classement des variances empiriques.

C-2 - Première approche : avec une boucle

La façon la plus intuitive de mettre en œuvre cette étude de simulation en R consiste à utiliser une boucle for. Avec cette approche, il est nécessaire d'initialiser une matrice de 3 lignes et colonnes (ou l'inverse) dans laquelle seront stockées les valeurs des trois estimateurs pour chaque simulation. Une fois la matrice remplie dans la boucle, il ne reste plus qu'à calculer la moyenne et la variance par ligne pour obtenir les résultats souhaités.

La figure C.1 présente un exemple de code adéquat pour réaliser la simulation à l'aide d'une boucle.

Si l'on souhaite pouvoir exécuter le code de la figure C.1 facilement à l'aide d'une seule expression, il suffit de placer l'ensemble du code dans une fonction. La fonction simul1 de la figure C.2 reprend le code de la figure C.1, sans les commentaires. On a alors :

1. > simul1(10000, 100, 0)

2. $biais

3. Moyenne       Mediane   Mi-etendue

4. 6.945790e-05 -5.293742e-04 9.199309e-05

6. $variances

7. Moyenne      Mediane   Mi-etendue

8. 8.128808e-04 2.396603e-03 4.918073e-05

C-3 - Seconde approche : avec sapply

On le sait, les boucles sont inefficaces en R. Il est en général plus efficace de déléguer les boucles aux fonctions lapply et sapply (section 7.3). On rappelle que la syntaxe de ces fonctions est

Fig. C.1 - Code pour la simulation utilisant une boucle for

1. ## Bonne habitude à prendre : stocker les constantes dans

2. ## des variables faciles à modifier au lieu de les écrire

3. ## explicitement dans le code.

4. size <- 100             # taille de chaque échantillon

5. nsimul <- 10000         # nombre de simulations

6. theta <- 0              # la valeur du paramètre

8. ## Les lignes ci-dessous éviteront de faire deux additions

9. ## 'nsimul' fois.

10. a <- theta - 0.5        # borne inférieure de l'uniforme

11. b <- theta + 0.5        # borne supérieure de l'uniforme

13. ## Initialisation de la matrice dans laquelle seront

14. ## stockées les valeurs des estimateurs. On donne également

15. ## des noms aux lignes de la matrice afin de facilement

16. ## identifier les estimateurs.

17. x <- matrix(0, nrow = 3, ncol = nsimul)

18. rownames(x) <- c("Moyenne", "Mediane", "Mi-etendue")

20. ## Simulation comme telle.

21. for (i in 1:nsimul)

22. {

23. u <- runif(size, a, b)

24. x[, i] <- c(mean(u),  # moyenne

25. median(u) # médiane

26. mean(range(u))) # mi-étendue

27. }

29. ## On peut maintenant calculer la moyenne et la variance

30. ## par ligne.

31. rowMeans(x) - theta     # vérification du biais

32. apply(x, 1, var)        # comparaison des variances

Fig. C.2 - Définition de la fonction simul1

1. simul1 <- function(nsimul, size, theta)

2. {

3. a <- theta - 0.5

4. b <- theta + 0.5

5. x <- matrix(0, nrow = 3, ncol = nsimul)

6. rownames(x) <- c("Moyenne", "Mediane", "Mi-etendue")

7. for (i in 1:nsimul)

8. {

9. u <- runif(size, a, b)

10. x[, i] <- c(mean(u), median(u), mean(range(u)))

11. }

12. list(biais = rowMeans(x) - theta,

13. variances = apply(x, 1, var))

14. }

1. lapply(x, FUN, ...)

2. sapply(x, FUN, ...)

Ces fonctions appliquent la fonction FUN à tous les éléments de la liste ou du vecteur x et retournent les résultats sous forme de liste (lapply) ou, lorsque c'est possible, de vecteur ou de matrice (sapply). Il est important de noter que les valeurs successives de x seront passées comme premier argument à la fonction FUN. Le cas échéant, les autres arguments de FUN sont spécifiés dans le champ ....

Pour pouvoir utiliser ces fonctions dans le cadre d'une simulation comme celle dont il est question ici, il s'agit de définir une fonction qui fera tous les calculs pour une simulation, puis de la passer à sapply pour obtenir les résultats de simulations. La figure C.3 présente une première version d'une telle fonction. On remarquera que l'argument i ne joue aucun rôle dans la fonction. Voici un exemple d'utilisation pour un petit nombre de simulations (4) :

Fig. C.3 - Définitions des fonctions fun1 et simul2

1. fun1 <- function(i, size, a, b)

2. {

3. u <- runif(size, a, b)

4. c(Moyenne = mean(u),

5. Mediane = median(u),

6. "Mi-etendue" = mean(range(u)))

7. }

8. simul2 <- function(nsimul, size, theta)

9. {

10. a <- theta - 0.5

11. b <- theta + 0.5

13. x <- sapply(1:nsimul, fun1, size, a, b)

15. list(biais = rowMeans(x) - theta,

16. variances = apply(x, 1, var))

17. }

1. > sapply(1:4, fun1, size = 10, a = -0.5, b = 0.5)

3. [,1]        [,2]         [,3]        [,4]

4. Moyenne -0.076212566 -0.07839662 0.0006719503  0.02844156

5. Mediane -0.117479511 -0.01346815 0.0403468686  0.02980248

6. Mi-etendue -0.005960669 -0.05199959 0.0395837893 -0.03937694

On remarque donc que les résultats de chaque simulation se trouvent dans les colonnes de la matrice obtenue avec sapply.

Pour compléter l'analyse, on englobe le tout dans une fonction simul2, dont le code se trouve à la figure C.3 :

1. > simul2(10000, 100, 0)

3. $biais

4. Moyenne      Mediane   Mi-etendue

5. 0.0003148615 0.0002347348 0.0001190661

7. $variances

8. Moyenne      Mediane   Mi-etendue

9. 8.272098e-04 2.401754e-03 4.875476e-05

Il est généralement plus facile de déboguer le code avec cette approche puisque l'on peut rapidement circonscrire un éventuel problème à fun1 ou simul2.

C-4 - Variante de la seconde approche

Une chose manque d'élégance dans la seconde approche : l'obligation d'inclure un argument factice dans la fonction fun1. La fonction replicate (section 6.5) permet toutefois de passer outre cette contrainte. En effet, cette fonction exécute un nombre donné de fois une expression quelconque.

Fig. C.4 - Définitions des fonctions fun2 et simul3

1. fun2 <- function(size, a, b)

2. {

3. u <- runif(size, a, b)

4. c(Moyenne = mean(u),

5. Mediane = median(u),

6. "Mi-etendue" = mean(range(u)))

7. }

9. simul3 <- function(nsimul, size, theta)

10. {

11. a <- theta - 0.5

12. b <- theta + 0.5

14. x <- replicate(nsimul, fun2(size, a, b))

16. list(biais = rowMeans(x) - theta,

17. variances = apply(x, 1, var))

18. }

Les fonctions fun2 et simul3 de la figure C.4 sont des versions légèrement modifiées de fun1 et simul2 pour utilisation avec replicate. On a alors

1. Moyenne      Mediane   Mi-etendue

2. 8.578549e-04 2.483565e-03 4.928925e-05

C-5 - Gestion des fichiers

Pour un petit projet comme celui utilisé en exemple ici, il est simple et pratique de placer tout le code informatique dans un seul fichier de script. Pour un plus gros projet, cependant, il vaut souvent mieux avoir recours à plusieurs fichiers différents. Le présent auteur utilise pour sa part un fichier par fonction.

À des fins d'illustration, supposons que l'on utilise l'approche de la section C.4 avec la fonction replicate et que le code des fonctions fun2 et simul3 est sauvegardé dans des fichiers fun2.R et simul3.R, dans l'ordre. Si l'on crée un autre fichier, disons go.R, ne contenant que des expressions source pour lire les autres fichiers, il est alors possible de démarrer des simulations en exécutant ce seul fichier. Dans notre exemple, le fichier go.R contiendrait les lignes suivantes :

source("fun2.R")
source("simul3.R")
simul3(10000, 100, 0)

Une simple commande

> source("go.R")

exécutera alors une simulation complète.

C-6 - Exécution en lot

Les utilisateurs plus avancés pourront vouloir exécuter leur simulation R en lot (batch) pour en accélérer le traitement. Dans ce mode, aucune interface graphique n'est démarrée et tous les résultats sont redirigés vers un fichier pour consultation ultérieure. Pour les simulations demandant un long temps de calcul, c'est très pratique.

On exécute R en lot depuis la ligne de commande (Invite de commande sous Windows, Terminal sous OS X ou Linux). Une fois placé dans le répertoire contenant les fichiers de script, il suffit d'entrer à la ligne de commande

R CMD BATCH go.R

La sortie de cette commande (et donc tous les résultats des expressions R du fichier go.R) sera placée par défaut dans le fichier go.Rout. Sous Windows, le dossier d'installation de R peut ne pas se trouver dans la variable d'environnement %PATH%, auquel cas il faut spécifier le chemin d'accès complet de l'exécutable à la ligne de commande :

"c:\Program Files\R\R-x.y.z\bin\R" CMD BATCH go.R

Remplacer R-x.y.z par le numéro de version courant de R.

C-7 - Conclusion

Le nombre de simulations, et la taille de l'échantillon, , ont tous deux un impact sur la qualité des résultats, mais de manière différente. Quand augmente, la précision des estimateurs augmente. Ainsi, dans l'exemple ci-dessus, le biais et la variance des estimateurs de seront plus faibles. D'autre part, l'augmentation du nombre de simulations diminue l'impact des échantillons aléatoires individuels et, de ce fait, améliore la fiabilité des conclusions de l'étude.

D'ailleurs, les conclusions de l'étude de simulation sur le biais et la variance des trois estimateurs de la moyenne d'une loi uniforme sont les suivantes : les trois estimateurs sont sans biais et la mi-étendue a la plus faible variance. En effet, on peut démontrer mathématiquement que, pour impair,

et donc

pour tout .

D - Installation de packages dans R

Un package R est un ensemble cohérent de fonctions, de jeux de données et de documentation permettant de compléter les fonctionnalités du système de base ou d'en ajouter de nouvelles. Les packages sont normalement installés depuis le site Comprehensive R Archive Network (CRAN ; http://cran.r-project.org).

Cette annexe explique comment configurer R pour faciliter l'installation et l'administration de packages externes.

Les instructions ci-dessous se basent sur la création d'une bibliothèque personnelle où seront installés les packages R téléchargés de CRAN. Il est fortement recommandé de créer une telle bibliothèque. Cela permet d'éviter les éventuels problèmes d'accès en écriture dans la bibliothèque principale et de conserver les packages intacts lors des mises à jour de R. On montre également comment spécifier le site miroir de CRAN pour éviter d'avoir à le répéter à chaque installation de package.

1. Identifier le dossier de départ de l'utilisateur. En cas d'incertitude, examiner la valeur de la variable d'environnement HOME(7), depuis R avec la commande

> Sys.getenv("HOME")

ou, pour les utilisateurs de Emacs, directement depuis l'éditeur avec

M-x getenv RET HOME RET

Nous référerons à ce dossier par le symbole ~.

2. Créer un dossier qui servira de bibliothèque de packages personnelle. Dans la suite, nous utiliserons ~/R/library.

3. La configuration de R se fait à l'aide de simples fichiers texte, comme pour GNU Emacs ; voir la section A.7. Dans un fichier nommé ~/.Renviron (donc situé dans le dossier de départ), enregistrer la ligne suivante :

R_LIBS_USER="~/R/library"

Au besoin, remplacer le chemin ~/R/library par celui du dossier créé à l'étape précédente. Utiliser la barre oblique avant (/) dans le chemin pour séparer les dossiers.

Sous OS X, ajouter dans le fichier ~/.Renviron la ligne suivante :

R_INTERACTIVE_DEVICE=quartz

Ainsi, R utilisera toujours l'interface Quartz native de OS X pour afficher les graphiques.

4. Dans un fichier nommé ~/.Rprofile, enregistrer l'option suivante :

options(repos = "http://cran.ca.r-project.org")

Si désiré, remplacer la valeur de l'option repos par l'URL d'un autre site miroir de CRAN.

Les utilisateurs de GNU Emacs voudront ajouter une option pour éviter que R ait recours aux menus graphiques Tcl/Tk. Le code à entrer dans le fichier ~/.Rprofile sera plutôt

options(repos = "http://cran.ca.r-project.org", menu.graphics = FALSE)

Consulter la rubrique d'aide de Startup pour les détails sur la syntaxe et l'emplacement des fichiers de configuration, celles de library et .libPaths pour la gestion des bibliothèques et celle de options pour les différentes options reconnues par R.

Après un redémarrage de R, la bibliothèque personnelle aura préséance sur la bibliothèque principale et il ne sera plus nécessaire de préciser le site miroir de CRAN lors de l'installation de packages. Ainsi, la simple commande

> install.packages("actuar")

téléchargera le package de fonctions actuarielles actuar depuis le miroir canadien de CRAN et l'installera dans le dossier ~/R/library. Pour charger le package en mémoire, on fera

> library("actuar")

On peut arriver au même résultat sans utiliser les fichiers de configuration .Renviron et .Rprofile. Il faut cependant recourir aux arguments lib et repos de la fonction install.packages et à l'argument lib.loc de la fonction library. Consulter les rubriques d'aide de ces deux fonctions pour de plus amples informations.

Réponses des exercices

Chapitre 2

Solution de l'exercice 2.1

1. Il y a plusieurs façons de créer les troisième et quatrième éléments de la liste. Le plus simple consiste à utiliser numeric() et logical() :

> x <- list(1:5, data = matrix(1:6, 2, 3), numeric(3),
+ test = logical(4))

2. > names(x)

3. >

mode(x$test)
> length(x$test)

4. > dim(x$data)

5. > x[[2]][c(2, 3)]

6. > x[[3]] <- 3:8

Solution de l'exercice 2.2

1. > x[2]

2. > x[1:5]

3. > x[x > 14]

4. > x[-c(6, 10, 12)]

Solution de l'exercice 2.3

1. > x[4, 3]

2. > x[6, ]

3. > x[, c(1, 4)]

4. > x[x[, 1] > 50, ]

Chapitre 3

Solution de l'exercice 3.1

1. > rep(c(0, 6), 3)

2. > seq(1, 10, by = 3)

3. > rep(1:3, 4)

4. > rep(1:3, 1:3)

5. > rep(1:3, 3:1)

6. > seq(1, 10, length = 3)

7. > rep(1:3, rep(4,3))

Solution de l'exercice 3.2

1. > 11:20 / 10

2. > 2 * 0:9 + 1

3. > rep(-2:2, 2)

4. > rep(-2:2, each = 2)

5. > 10 * 1:10

Solution de l'exercice 3.3 Soit mat une matrice.

1. > apply(mat, 1, sum)

2. > apply(mat, 2, sum)

3. > apply(mat, 1, mean)

4. > apply(mat, 2, mean)

Solution de l'exercice 3.4> cumprod(1:10)

Solution de l'exercice 3.5x == (x %% y) + y * (x %/% y)

Solution de l'exercice 3.6

> x[1:5]
> head(x, 5)

5. > max(x)

> mean(x[1:5])
> mean(head(x, 5))

6. S

> mean(x[16:20])
> mean(x[(length(x) - 4):length(x)])  # plus général
> mean(tail(x, 5))                    # plus lisible!

Solution de l'exercice 3.7

1. (j - 1)*I + i

2. ((k - 1)*J + j - 1)*I + i

Solution de l'exercice 3.8

1. > rowSums(mat)

2. > colMeans(mat)

3. > max(mat[1:3, 1:3])

4. > mat[rowMeans(mat) > 7,]

Solution de l'exercice 3.9> temps[match(unique(cummin(tps)), temps)]

Chapitre 4

Solution de l'exercice 4.1> sum(P / cumprod(1 + i))

Solution de l'exercice 4.2

> x <- c(7, 13, 3, 8, 12, 12, 20, 11)
> w <- c(0.15, 0.04, 0.05, 0.06, 0.17, 0.16, 0.11, 0.09)
> sum(x * w)/sum(w)

Solution de l'exercice 4.3> 1/mean(1/x)

Solution de l'exercice 4.4

> lambda <- 2
> x <- 5
> exp(-lambda) * sum(lambda^(0:x)/gamma(1 + 0:x))

Solution de l'exercice 4.5

> x <- 10^(0:6)
> probs <- (1:7)/28

1. > sum(x^2 * probs) - (sum(x * probs))^2

Solution de l'exercice 4.6

> i <- 0.06
> 4 * ((1 + i)^0.25 - 1)

Solution de l'exercice 4.7

> n <- 1:10
> i <- seq(0.05, 0.1, by = 0.01)
> (1 - outer((1 + i), -n, "^"))/i

ou

> n <- 1:10
> i <- (5:10)/100
> apply(outer(1/(1+i), n, "^"), 1, cumsum)

Solution de l'exercice 4.8

> v <- 1/1.06
> k <- 1:10
> sum(k * v^(k - 1))

Solution de l'exercice 4.9

> pmts <- rep(1:4, 1:4)
> v <- 1/1.07
> k <- 1:10
> sum(pmts * v^k)

Solution de l'exercice 4.10

> v <- cumprod(1/(1 + rep(c(0.05, 0.08), 5)))
> pmts <- rep(1:4, 1:4)
> sum(pmts * v)

Chapitre 5

Solution de l'exercice 5.1

1. variance <- function(x, biased = FALSE)

2. {

3. if (biased)

4. {

5. n <- length(x)

6. (n - 1)/n * var(x)

7. }

8. else

9. var(x)

10. }

Solution de l'exercice 5.2 Une première solution utilise la transposée. La première expression de la fonction s'assure que la longueur de data est compatible avec le nombre de lignes et de colonnes de la matrice demandée.

1. matrix2 <- function(data = NA, nrow = 1, ncol = 1,

2. bycol = FALSE, dimnames = NULL)

3. {

4. data <- rep(data, length = nrow * ncol)

6. if (bycol)

7. dim(data) <- c(nrow, ncol)

8. else

9. {

10. dim(data) <- c(ncol, nrow)

11. data <- t(data)

12. }

14. dimnames(data) <- dimnames

15. data

16. }

La seconde solution n'a pas recours à la transposée. Pour remplir la matrice par ligne, il suffit de réordonner les éléments du vecteur data en utilisant la formule obtenue à l'exercice 3.7.

1. matrix2 <- function(data = NA, nrow = 1, ncol = 1,

2. bycol = FALSE, dimnames = NULL)

3. {

4. data <- rep(data, length = nrow * ncol)

6. if (!bycol)

7. {

8. i <- 1:nrow

9. j <- rep(1:ncol, each = nrow)

10. data <- data[(i - 1)*ncol + j]

11. }

12. dim(data) <- c(nrow, ncol)

13. dimnames(data) <- dimnames

14. data

15. }

Solution de l'exercice 5.3

1. phi <- function(x)

2. {

3. exp(-x^2/2) / sqrt(2 * pi)

4. }

Solution de l'exercice 5.4

1. Phi <- function(x)

2. {

3. n <- 1 + 2 * 0:50

4. 0.5 + phi(x) * sum(x^n / cumprod(n))

5. }

Solution de l'exercice 5.5 La première solution utilise une fonction interne et une structure de contrôle if ... else.

1. Phi <- function(x)

2. {

3. fun <- function(x)

4. {

5. n <- 1 + 2 * 0:50

6. 0.5 + phi(x) * sum(x^n / cumprod(n))

7. }

9. if (x < 0)

10. 1 - fun(-x)

11. else

12. fun(x)

13. }

Seconde solution sans structure de contrôle if ... else. Rappelons que dans des calculs algébriques, FALSE vaut 0 et TRUE vaut 1.

1. Phi <- function(x)

2. {

3. n <- 1 + 2 * 0:50

4. neg <- x < 0

5. x <- abs(x)

6. neg + (-1)^neg * (0.5 + phi(x) * sum(x^n / cumprod(n)))

7. }

Solution de l'exercice 5.6

1. Phi <- function(x)

2. {

3. n <- 1 + 2 * 0:30

4. 0.5 + phi(x) * colSums(t(outer(x, n, "^")) / cumprod(n))

5. }

Solution de l'exercice 5.7

1. prod.mat <- function(mat1, mat2)

2. {

3. if (ncol(mat1) == nrow(mat2))

4. {

5. res <- matrix(0, nrow = nrow(mat1),

6. ncol = ncol(mat2))

7. for (i in 1:nrow(mat1))

8. {

9. for (j in 1:ncol(mat2))

10. {

11. res[i, j] <- sum(mat1[i,] * mat2[,j])

12. }

13. }

14. res

15. }

16. else

17. stop("Les dimensions des matrices ne permettent pas le produit matriciel.")

18. }

1. prod.mat<-function(mat1, mat2)

2. {

3. if (ncol(mat1) == nrow(mat2))

4. {

5. res <- matrix(0, nrow = nrow(mat1),

6. ncol = ncol(mat2))

7. for (i in 1:nrow(mat1))

8. res[i,] <- colSums(mat1[i,] * mat2)

9. res

10. }

11. else

12. stop("Les dimensions des matrices ne permettent pas le produit matriciel.")

13. }

Solutions bonus : deux façons de faire équivalentes qui cachent la boucle dans un sapply.

1. prod.mat<-function(mat1, mat2)

2. {

3. if (ncol(mat1) == nrow(mat2))

4. t(sapply(1:nrow(mat1),

5. function(i) colSums(mat1[i,] * mat2)))

6. else

7. stop("Les dimensions des matrices ne permettent pas le produit matriciel.")

8. }

10. prod.mat<-function(mat1, mat2)

11. {

12. if (ncol(mat1) == nrow(mat2))

13. sapply(1:ncol(mat2),

14. function(j) colSums(t(mat1) * mat2[,j]))

15. else

16. stop("Les dimensions des matrices ne permettent pas le produit matriciel.")

17. }

Solution de l'exercice 5.8 Le calcul à faire n'est qu'un simple produit matriciel, donc :

notes.finales <- function(notes, p) notes %*% p

Solution de l'exercice 5.10

1. param <- function (moyenne, variance, loi)

2. {

3. loi <- tolower(loi)

4. if (loi == "normale")

5. {

6. param1 <- moyenne

7. param2 <- sqrt(variance)

8. return(list(mean = param1, sd = param2))

9. }

10. if (loi == "gamma")

11. {

12. param2 <- moyenne/variance

13. param1 <- moyenne * param2

14. return(list(shape = param1, scale = param2))

15. }

16. if (loi == "pareto")

17. {

18. cte <- variance/moyenne^2

19. param1 <- 2 * cte/(cte-1)

20. param2 <- moyenne * (param1 - 1)

21. return(list(alpha = param1, lambda = param2))

22. }

23. stop("La loi doit etre une de \"normale\", \"gamma\" ou \"pareto\"")

24. }

Chapitre 6

Solution de l'exercice 6.1 Soit Xij et wij des matrices, et Xijk et wijk des tableaux à trois dimensions.

1. > rowSums(Xij * wij)/rowSums(wij)

2. > colSums(Xij * wij)/colSums(wij)

3. > sum(Xij * wij)/sum(wij)

4. > apply(Xijk * wijk, c(1, 2), sum) / apply(wijk, c(1, 2), sum)

5. > apply(Xijk * wijk, 1, sum)/apply(wijk, 1, sum)

6. > apply(Xijk * wijk, 2, sum)/apply(wijk, 2, sum)

7. > sum(Xijk * wijk)/sum(wijk)

Solution de l'exercice 6.2

1. > unlist(lapply(0:10, seq, from = 0))

2. > unlist(lapply(1:10, seq, from = 10))

3. > unlist(lapply(10:1, seq, to = 1))

Solution de l'exercice 6.3

1. > x <- lapply(seq(100, 300, by = 50), rpareto, shape = 2, scale = 5000)

2. > names(x) <- paste("sample", 1:5, sep = "")

3. > sapply(x, mean)

> lapply(x, function(x) sort(ppareto(x, 2, 5000)))
> lapply(lapply(x, sort), ppareto,  shape = 2, scale = 5000)

4. > hist(x$sample5)

5. > lapply(x, "+", 1000)

Solution de l'exercice 6.4

1. > mean(sapply(x, function(liste) liste$franchise))

Les crochets utilisés pour l'indiçage constituent en fait un opérateur dont le « nom » est [[. On peut donc utiliser cet opérateur dans la fonction sapply :

> mean(sapply(x, "[[", "franchise"))

2. > sapply(x, function(x) mean(x$nb.acc))

3. > sum(sapply(x, function(x) sum(x$nb.acc))) ou > sum(unlist(sapply(x, "[[", "nb.acc")))

4. > mean(unlist(lapply(x, "[[", "montants")))

5. > sum(sapply(x, function(x) sum(x$nb.acc) == 0))

6. > sum(sapply(x, function(x) x$nb.acc[1] == 1))

7. > var(unlist(lapply(x, function(x) sum(x$nb.acc))))

8. > sapply(x, function(x) var(x$nb.acc))

> y <- unlist(lapply(x, "[[", "montants"))
> sum(y <= x)/length(y)

La fonction ecdf retourne une fonction permettant de calculer la fonction de répartition empirique en tout point :

> ecdf(unlist(lapply(x, "[[", "montants")))(x)

> y <- unlist(lapply(x, "[[", "montants"))
> colSums(outer(y, x, "<="))/length(y)

La fonction retournée par ecdf accepte un vecteur de points en argument :

> ecdf(unlist(lapply(x, "[[", "montants")))(x)

Chapitre 7

Solution de l'exercice 7.1

> f <- function(x) x^3 - 2 * x^2 - 5
> uniroot(f, lower = 1, upper = 4)

1. Comme un simple graphique le démontre, il y a deux racines dans l'intervalle.

> f <- function(x) x^3 + 3 * x^2 - 1
> curve(f, xlim = c(-4, 0))
> uniroot(f, lower = -4, upper = -1)
> uniroot(f, lower = -1, upper = 0)

> f <- function(x) x - 2^(-x)
> uniroot(f, lower = 0, upper = 1)

> f <- function(x) exp(x) + 2^(-x) + 2 * cos(x) - 6
> uniroot(f, lower = 1, upper = 2)

> f <- function(x) exp(x) - x^2 + 3 * x - 2
> uniroot(f, lower = 0, upper = 1)

Solution de l'exercice 7.2

1. > X <- c(2061, 1511, 1806, 1353, 1600)

2. > w <- c(100155, 19895, 13735, 4152, 36110)

3. > g <- function(a, X, w, s2)

4. + {

5. +     z <- 1/(1 + s2/(a * w))

6. +     Xz <- sum(z * X)/sum(z)

7. +     sum(z * (X - Xz)^2)/(length(X) - 1)

8. + }

9. > uniroot(function(x) g(x, X, w, 140E6) - x, c(50000, 80000))

Solution de l'exercice 7.3

1. > dpareto <- function(x, alpha, lambda)

2. + {

3. +     (alpha * lambda^alpha)/(x + lambda)^(alpha+1)

4. + }

5. > f <- function(par, x) -sum(log(dpareto(x, par[1], par[2])))

6. > optim(c(1, 1000), f, x = x)

ou, en utilisant le truc du logarithme des paramètres expliqué dans le code informatique de la section 7.4 pour éviter les soucis de convergence :

1. > dpareto <- function(x, logAlpha, logLambda)

2. + {

3. +     alpha <- exp(logAlpha)

4. +     lambda <- exp(logLambda)

5. +     (alpha * lambda^alpha)/(x + lambda)^(alpha+1)

6. + }

7. > optim(c(log(2), log(1000)), f, x = x)

8. > exp(optim(c(log(2), log(1000)), f, x = x)$par)

Chapitre 8

Solution de l'exercice 8.1

1. > x <- rlnorm(1000, meanlog = log(5000) - 0.5, sdlog = 1)

2. > hist(x)

Solution de l'exercice 8.2

1. > x <- rpois(10000, lambda = rgamma(10000, shape = 5, rate = 4))

2. > xx <- seq(min(x), max(x))

3. > px <- table(x)

4. > plot(xx, dnbinom(xx, size = 5, prob = 0.8),

5. + type = "h", lwd = 5, col = "blue")

6. > points(xx, px/length(x), pch = 16)

Solution de l'exercice 8.3

1. > w <- rbinom(1, 10000, 0.55)

2. > x <- c(rlnorm(w, 3.5, 0.6), rlnorm(10000 - w, 4.6, 0.3))

3. > hist(x)

Code source

Code source

Code informatique des sections d'exemples

Télécharger

Code source du document

Accéder au dépôt

Bibliographie

Abelson, H., G. J. Sussman et J. Sussman. 1996, Structure and Interpretation of Computer Programs, 2e éd., MIT Press, ISBN 0-26201153-0.

Becker, R. A. 1994, «A brief history of S», cahier de recherche, AT&T Bell Laboratories. URL http://cm.bell-labs.com/cm/ms/departments/sia/doc/94.11.ps.

Becker, R. A. et J. M. Chambers. 1984, S: An Interactive Environment for Data Analysis and Graphics, Wadsworth, ISBN 0-53403313-X.

Becker, R. A., J. M. Chambers et A. R. Wilks. 1988, The New S Language: A Programming Environment for Data Analysis and Graphics, Wadsworth & Brooks/Cole, ISBN 0-53409192-X.

Braun, W. J. et D. J. Murdoch. 2007, A First Course in Statistical Programming with R, Cambridge University Press, ISBN 978-0-52169424-7.

Cameron, D., J. Elliott, M. Loy, E. S. Raymond et B. Rosenblatt. 2004, Leaning GNU Emacs, 3e éd., O'Reilly, Sebastopol, CA, ISBN 0-59600648-9.

Chambers, J. M. 1998, Programming with Data: A Guide to the S Language, Springer, ISBN 0-38798503-4.

Chambers, J. M. 2000, «Stages in the evolution of S», URL http://cm.bell-labs.com/cm/ms/departments/sia/S/history.html.

Chambers, J. M. 2008, Software for Data Analysis: Programming with R, Springer, ISBN 978-0-38775935-7.

Chambers, J. M. et T. J. Hastie. 1992, Statistical Models in S, Wadsworth & Brooks/Cole, ISBN 0-53416765-9.

Dutang, C., V. Goulet et M. Pigeon. 2008, «actuar: An R package for actuarial science», Journal of Statistical Software, vol. 25, no 7. URL http://www.jstatsoft.org/v25/i07.

Hornik, K. 2013, « The R FAQ », URL http://cran.r-project.org/doc/FAQ/R-FAQ.html.

Iacus, S. M., S. Urbanek, R. J. Goedman et B. D. Ripley. 2013, « R for Mac OS X FAQ », URL http://cran.r-project.org/bin/macosx/RMacOSX-FAQ.html.

IEEE. 2003, 754-1985 IEEE Standard for Binary Floating-Point Arithmetic, IEEE, Piscataway, NJ.

Ihaka, R. et R. Gentleman. 1996, «R: A language for data analysis and graphics», Journal of Computational and Graphical Statistics, vol. 5, no 3, p. 299-314.

Ligges, U. 2003, « R-winedt », dans Proceedings of the 3rd International Workshop on Distributed Statistical Computing (DSC 2003), édité par K. Hornik, F. Leisch et A. Zeileis, TU Wien, Vienna, Austria, ISSN 1609-395X. URL http://www.ci.tuwien.ac.at/Conferences/DSC-2003/Proceedings/.

Redd, A. 2010, «Introducing NppToR: R interaction for Notepad++», R Journal, vol. 2, no 1, p. 62-63. URL http://journal.r-project.org/archive/2010-1/RJournal_2010-1.pdf.

Ripley, B. D. et D. J. Murdoch. 2013, « R for Windows FAQ », URL http://cran.r-project.org/bin/windows/base/rw-FAQ.html.

Venables, W. N. et B. D. Ripley. 2000, S Programming, Springer, New York, ISBN 0-38798966-8.

Venables, W. N. et B. D. Ripley. 2002, Modern Applied Statistics with S, 4e éd., Springer, New York, ISBN 0-38795457-0.

Venables, W. N., D. M. Smith et R Core Team. 2013, An Introduction to R, R Foundation for Statistical Computing. URL http://cran.r-project.org/doc/manuals/R-intro.html.

Wheeler, B. 2013, SuppDists : Supplementary Distributions. URL http://cran.r-project.org/package=SuppDists, R package version 1.1-9.

Note de la rédaction Developpez.com

Tous nos remerciements à Vincent Goulet pour l'autorisation de publication de ce cours. Nos remerciements également à Winjerome pour la mise au gabarit Developpez.com, et à Claude Leloup pour la relecture orthographique.

Index

Les éléments avec une étoile indiquent les pages où les concepts sont introduits, définis ou expliqués.

- Symboles -

distribution

liste (1)

racine

! (*1)

Bêta (1)

lm (1)

d'un polynôme (1)

!= (*1)

binomiale (1)

log (1)

d'une fonction (1)

$ (*1), (2), (3)

Binomiale (1)

logical (*1), (2)

rang (1)

$<- (1)

Binomiale négative (1)

logical (mode) (*1), (2), (3)

range (*1), (1)

%% (*1)

Cauchy (1)

longueur (*1), (2)

rank (*1), (1)

%*% (*1), (2)

Exponentielle (1)

lower (1)

rbind (*1), (2), (3)

%/% (*1)

F (1)

ls (1), (2), (3)

rbinom (1)

%in% (*1), (2)

gamma (1), (2), (3)

- M -

Re (1)

%o% (*1), (2)

Géométrique (1)

mapply (1), (*2), (3)

renverser un vecteur (1)

& (*1)

Hypergéométrique (1)

match (1), (2)

rep (1), (*2), (3), (4), (5)

&& (**1)

Khi carré (1)

matrice (1), (2), (3), (4)

repeat (*1), (2), (3), (4), (5), (6)

* (*1)

log-normale (1), (1)

diagonale (1)

répertoire de travail (1)

+ (*1)

logistique (1)

identité (1)

répétition de valeurs (1)

- (*1)

mélange discret (1)

inverse (1)

replace (1), (2), (3), (4)

-> (1), (2)

mélange Poisson/gamma (1)

moyennes par colonne (1)

replicate (*1), (2), (3)

->> (*1)

normale (1)-(2), (3)

moyennes par ligne (1)

return (*1)

-Inf (*1)

Pareto (1), (2)

sommes par colonne (1)

rev (*1), (2), (3)

... (*1), (2)

Poisson (1), (2), (3)

sommes par ligne (1)

rgamma (1), (2)

/ (*1)

t (1)

transposée (1)

rm (1), (2)

: (1), (2)

Uniforme (1)

matrix (1), (*2), (3), (4), (5), (6)

RNGkind (1)

; (1)

Weibull (1)

matrix (classe) (2)

rnorm (1), (2), (3)

< (*1)

Wilcoxon (1)

max (1), (*2), (3), (4)

round (1), (*2), (3)

<- (*1), (2)

dnorm (1)

maximum

row.names (1)

<<- (*1), (2)

dossier de travail, voir répertoire de travail

cumulatif (1)

rowMeans (*1), (2), (3)

<= (*1)

dpois (1)

d'un vecteur (1)

rowSums (*1), (2), (3), (4), (5)

= (1)

- E -

local (1)

rpois (1)

== (*1)

écart type (1)

parallèle (1)

RStudio (1), (2)-(3)

> (*1)

ecdf (1)

position dans un vecteur (1)

configuration (1)

>= (*1)

else (*1), (2), (3), (4)

mean (1), (2), (*3), (4), (5), (6)

projets (1)

[ (*1)

Emacs (1), (2), (3), (4)

median (*1), (2)

sous-fenêtres (1)

[ ] (1), (2), (3), (*4)

déplacement du curseur (1)

médiane (1)

symbole d'assignation (1), (2)

[<- (*1)

mode ESS (1), (2), (3)

methods (*1)

runif (1), (2)

[[ (1)

nouveau fichier (1)

min (1), (*2), (3)

- S -

[[ ]] (*1), (2)

rechercher et remplacer (1)

minimum

S (1), (2)

^ (*1)

sauvegarder (1)

cumulatif (1)

S+ (1)

{ } (*1)

sauvegarder sous (1)

d'un vecteur (1)

S-PLUS (1)

| (*1)

sélection (1)

fonction non linéaire (1)

sample (1), (2), (3), (4), (5), (*6), (7)

|| (*1)

étiquette (1), (2)

local (1)

sapply (1), (2), (*3), (4)-(5), (6)-(7), (8), (9), (10)

- A -

eval (1), (2)

parallèle (1)

save.image (1), (2), (3), (4)

abs (1), (2)

exists (1)

position dans un vecteur (1)

Scheme (1)

add (1)

exp (1), (2)-(3), (4)

mode (*1), (*2), (3)

sd (*1), (2), (3)

affectation (1)

expression (1), (2), (3), (4)

(mode) (*1), (2)

search (*1), (2)

apply (1), (2), (3), (*4), (5), (6), (7), (8)

expression (mode) (*1)

moyenne

seq (1), (2), (*3), (4), (5), (6)

array (*1), (2), (3)

extraction, voir aussi indiçage

arithmétique (1)

seq_len (*1)

array (classe) (1)

derniers éléments (1)

harmonique (1)

simulation

arrondi (1)

éléments différents (1)

pondérée (1), (2)

nombres uniformes (1)

as.data.frame (*1)

premiers éléments (1)

tronquée (1)

planification (1)-(2)

attach (*1), (2)

- F -

- N -

variables aléatoires (1)

attr (*1), (2), (3)

F, voir FALSE

NA (*1), (2)

sin (1)

attribut (*1)

factorial (1), (*2)

na.rm (*1), (2), (3)

solve (1), (*2), (3)

attributes (*1), (2)

FALSE (1), (2)

names (1)

somme (1)

- B -

floor (*1), (2)

names (attribut) (*1)

cummulative (1)

biscuits, voir Syndrôme de la plaque à biscuits

fonction

NaN (*1)

sort (*1), (2)

boucle (1), (2), (3)

anonyme (1)

nchar (*1), (2)

source (1)

break (*1), (2), (3), (4)

appel (1)

ncol (1), (2), (3), (*4), (5), (6)

stop (1)

by (1), (2)

débogage (1)

next (*1)

structure (1)

byrow (1), (2)

définie par l'usager (1)

nlm (*1), (2), (3)

style (1)

- C -

générique (1)

nlminb (*1)

suite des nombres (1), (2)

c (*1)

maximum local (1)

noms d'objets

sum (1), (*2), (3), (4), (5)

cat (1)

minimum (1)

conventions (1)

summary (*1), (2), (3), (4)

cbind (*1), (2), (3), (4)

minimum local (1)

réservés (1)

switch (*1)

ceiling (*1), (2)

optimisation (1)

Notepad++ (1)

Syndrôme de la plaque à biscuits (1)

character (*1), (2)

racine (1)

nrow (1), (2), (3), (*4), (5), (6)

system.time (1)

character (mode) (*1), (2)

résultat (1)

NULL (*1), (2)

- T -

choose (1)

for (*1), (2), (3), (4), (5)

NULL (mode) (1)

T, voir TRUE

class (1), (2), (3)

function (*1), (2)-(3), (4), (5)

numeric (*1), (2), (3), (4)

t (1), (*2), (3)

class (attribut) (*1)

function (mode) (*1)

numeric (mode) (*1), (2)

table (1)

colMeans (*1), (2), (3), (4)

- G -

- O -

tableau (1), (2)

colSums (*1), (2), (3), (4), (5)

gamma (1), (2), (*3)

optim (*1), (2)

tail (*1), (2)

compilé (langage) (1)

- H -

optimize (*1), (2)

tri (1)

complex (1)

head (*1), (2)

order (*1), (2)

TRUE (1), (2)

complex (mode) (*1)

hist (1), (2)

ordre (1)

trunc (*1), (2)

cos (1)

- I -

outer (*1), (2)-(3), (4), (5), (6), (7)

typeof (*1)

cummax (*1), (2)

if (*1), (2), (3), (4), (5)

- P -

- U -

cummin (*1), (2)

ifelse (*1), (2), (3)

package (*1)

unique (*1), (2)

cumprod (*1), (2), (3)

Im (1)

bibliothèque personnelle (1)

uniroot (1), (2)

cumsum (*1), (2)

indiçage

installation (1), (2)

unlist (1), (*2), (3), (4)

curve (1)

liste (*1), (2)

paste (1)

upper (1)

- D -

matrice (1), (*2), (3)

pgamma (1), (2)

- V -

data (1), (2), (3)

vecteur (*1), (2)

plot (1), (2), (3), (4), (5)

valeur actuelle (1), (2)-(3)

data.frame (*1)

Inf (*1)

pmax (*1), (2)

var (*1), (2), (3)

data.frame (classe) (1)

install.packages (*1)

pmin (*1), (2)

variable

dbeta (1)

interprété (langage) (1)

pnorm (1)

globale (1)

dbinom (1)

is.finite (*1)

point fixe (1), (2)

locale (1)

density (1)

is.infinite (*1)

points (1)

variance (1)

det (*1)

is.na (*1), (2), (3)

polyroot (*1), (2)

vecteur (1), (2)

detach (*1), (2)

is.nan (*1)

print (1), (2), (3), (4), (5), (6)

vector (1)

dgamma (1)

is.null (*1)

prod (*1), (2), (3)

vide, voir NULL

diag (1), (*2), (3)

- L -

produit (1)

- W -

diff (*1), (2)

lapply (1), (2), (*3), (4)-(5), (6), (7), (8)

cumulatif (1)

which (*1), (2)

différences (1)

length (1), (*2), (3), (4)

extérieur (1)

which.max (*1), (2)

dim (1), (2), (3)

lfactorial (1)

- Q -

which.min (*1), (2)

dim (attribut) (*1), (2)

lgamma (1)

q (1), (2)

while (*1), (2), (3)

dimension (1), (2)

library (*1), (2)

quantile (*1), (2)

WinEdt (1)

dimnames (1), (2)

list (*1), (2), (3)

Quartz (1)

dimnames (attribut) (*1)

list (mode) (*1), (2)

- R -

(1)http://www.gnu.org

(2)Merci à Markus Gesmann pour cette observation.

(3)http://www.r-project.org

(4)En fait, cela devient un appel à l'opérateur « ( » qui ne fait que retourner son argument.

(5)À ce compte, GNU Emacs constitue un environnement de développement intégré. Seulement, nous avons davantage insisté sur ses fonctionnalités d'éditeur de texte dans le présent document.

(6)Les sous-fenêtres sont appelées panes (en anglais) dans l'application.

(7)Pour les utilisateurs de GNU Emacs sous Windows, la variable est créée par l'assistant d'installation de Emacs lorsqu'elle n'existe pas déjà.