Modélisation stochastique par R

L3 M1A

Faculté des Sciences et Techniques (FAST)

Université d’Abomey-Calavi (UAC)
Qu’est ce que R ?

I langage de programmation et logiciel d’analyse statistique et

graphique ;
I logiciel gratuit disponible sur le site
http ://www.r-project.org
I avec
I les fichiers d’installation, de mise à jour, la documentation ;
I les librairies (ou packages) de base qui regroupent les fonctions
usuelles ;
I les librairies (ou packages) complémentaires.
Au démarrage de R

> est le prompt qui attend les commandes et qui apparaît en début
de chaque ligne de commandes.

# permet d’insérer un commentaire ;

+ apparaît au début d’une ligne si la ligne de commande

précédente est incomplète ;

; permet de séparer plusieurs commandes écrites sur la même ligne ;

q() permet d’arrêter le logiciel.

Au démarrage

Toutes les librairies ne sont pas chargées au démarrage ;

library() affiche la liste des librairies installées ;

library(lib) charge la librairie lib ; les fonctions de cette librairie

peuvent alors être utilisées ;

library(help = lib) afiche la liste des fonctions de la librairie lib.

Aide sur une fonction R

I help(FUN) ou ?FUN (sous windows) ? ?FUN (sous Linux)

affiche l’aide sur la fonction FUN

I example(FUN) exécute les exemples généralement présentés à

la fin du fichier d’aide de la fonction FUN

I demo() affiche la liste des démonstrations des fonctionnalités

de R.
Objets

Les objets de base sont

I vecteurs

I matrices

I tableaux de données (data.frame) : matrices dont toutes les

colonnes ne sont pas nécessairement du même type
Objets

Les objets sont caractérisés par :

I leur nom,
I leur contenu
I des attributs qui précisent le type de données représenté et les
traitements statistiques possibles.

ls() # affiche la liste des objets en mémoire.

rm(o) # supprime l’objet o.

rm(list=ls()) # supprime tous les objets en mémoire.

Attributs des objets

I Le mode : c’est le type des éléments contenus dans un objet ;

les quatre principaux sont :
I numérique
I caractère
I complexe
I logique (FALSE ou TRUE).

I la longueur : c’est le nombre d’éléments de l’objet.

mode et length affichent respectivement le mode et la longueur

d’un objet.
Opérateurs et connecteurs

+ Addition
- Soustraction
Opérateurs arithmétiques
* Produit
/ Division

> Supérieur
< Inférieur
Opérateurs logiques >= Supérieur ou égal
<= Inférieur ou égal
== Egal
!= Différent

& Et
Connecteurs
| Ou
Opérations de base sur les nombres

2+3 # addition
2-3 # soustraction
2*3 # produit
2/3 # division
2ˆ3 # puissance
sqrt(2) # racine carrée
sqrt(as.complex(-1)) # racine carrée
sqrt(-1+2i) # racine carrée
abs(-1+2i) # module
exp(1) # constante d’Euler e
exp(2) # e2
Opérations de base

log(5) # ln(5)
log10(5) # ln(5)/ln(10)
1.e4 # 104
1.e-4 # 10−4
Inf # +∞
-Inf # −∞
Na # valeur manquante
NaN # valeurs indéfinies ; ex −∞ + ∞
x=1 # affectation
x<-1 # affectation
Fonctions trigonométriques

sin() # sinus
cos() # cosinus
tan() # tangente
acos() # arccosinus
asin() # arcsinus
atan() # arctangente
cosh() # cosinus hyperbolique
sinh() # sinus hyperbolique
tanh() # tangente hyperbolique
acosh() # argument cos hyperbolique
asinh() # argument sin hyperbolique
atanh() # argument tan hyperbolique
Opération de base

Re(z) # partie réelle de z

Im(z) # partie imaginaire de z
Mod(z) # module de z
Arg(z) # argument de z
factorial(n) # n!
choose(n,k) # Cnk
Opérations sur les vecteurs

√
v=c(1,sqrt(2),1/3,4,-3.235,5) # Création du vecteur (1, 2,1/3,4,-3
v # Affichage du vecteur
v-1 # Soustraction
2*v # Multiplication par 2
v/2 # Division par 2
v[2] # Affichage de la 2ème coordonnée d
floor(v) # Partie entière des composantes de v
round(v,2) ) # Arrondi à la 2ème décimale des com
u=v[3 :5] # Création du vecteur u = (1, 4, 5)
w=v[c(1,3,6)] # Création du vecteur w = (1, 1, 5)
x=c(0,w,u) # Création du vecteur x = (0, 1, 1, 5,
Opérations sur les vecteurs

v [v < 4] # Affichage des coordonnées de v < 4

v [which(v >= 3)] # Affichage des coordonnées de v ≥3
u+w # Somme de u et w
u*w # Produit terme à terme
u/w # Division terme à terme
t(u) # Transposée de u
sum(u) # Somme des coordonnées de u
prod(u) # Produit des coordonnées de u
length(u) # Affichage de la dimension de u
unique(u) #Affichage des composantes qui ne se répètent p
duplicated(u) #Affichage des composantes qui se répètent dans
Opérations sur les vecteurs

t(u)%*%w # Produit scalaire de t(u) et w

sort(v) ou sort(v,F) # Tri croissant des coordonnées de v
sort(v,T) # Tri décroissant des coordonnées de v
min(v) # Affichage du plus petit élément de v
max(v) # Affichage du plus grand élément de v
which.min(v) # Affichage de la position du minimum de v
which.max(v) # Affichage de la position du maximum de v
cumsum(v) # Sommes cumulées des éléments de v
Générer une séquence de données

I rep(v,n) permet de répéter n fois les éléments d’un vecteur v ;

Exemple
rep(1, 10)
rep(c(2,3),4)
v=c(0,-1,2,1.25) ; rep(v,3)
I seq(a,b,by=r) génère des termes d’une suite arithmétique de
raison r, où a est le 1er terme et le dernier est supérieur ou
égal à b

Exemple
seq(10,1,-0.5)
seq(2,12,8)
Générer une séquence régulière de données

I a :b permet de générer des termes d’une suite arithmétique de

raison 1, où le 1er terme est a et le dernier est inférieur ou égal
à b;

Exemple
1 :10
0.5 :11
I b :a génère des termes d’une suite arithmétique de raison -1,
où le 1er terme est b et le dernier est supérieur ou égal à a.

Exemple
10 :0.5
8.3 :1.2
Matrices

Pour créer une matrice à n lignes et p colonnes :

matrix(data = x, nrow = n, ncol = p, byrow = FALSE)

où x est le vecteur contenant les éléments de la matrice qui seront

rangés par colonne.
Si byrow=TRUE dans la syntaxe, les éléments de la matrice seront
rangés par ligne.
Opérations sur les matrices

x=1 :20
A=matrix(x,nrow=5,byrow=T) # Création de la matrice 5 × 4
# dont la 1ère ligne est (1,2,3,4)
B=matrix(a,nrow=5) # Création de la matrice 5 × 4
# dont la 1ère colonne est (1,2,3,4,5)
t(A) #Transposition de la matrice A
A[1,2] # Sélection de l’élément à l’intersection
# de la 1ère ligne et de la 3ème colonne de A
A[,3] # Sélection de la 3ème colonne de A
A[c(3,4),] # Sélection des 3ème et 4ème lignes de A
A[,c(3,4)] # Sélection des 3ème et 4ème colonnes de A
A[-2,] # Suppression de la 2ème ligne de A
Opérations sur les matrices

A[,-c(1,4)] # Suppression des 1ère et 4ème colonnes de A

A+1 # Addition de 1 à tous les termes de A
A*2 # Multiplication par 2 à tous les termes de A
A+1 :5 # Addition du vecteur 1 :5 à toutes les colonnes de A
A+B # Somme de A et B
A*B # Produit terme à terme de A et B
A[A>9] # Affichage des éléments de A supérieurs à 9
diag(1 :4) # Création d’une matrice diagonale de diagonale (1,2,3,4
diag(A) # Affichage de la diagonale de la matrice A
Opérations sur les matrices

rbind(A,B) # Concaténation verticale des matrices A et B

cbind(x1,x4) # Concaténation horizontale des matrices A et B
colSums(A) # Calcul des sommes par colonne de A
rowSums(A) # Calcul des sommes par ligne de A
colMeans(A) # Calcul des moyennes par colonne de A
rowMeans(A) # Calcul des moyennes par ligne de A
C=A%*%t(B) # Produit matriciel de A par la transposée de B
det(C) # Calcul du déterminant de C
eigen(A) # Calcul des valeurs/vecteurs propres de A
solve(A) # Calcul de l’inverse de A
solve(A,b) # Calcul de la solution x de l’équation Ax=b
Statistique descriptive ; caractéristiques numériques

table(Nile) # Affichage du tableau des modalités de x av

range(Nile) # Affichage du maximum et du minimum de
mean(Nile) # Calcul de la moyenne de Nile
var(Nile) # Calcul de la variance empirique
sd(Nile) # Calcul de l’écart-type Nile
median(Nile) # Calcul de la médiane empirique de Nile
quantile(Nile) # Calcul des quantiles empiriques de Nile
summary(Nile) # Résumé statistique de Nile
x=c(0.5,1,2.2) ; y=x-1
cov(x,y) #covariance entre x et y
corr(x,y) #coefficient de corrélation entre x et y
plot(ecdf(x)) # tracé de la fonction de répartition de x
Statistique descriptive : représentations graphiques

plot(x) # Tracé du nuage de points

pie(x) # Diagramme en secteurs
barplot(x) # Diagramme en barres
hist(Nile,nclass=10) # Histogramme de Nile constitué de 10 classes
boxplot(Nile) # Boîte à moustâche
Lancer de dés et tirage dans une urne

I Un résultat du lancer d’un dé bien équilibré, n fois

sample(1 :6,n,replace=T)

I Une urne contenant 10 boules dont 3 noires, 5 blanches et

deux rouges

urne=c(rep("black",3),rep("white",5), rep("red",2))

I Un tirage de quatre boules sans remise dans cette urne

sample(urne,4,replace=F)
I Un tirage de quatre boules avec remise dans cette urne
sample(urne,4,replace=T)
Echantillonnage

I sample(v,n,rep=T) tire avec remise n éléments parmi les

composantes d’un vecteur v

sample(10 :15,8,rep=T)
I sample(v,n,rep=F) tire sans remise de n éléments parmi les
composantes d’un vecteur v

sample(1 :20,8,rep=F)
I sample(v) permute les composantes du vecteur v

sample(c(0.5,1,-0.2,3))
Echantillon d’une loi discrète finie

Pour générer un échantillon de taille n d’une v.a.r.d de loi de

probabilité {(xi , pi ), i = 1, · · · , k},

sample(x,n,replace=T,prob)

où x = (x1 , x2 , · · · , xk )
prob = (p1 , p2 , · · · , pk )

I Si replace=F, c’est un tirage sans remise et l’échantillon n’est

pas iid ;
I Si prob n’est pas spécifié, il y a équiprobabilité ;
Lois de probabilité usuelles

Sous R, les lois de probabilité usuelles sont notées comme suit :

Loi Notation Paramètres

uniforme unif min=0, max=1
binomiale binom size, prob
Poisson pois lambda
géométrique geom prob
normale norm mean=0, sd=1
exponentielle exp rate=1
Student t df
Fisher Snedecor f df1, df2
Khi-deux chisq df
Calcul de probabilité, fonction de répartition, quantile,
échantillon des lois usuelles
I dnotationdist :c’est la fonction de densité pour une loi continue
et la fonction de probabilité P(X = k) pour une loi discrète.

I rnotationdist : c’est la fonction qui génère des réalisations

aléatoires indépendantes de la loi ;

I pnotationdist : (avec l’option lower.tail=TRUE (par défaut))

c’est la fonction de répartition P(X ≤ x) ;

I pnotationdist : (avec l’option lower.tail=FALSE ) c’est la

fonction P(X > x) ;

I qnotationdist : c’est la fonction qui permet d’obtenir les

quantiles ;
Loi binomiale

Exemple
Soit X une v.a de loi B(35, 0.45)

dbinom(10, size = 35, prob = 0.45)# P(X = 10)

pbinom(15, size = 35, prob = 0.45)# P(X ≤ 15)

pbinom(11, size = 35, prob = 0.45, lower.tail = FALSE)#

P(X > 11)

sum(dbinom(16 :19, size = 35, prob = 0.457))#

P(16 ≤ X ≤ 19)pbinom(19, size = 35, prob = 0.45)-pbinom(15,

size = 35, prob = 0.45)# P(16 ≤ X ≤ 19)

Loi de Poisson

Exercice
Identifier ce que chacune des commandes suivantes effectue :

1. dpois(6, lambda=5)
2. dpois(0 :5,5)
3. dpois(c(2,4,8,10), 5)
4. ppois(6, 10, lower.tail = TRUE)
5. ppois(6, 10, F)
6. qpois(0.25, 10, T)
7. qpois(0.25, 10, F)
8. rpois(50,5).
Loi géométrique

Exercice
Quelle est la définition de la loi géométrique considérée par R ?

Identifier ce que chacune des commandes suivantes effectue :

1. dgeom(1,1/3)
2. dgeom(0 :9,1/3)
3. dgeom(1,c(1/3,2/3))
4. pgeom(9, 1/3,)
5. pgeom(9, 1/3, F)
6. qgeom(0.35, 1/3, T)
7. qgeom(0.35, 1/3, F)
8. rgeom(1,1/3)
Loi uniforme

Exercice
Identifier ce que chacune des commandes suivantes effectue :

1. runif(12,0,1)
2. runif(15)
3. dunif(0.5, min = 0, max = 1) ; dunif(1.5)
√
4. dunif(- 2,-3,5)
5. punif(-4,-3,5,F) ; punif(-3,-3,5,F)
6. punif(5,T)
7. punif(5.25,-3,5)
8. qunif(0.75, -3,5)
9. qunif(0.6, -3,5,F)
Loi exponentielle

Exercice
Identifier ce que chacune des commandes suivantes effectue :
1. dexp(3.2, rate = 1)
2. pexp(10, 1.5, T)
3. pexp(10, 1.5, F)
4. qexp(0.5, 1,T)
5. qexp(0.5, 1, F)
6. qexp(c(0.25,0.5,0.75), 1,T)
7. qexp(0.5, 2 :5, F)
8. rexp(20, 17.5)
Loi normale

Exercice
Identifier ce que chacune des commandes suivantes effectue :

1. dnorm(0, mean = 0, sd = 1) ; dnorm(0, 0, 1)

2. dnorm(2)
3. pnorm(5, 0, 1, T)
4. pnorm(5,T)
5. pnorm(5, 5, 1, F)
6. pnorm(5, 5, 1, T)
7. qnorm(0.95,10, 10, T)
8. qnorm(0.75)
9. qnorm(0.95,10, 10, F)
10. rnorm(10, 77, 25) ; rnorm(10)
Loi de Student

Exercice
Identifier ce que chacune des commandes suivantes effectue :

1. dt(0.5, 5)
2. pt(0.5,5,T)
3. pt(0.5,5,F)
4. pt(0.5,5,c(F,T)) ; pt(0.5,5 :10,T)
5. pt(0.5,5 :10,F)
6. qt(0.5,5,T) ; qt(0.5,5,F)
7. qt(0.5,5,c(F,T))
8. qt(0.5,5 :10,T)
9. qt(0.5,5 :10,F)
10. rt(25,5)
Représentations graphiques

plot qui réalise des représentations graphiques.

Exemple
plot(1 :20)
x=rnorm(50)
print(x)
summary(x)
y=2*x+7+rnorm(50,sd=0.6)
plot(x,y)
Quelques paramétrages de la fonction plot

I xlab (resp ylab) pour spécifier le titre de en abcisse (resp en

ordonnée) ;
I col pour spécifier la couleur ;
I xlim (resp ylim) pour spécifier les bornes de l’axe des abcisses
(resp des ordonnées) ;

Exemples
x<- seq(-pi,pi,0.05)
y <- cos(x)
plot(x,y,xlab="x",ylab="cos x")
plot(x,y,type="l", main="trait continu")
plot(x,y,type="s")
plot(x,y,col="blue") plot(x,y,xlim=c(-2*pi,2*pi))
plot(x,y,xlim=c(-2*pi,2*pi),ylim=c(-2,2))
Représentation d’une fonction

Pour représenter la fonction fun :

I plot(x,fun(x),type="l")
I plot(fun,min,max)
I curve(fun,min,max).

Exemples
x<- 0 :30
plot(x,log(xˆ 3+2*x+1),type="l")
curve(log(xˆ 3+2*x+1),0 :30)
plot(sin,-pi,pi)
Plusieurs représentations sur le même graphe
I matplot permet de représenter dans un même repère les
colonnes d’une matrice
I points permet d’ajouter des points à un graphique existant
I segments permet d’ajouter un segment joignant deux points à
un graphique.

Exemples
A=matrix(5 :16, 4)
matplot(A, type="o")
matplot(A, type="c")
matplot(A, type="l")
matplot(A, type="s")
points(1 :4,1 :4+5,col="blue")
points(1 :4,1 :4+5,col="blue",type="s")
curve(sin,-pi,pi)
segments(0,-1,0,0,col="red")
segments(-2,0.5,-1,-0.5,col="blue")
Plusieurs représentations sur le même graphe

I lines permet d’ajouter une ligne à un graphique existant.

I abline permet d’ajouter une ou plusieurs droites verticales ou
horizontales à un graphique.

Exemples
curve(tan(x),-5,5)
lines(-5 :5,10*(-5 :5))
abline(h=0,v=0,col="green")
Plusieurs représentations sur une même page
Pour représenter n × m graphiques sur n lignes et m colonnes dans
la même fenêtre, il suffit d’écrire :

par(mfrow=c(n,m))

avant les commandes permettant de représenter les n × m

graphiques.
Exemple
par(mfcol=c(1,2))
plot(dnorm,-5,5)
plot(dexp,0,5)

Exemple
par(mfrow=c(1,2))
plot(dnorm,-5,5)
plot(dexp,0,5)
Exporter un graphique au format jpeg

I Pour exporter un graphique en jpeg , il faut

1. écrire jpeg(filename="nomdufichier%d.jpeg")
2. réaliser le graphique
3. écrire dev.off() pour enregistrer le graphique dans le fichier
nomdufichier.jpeg.

I Pour exporter un graphique en bmp ou png , il faut procéder

de même en remplaçant jpeg par bmp ou par png.
I Pour exporter les graphiques en pdf, il faut procéder comme
suit
1. écrire pdf("nomdufichier.pdf", height=, width=)
2. réaliser le graphique
3. écrire dev.off() pour enregistrer le graphique le graphique dans
le fichier nomdufichier.pdf.
Importer/Exporter des données

I scan() permet de créer un vecteur à partir de données stockées

dans un fichier ;

I read.table() permet d’importer un fichier de données ;

I read.csv() et read.csv2() sont analogues à read.table() ;

I write.table() permet d’exporter des données.

Structures de contrôle : condition
I Pour exécuter des instructions si une condition est vraie et des
instructions alternatives sinon :

if (condition) instructions else instructions alternatives

I S’il n’y a pas d’instructions alternatives

if (condition) instructions

I Si les instructions se limitent à un seul calcul, on peut utiliser

ifelse

ifelse(condition, instructions, instructions alternatives)

Exemple
if (x>0) y=log(x) else y=0
y=ifelse(x>0,log(x),0)
Structures de contrôle : itération avec for

Pour exécuter des instructions pour les valeurs du compteur

appartenant à la séquence

for(compteur in séquence) instructions

Exemple
somme<-0
for(i in 1 :10)
{
somme<-somme+i}
somme
Structures de contrôle : itération avec while

Pour exécuter des instructions tant qu’une condition est vraie

while(condition) instructions

Exemple
i=1
somme<-0
while ( i <= 10 )
{
somme<-somme+i
i<-i+1
} somme
Structures de contrôle : itération avec repeat
Pour exécuter de façon répétée les instructions jusqu’à ce qu’une
condition soit vérifiée :

repeat (instructions) if(condition) break

Exemple
somme< −0
i=1
repeat
{
somme<-somme+i
i<-i+1
if(i>10) break
}
somme
Créer une fonction

mafonction=function (liste des paramètres)

{
instructions
return (l’objet contenant les résultats )
}
Pour exécuter la fonction :
mafonction(valeurs des paramètres)
Exemple
ecart-type <- function(x) sqrt(mean(xˆ 2)-mean(x)ˆ 2)
x<-c(9,5,2,3,7)
ecart-type(x)
Exercice 1
√
1. Créer un vecteur x de 50 réels choisis au hasard entre 2 et
25.
2. Arrondir les éléments de x à la 3ème décimale.
3. Extraire les éléments de x d’indice impair.
4. Afficher les éléments de x compris strictement entre 10 et 15.

Exercice 2
On lance 10 fois une pièce de monnaie bien équilibrée. Soit X le
nombre de Pile obtenus.
1. Calculer la probabilité d’obtenir exactement 4 fois Pile.
2. Calculer la probabilité d’obtenir au moins 4 fois Face.
3. Représenter la fonction de répartition de X .
4. Calculer la médiane.
Exercice 3
On lance une pièce de monnaie bien équilibrée jusqu’à obtenir
Face.Calculer
1. la probabilité d’obtenir Face au 5ème lancer.
2. la probabilité d’effectuer plus de 5 lancers pour obtenir Face.
3. le nombre moyen de lancers à effectuer pour obtenir Face.

Exercice 4
P500 3;
1. Calculer i=1 i
2. Effectuer 20 tirages successifs avec remise dans une population
de 100 personnes dont 45 hommes et 55 femmes.
3. Choisir au hasard et sans remise 15 observations dans un
ensemble 20 observations d’une loi uniforme sur [0, 1] et 15
observations d’une loi de Poisson P(5).
Exercice 5
1. Simuler 100 lancers indépendants d’un dé cubique équilibré
dont les faces sont numérotées de 6 à 12.
2. Former la distribution des fréquences de la série statistique
ainsi obtenue ;
3. Déterminer le mode de cette série.
4. Déterminer les déciles de cette série.

Exercice 6
1. Simuler l’épreuve qui consiste à lancer 25 fois une pièce de
monnaie équilibrée.
2. Calculer la probabilité d’obtenir Pile 5 fois à l’issue de cette
épreuve.
3. Calculer la probabilité d’obtenir Face au moins 2 fois à l’issue
de cette épreuve.
4. Calculer la probabilité d’obtenir Face plus de 5 fois et au plus
13 fois à l’issue de cette épreuve.
Exercice 7 √
Soit une variable aléatoire X ∼ N (2, 5, 0, 5).
1. Calculer P(|X | ≥ 4, 75);
2. Déterminer x tel que P(X ≤ x) = 0, 95 et y tel que
P(X > y ) = 0, 95.
3. Déterminer z tel que P(|X | ≤ z) = 0.975.
4. Représenter la densité de probabilité de X sur l’intervalle
[−20, 20] avec une subdivision régulière de pas 0,1.
5. Ajouter un segment vertical joignant le point (x,0) à la courbe.
Que représentent les aires situées à gauche et à droite du
segment, entre la courbe et l’axe des abcisses ?
Exercice 8
1. Simuler un échantillon (x1 , x2 , · · · , x1 000) d’une v.a. X de loi
de Bernoulli de paramètre 0.75.
2. Représenter l’histogramme de série en considérant 20 classes.
3. Tracer la courbe de m 7→ x̄m , ∀m = 1, 2, · · · , 1000.
4. Ajouter au graphique la droite d’équation y = 0.75.
5. Que constatez-vous ?
6. Quel résultat important de probabilité cette expérience illustre
t-elle ?
Exercice 9
airquality désigne un jeu de données disponible sous R.
1. Afficher airquality. De quel type d’objet s’agit-il ?
2. Déterminer le résumé des statistiques descriptive des airquality.
3. Représenter les boîtes à moustache de airquality.
4. Créer la matrice A des trois dernières colonnes de airquality.
5. Créer la matrice V des covariances des trois variables en
colonne dans A puis calculer le déterminant de V .
6. Afficher toutes les observations supérieures ou égales au 3ème
quartile de la 3ème colonne de A.
7. Afficher les observations des quatre premières variables
obtenues le 20 juillet.
8. Afficher les observations des quatre premières variables
obtenues le 20 de chaque mois considéré.
9. Retirer les lignes de airquality contenant des données
manquantes.