Programmation avec R

Nous rangeons cet ensemble de valeurs dans un objet appelé dé :

1:6 -> dé

Pour examiner le contenu d’un objet il suffit d’entrer son nom dans la console :

dé

[1] 1 2 3 4 5 6

Nous utilisons maintenant cet objet dé comme une urne dans laquelle nous faisons un

tirage aléatoire, c’est-à-dire au hasard, grâce à la fonction sample() :

1:6

sample(dé, 1) [1] 6

Le 1 après la virgule dans sample(dé, 1) signifie simplement que l’on ne veut faire qu’un seul tirage. Pour faire plusieurs lancers de dé sans avoir à retaper la commande en entier il suffit d’utiliser la touche ↑ pour invoquer une instruction antérieure.

sample(dé, 1)

[1] 1 sample(dé, 1) [1] 2 sample(dé, 1) [1] 6

On peut également faire plusieurs tirages d’un coup :

sample(dé, 3)

[1] 3 6 2

3 R est-il un bon simulateur ?

Nous venons de simuler trois jeux de hasard avec un ordinateur, c’est-à-dire un système au comportement complètement déterministe. Il y a quelque chose d’un peu paradoxal ici. On peut se convaincre que les résultats donnés par R ne doivent rien au hasard en contrôlant la condition initiale du générateur de nombres aléatoires grâce à la fonction set.seed() :

set.seed(1)

matrix(sample(dé), ncol = 2) [,1] [,2] [1,] 2 4

| Date : 14-11-2021.

[2,] 6 1 [3,] 3 5

set.seed(1)

matrix(sample(dé), ncol = 2) [,1] [,2] [1,] 2 4 [2,] 6 1 [3,] 3 5

On constate que l’on a obtenu exactement le même tirage. C’est pourquoi on parle de

tirage pseudo-aléatoires. Si on change la condition initiale on obtient des tirages différents :

set.seed(2)

matrix(sample(dé), ncol = 2) [,1] [,2] [1,] 2 1 [2,] 4 6 [3,] 3 5

set.seed(2)

matrix(sample(dé), ncol = 2) [,1] [,2] [1,] 2 1 [2,] 4 6 [3,] 3 5

! ! Propriétés souhaitables d’un tirage pseudo-aléatoire ! !

Ce que l’on souhaite c’est que les tirages pseudo-aléatoires se comportent comme si on avait un tirage réellement aléatoire. Reprenons le cas du lancer de dé. Il est raisonnable de souhaiter que le dé ne soit pas truqué et que toutes les faces aient la même probabilité, 1/6 , d’être obtenues. Est-ce le cas ? Faisons une expérience, reproductible en invoquant set.seed(1) au préalable, en lançant le dé six-cent- mille (600000) fois :

set.seed(1)

table(sample(dé, 600000, replace = TRUE)) -> obs obs 1 2 3 4 5 6 100129 99924 99734 100389 100207 99617

Dans cette expérience nous avons donc obtenu : (cid:70) La face 1 cent-mille-cent-vingt-neuf (100129) fois. (cid:70) La face 2 quatre-vingt-dix-neuf-mille-neuf-cent-vingt-quatre (99924) fois. (cid:70) La face 3 quatre-vingt-dix-neuf-mille-sept-cent-trente-quatre (99734) fois. (cid:70) La face 4 cent-mille-trois-cent-quatre-vingt-neuf (100389) fois. (cid:70) La face 5 cent-mille-deux-cent-sept (100207) fois. (cid:70) La face 6 quatre-vingt-dix-neuf-mille-six-cent-dix-sept (99617) fois.

| Date : 14-11-2021.

Une représentation graphique valant généralement mieux qu’un long discours on peut éga-

lement résumer les résultats de notre expérience avec un diagramme en bâtons :

barplot(obs, main = "Résultat de l’expérience", xlab = "Face du dé", ylab ="Nombre

de tirages obtenus")

Nous pouvons donc constater que le dé n’est pas truqué puisque nous avons obtenu à peu près chaque face cent-mille (100000) fois, soit avec une probabilité voisine de 1/6 , puisque nous avons fait six-cent-mille (600000) lancers en tout.

4 Loi, fonction de répartition et quantiles

Les fonctions usuelles pour les lois de probabilité sont composées de deux parties : (cid:70) un préfixe : d pour la densité, p pour la fonction de répartition et q pour la fonction

quantile.

(cid:70) un suffixe spécifiant la loi à laquelle on s’intéresse : binom pour la loi binomiale, norm

pour la loi normale, etc.

Les paramètres des lois sont précisés à l’intérieur de la fonction : (cid:70) size et prob pour la loi binomiale ; (cid:70) mean et sd pour la loi normale.

La simulation de données consiste typiquement à tirer un échantillon (X1, . . . , Xn) de va- riables indépendantes et identiquement distribuées (en abrégé i.i.d.) suivant la même loi qu’une variable aléatoire X. Les fonctions utilisées sur R pour simuler des données suivant une loi s’écrivent de la même manière que pour la densité, la fonction de répartition ou le quantile, mais cette fois avec le préfixe r.

Le tableau suivant récapitule les suffixes et les fonctions d, p, q et r des lois les plus utilisées

en statistique :

4.1 Génération d’un échantillon gaussien

(cid:70) retourne un vecteur de 100 valeurs tirées d’une distribution gaussienne de moyenne 0 et

d’écart-type 1.

rnorm(100,mean=0,sd=1) rnorm(100,0,1) rnorm(100)

(cid:70) obtenir une distribution de valeurs bimodale.

rnorm(1000,mean=c(0,3),sd=c(0.5,1))

4.2 Génération d’un échantillon suivant une loi uniforme

| Date : 14-11-2021.

suffixe

masse ou densité

répartition

quantile

Binomiale Binomiale négative Géométrique Hypergéométrique Poisson Uniforme Exponentielle Normale Log-normale Beta Gamma Cauchy Student Fisher Weibull χ2 Logistique

binom nbinom geom hyper pois unif exp norm lnorm beta gamma cauchy t f weibull chisq logis

dbinom dnbinom dgeom dhyper dpois dunif dexp dnorm dlnorm dbeta dgamma dcauchy dt df dweibull dchisq dlogis

pbinom pnbinom pgeom phyper ppois punif pexp pnorm plnorm pbeta pgamma pcauchy pt pf pweibull pchisq plogis

qbinom qnbinom qgeom qhyper qpois qunif qexp qnorm qlnorm qbeta qgamma qcauchy qt qf pweibull qchisq qlogis

génération de nombres aléatoires rbinom rnbinom rgeom rhyper rpois runif rexp rnorm rlnorm rbeta rgamma rcauchy rt rf rweibull rchisq rlogis

4.2 Génération d’un échantillon suivant une loi uniforme

(cid:70) retourne 100 valeurs issues d’une distribution uniforme entre 0 et 5.

runif(100,0,5) runif(100,min=0,max=5)

(cid:70) obtenir une distribution de valeurs bimodale.

runif(100,c(0,1),c(5,7))

4.3 Exemples

Donner toutes les valeurs de la loi binomiale de paramètres n = 10 et p = 1/3 (B(10, 1/3)).

options(digits = 7) # par défaut

dbinom(0:10, 10, 1/3) [1] 1.734153e-02 8.670765e-02 1.950922e-01 2.601229e-01 2.276076e-01 1.365645e-01 [7] 5.690190e-02 1.625768e-02 3.048316e-03 3.387018e-04 1.693509e-05

options(digits = 4)

dbinom(0:10, 10, 1/3) [1] 1.734e-02 8.671e-02 1.951e-01 2.601e-01 2.276e-01 1.366e-01 5.690e-02 1.626e-02 [9] 3.048e-03 3.387e-04 1.694e-05

Quelle est la probabilité d’obtenir 1 avec une loi binomiale B(10, 1/3) ?

dbinom(1, 10, 1/3)

[1] 0.08671

4.4 Exercices

| Date : 14-11-2021.

Quelle est la probabilité d’obtenir plus de 45 et moins de 55 avec une loi binomiale

B(100, 1/2) ?

Calculer P(45 < X < 55) pour X ∼ B(100, 1/2)

sum(dbinom(46:54, 100, 0.5)) [1] 0.6318

Quelle est la probabilité d’obtenir au plus 1 avec une loi binomiale pour B(10, 1/3) ?

pbinom(1, 10, 1/3)

[1] 0.104

dbinom(0, 10,1/3)

[1] 0.01734 dbinom(1, 10,1/3) [1] 0.08671 dbinom(0, 10,1/3) + dbinom(1, 10,1/3) [1] 0.104

N.B. : Noter que d donne les valeurs P(X = k) et que p donne les valeurs P(X (cid:54) x).

4.4 Exercices

4.4.1 Exercice 1

(cid:70) Quelle est la probabilité de dépasser strictement 4 pour une loi de Poisson de paramètre

2.7 (P(2.7)) ?

(cid:70) Quelle est la probabilité de dépasser 1.96 pour une loi normale centrée réduite (N (0, 1)) ? (cid:70) Quelle est la valeur x telle que P (X (cid:54) x) = 0.975 pour une loi normale (quantile) ? (cid:70) Quel est le quantile 1 % pour une loi T de Student à 5 ddl ? (cid:70) Donner un échantillon aléatoire simple de 10 valeurs d’une loi de Poisson de paramètre

égale à 2.7.

(cid:70) Donner un échantillon aléatoire simple de 10 valeurs d’une loi normale réduite. (cid:70) Donner un échantillon aléatoire simple de 10 valeurs d’une loi du Khi2 à 2 ddl. (cid:70) Donner un échantillon aléatoire simple de 10 valeurs d’une loi binomiale n = 100 et

p = 1/2.

4.4.2 Exercice 2

1. Soit la variable aléatoire X qui suit la loi normale standard.

a) Calculer les probabilités suivantes :

i) P(X < 1.45).

ii) P(X > 1.45).

iii) P(|X| < 2.05).

4.4 Exercices

| Date : 14-11-2021.

b) Trouver la valeur de u dans les cas suivants :

i) P(X < u) = 0.6.

ii) P(X > u) = 0.6.

2. Soit la variable aléatoire X qui suit la loi normale de paramètres µ = 30 et σ = 1.

Calculer P(X < 1.45) et P(X < 31.45).

3. Soit la variable aléatoire X qui suit la loi normale de paramètres µ = 0 et σ = 10.

Calculer P(X < 1.45) et P(X < 14.45).

#–––––––––––––––––––––––––––––––––––

Solution Exercise 2 #–––––––––––––––––––––––––––––––––––

#–––––––––––- question 1a –––––––––––––––––-

pnorm(1.45) pnorm(1.45, lower.tail=FALSE) pnorm(2.05)- pnorm(-2.05) pnorm(2.05)- pnorm(2.05, lower.tail=FALSE) #–––––––––––- question 1b –––––––––––––––––- qnorm(0.60) qnorm(0.60, lower.tail=FALSE) #–––––––––––- question 2 –––––––––––––––––- pnorm(1.45, mean=30) pnorm(31.45, mean=30) #–––––––––––- question 3 –––––––––––––––––- pnorm(1.45, mean=0,sd=10) pnorm(14.5, mean=0,sd=10)

4.4.3 Exercice 3

On jette dix fois de suite une pièce équilibrée et on s’intéresse à la variable aléatoire X

correspondant au nombre de Pile obtenus sur les 10 lancers.

1. Quelle est la loi de X ?

2. Déterminer P(X = 7).

3. Représenter la loi de la variable X (argument type="h" de la fonction plot).

4. Soit Y le nombre de Pile obtenus après 10 lancers lorsque la pièce est truquée avec probabilité d’obtenir Pile égale à p = 0.75. Représenter sur deux graphiques différents mais dans une même fenêtre graphique les lois des variables X et Y .

5. Calculer la valeur de la fonction de répartition de la variable Y au point 7.

6. Représenter la fonction de répartition de Y (plot avec l’argument type="s").

7. Lire sur le graphique la médiane de Y .

8. Retrouver la valeur de cette médiane ainsi que les valeurs des quartiles.

4.4 Exercices

| Date : 14-11-2021.

#–––––––––––––––––––––––––––––––––––

Solution Exercise 3 #––––––––––––––––––––––––––––––––––– # Question 1 : X suit une loi binomiale B(10,0.5) # Question 2 dbinom(7,size=10,prob=0.5) ou directement dbinom(7,10,0.5) # Question 3 k=0:10 p=dbinom(k,10,0.5) plot(k,p,type="h") # Question 4 k=0:10 px=dbinom(k,10,0.5) py=dbinom(k,10,0.75) par(mfrow=c(1,2)) plot(k,px,type="h",ylab="Loi de X") plot(k,py,type="h",ylab="Loi de Y") # X est stochastiquement inferieure a Y par(mfrow=c(1,1)) # Question 5 pbinom(7,10,0.75) # Question 6 k=0:10 F=pbinom(k,10,0.75) plot(k,F,type="s",ylab="cdf de Y") # Question 7 abline(h=0.5,lty=3) # mediane de Y = 8 # Question 8 qbinom(0.5,10,0.75) qbinom(c(0.25,0.5,0.75),10,0.75)