silico.biotoul.fr
 

M1 Traitement de Donnees Biologiques - TP 2 R

From silico.biotoul.fr

(Difference between revisions)
Jump to: navigation, search
m (Lecture et exploration du tableau de données)
 
(35 intermediate revisions not shown)
Line 6: Line 6:
* Estimation d'un intervalle de confiance d'un paramètre (moyenne, proportion)
* Estimation d'un intervalle de confiance d'un paramètre (moyenne, proportion)
-
<span style='color: #990000'>Créer un répertoire de travail sur le bureau (par exemple TDB-TP2) et commencez par télécharger le fichier source que vous allez utiliser et compléter pour générer le compte rendu de TP : [[Media:M1.TDB.TP_reglin_proba_ic_R.Rmd|M1.TDB.TP_reglin_proba_ic_R.Rmd]] (click droit de la souris -- enregistrer la cible sous...).</span> Ouvrez le logiciel RStudio et chargez ce fichier puis lancez sa compilation pour voir le compte rendu. Pour cela cliquez sur le bouton Knit HTML ou bien utilisez la combinaison de touches Ctrl + shift + K.  
+
<span style='color: #990000'>
 +
Créer un répertoire de travail sur le bureau (par exemple TDB-TP2) et commencez par télécharger le fichier source que vous allez utiliser et compléter pour générer le compte rendu de TP : [[Media:M1.TDB.TP_reglin_proba_ic_R.Rmd|M1.TDB.TP_reglin_proba_ic_R.Rmd]] (click droit de la souris -- enregistrer la cible sous...).
 +
</span>
 +
Ouvrez le logiciel RStudio et chargez ce fichier puis lancez sa compilation pour voir le compte rendu. Pour cela cliquez sur le bouton Knit HTML ou bien utilisez la combinaison de touches Ctrl + shift + K.
= Régression linéaire =
= Régression linéaire =
-
<span style='color: #990000'>On cherche à détecter une corrélation entre deux variables quantitatives mesurées sur les mêmes individus.</span>
+
On cherche à détecter une corrélation entre deux variables quantitatives mesurées sur les mêmes individus.
Notre exemple: pour différentes espèces de bactéries, nous connaissons (i) la taille du génome grâce au séquençage, (ii) le nombre de séquences codantes prédites par un algorithme bioinformatique.
Notre exemple: pour différentes espèces de bactéries, nous connaissons (i) la taille du génome grâce au séquençage, (ii) le nombre de séquences codantes prédites par un algorithme bioinformatique.
Téléchargez et sauvegardez le fichier [[Media:bacterial_genomes.txt|bacterial_genomes.txt]] dans le répertoire crée pour ce TP (click droit de la souris -- enregistrer la cible sous...).
Téléchargez et sauvegardez le fichier [[Media:bacterial_genomes.txt|bacterial_genomes.txt]] dans le répertoire crée pour ce TP (click droit de la souris -- enregistrer la cible sous...).
-
Quelle question peut-on se poser ?
+
<span style='color: #990000'>
 +
Quelle(s) question(s) peut-on se poser ?
 +
</span>
== Lecture et exploration du tableau de données ==
== Lecture et exploration du tableau de données ==
Line 22: Line 27:
</source>
</source>
-
Accéder directement aux variables simplement en donnant leurs noms:
+
Accéder directement aux colonnes du data.frame <tt>genomes</tt> simplement à partir de leurs noms (fonction <tt>attach</tt>):
<source lang='rsplus'>
<source lang='rsplus'>
attach(genomes)
attach(genomes)
Line 29: Line 34:
ORF_number
ORF_number
</source>
</source>
 +
 +
<span style='color: #990000'>
 +
Quelle est l'unité de mesure pour <tt>Genome_size</tt> ?
 +
</span>
Faire un graphique:
Faire un graphique:
Line 40: Line 49:
cov(Genome_size,ORF_number)
cov(Genome_size,ORF_number)
</source>
</source>
-
* coefficient de corrélation r de Pearson
+
* coefficient de corrélation ''r'' de Pearson
<source lang='rsplus'>
<source lang='rsplus'>
cor(Genome_size,ORF_number)
cor(Genome_size,ORF_number)
</source>
</source>
-
== Régression linéaire ==
+
== Régression ==
-
Calculer une fonction linéaire qui relie les 2 variables, avec la commande <span style='color: #990000;'>lm()</span>:
+
<span style='color: #990000;'>
 +
Calculer une fonction linéaire qui relie les 2 variables, avec la commande <tt>lm()</tt>
 +
</span>:
<source lang='rsplus'>
<source lang='rsplus'>
lm(ORF_number ~ Genome_size) # remarque: c'est une régression de "y" sur "x", d'où lm(y~x)
lm(ORF_number ~ Genome_size) # remarque: c'est une régression de "y" sur "x", d'où lm(y~x)
Line 62: Line 73:
</source>
</source>
 +
<span style='color: #990000;'>
Quelle est l'équation de la droite de régression ?
Quelle est l'équation de la droite de régression ?
 +
</span>
Calculez le coefficient de détermination R2 (% de variance expliquée par le modèle linéaire ==> bien si > 70%):
Calculez le coefficient de détermination R2 (% de variance expliquée par le modèle linéaire ==> bien si > 70%):
Line 82: Line 95:
== Calculs basés sur des données ==
== Calculs basés sur des données ==
-
L'exemple porte sur la taille des arbres dans une forêt de sequoia.
+
L'exemple porte sur la taille des arbres dans une forêt de séquoias.
Téléchargez et sauvegardez le fichier [[Media:sequoia.txt|sequoia.txt]] (click droit de la souris -- enregistrer la cible sous...).
Téléchargez et sauvegardez le fichier [[Media:sequoia.txt|sequoia.txt]] (click droit de la souris -- enregistrer la cible sous...).
 +
<span style='color: #990000;'>
 +
Chargez le fichier avec la commande <tt>read.table</tt> dans une variable <tt>sequoia</tt>, puis utilisez la commande <tt>attach</tt> pour accéder directement aux noms des colonnes. Utilisez également la commande <tt>names</tt> pour connaître le nom des colonnes dans le fichier.
 +
</span>
 +
 +
<!--
Lecture du tableau de données:
Lecture du tableau de données:
<source lang='rsplus'>
<source lang='rsplus'>
Line 90: Line 108:
attach(sequoia);names(sequoia)
attach(sequoia);names(sequoia)
</source>
</source>
 +
-->
Représentation graphique des données:  
Représentation graphique des données:  
Line 96: Line 115:
</source>
</source>
 +
<span style='color: #990000;'>
Quelle est la probabilité qu'un arbre mesure 80m ?
Quelle est la probabilité qu'un arbre mesure 80m ?
 +
</span>
 +
<source lang='rsplus'>
<source lang='rsplus'>
length(taille_arbre[taille_arbre==80])/length(taille_arbre)
length(taille_arbre[taille_arbre==80])/length(taille_arbre)
</source>
</source>
 +
ou mieux :
 +
<source lang='rsplus'>
 +
length( which(taille_arbre==80) ) / length(taille_arbre)
 +
</source>
 +
 +
<span style='color: #990000;'>
Quelle est la probabilité qu'un arbre mesure plus de 100 ?
Quelle est la probabilité qu'un arbre mesure plus de 100 ?
 +
</span>
 +
<source lang='rsplus'>
<source lang='rsplus'>
-
length(taille_arbre[taille_arbre>100])/length(taille_arbre)
+
length( which(taille_arbre>100) ) / length(taille_arbre)
</source>
</source>
== Calculs basés sur des lois (ou distributions) de probabilité  ==
== Calculs basés sur des lois (ou distributions) de probabilité  ==
-
Différentes lois sont programmées (binomiale,Poisson,normale,chi2...)
+
Différentes lois sont programmées (binomiale, Poisson, normale, chi2 ou <math>\chi^2</math>, ...)
Il y a plusieurs fonctions pour chaque loi. Par exemple, pour la loi normale:
Il y a plusieurs fonctions pour chaque loi. Par exemple, pour la loi normale:
-
* dnorm() : fonction densité (density)
+
* <tt>dnorm()</tt> : fonction densité ('''d'''ensity)
-
* pnorm() : fonction de répartition (probability)
+
* <tt>pnorm()</tt> : fonction de répartition ('''p'''robability)
-
* qnorm() : fonction quantile (quantile)
+
* <tt>qnorm()</tt> : fonction quantile ('''q'''uantile)
-
* rnorm() : générateur aléatoire (random)
+
* <tt>rnorm()</tt> : générateur aléatoire ('''r'''andom)
 +
 
 +
=== Calculs basés sur une loi normale centrée réduite <math>N(0,1)</math> ===
-
=== Calculs basés sur un loi normale centrée réduite N(0,1) (0 et 1 = moyenne et écart-type de la variable X) ===
+
'''Rappels :'''
 +
* Une '''loi normale''' est définie par '''2 paramètres''' : la '''moyenne''' et l''''écart-type'''.
 +
* La '''loi normale centrée réduite''' est la loi normale dont la moyenne est '''0''' et l'écart-type '''1'''. Elle est notée généralement <math>N(0,1)</math> (0 pour la moyenne et 1 pour l'écart-type de la variable X).
-
[[Image:loi_normale.jpeg]]
+
[[Image:loi_normale.jpeg|thumb|400px|link=|Fonction de densité de la loi normale centrée réduite <math>N(0,1)</math>]]
-
Prob[X=-1]:
+
Probabilité de X = 0, notée <math>P(X = 0)</math>
<source lang='rsplus'>
<source lang='rsplus'>
-
dnorm(-1,0,1)
+
dnorm(0,0,1)
</source>
</source>
-
Prob[X<-1]:
+
<math>P(X \le 0)</math> :
<source lang='rsplus'>
<source lang='rsplus'>
-
pnorm(-1,0,1)
+
pnorm(0,0,1)
</source>
</source>
-
Prob[X>-1]:
+
<math>P(X > 0)</math> :
<source lang='rsplus'>
<source lang='rsplus'>
-
1-pnorm(-1,0,1) # ou pnorm(-1,0,1, lower.tail=F)
+
1-pnorm(0,0,1) # ou pnorm(0,0,1, lower.tail=F)
</source>
</source>
-
Echantillonage + histogramme de 100 valeurs dans une N(0,1):
+
Histogramme pour 100 valeurs échantillonnées aléatoirement suivant une <math>N(0,1)</math> :
<source lang='rsplus'>
<source lang='rsplus'>
hist(rnorm(100,0,1))
hist(rnorm(100,0,1))
</source>
</source>
-
=== Calculs basés sur un loi Binomiale B(n,p) ===
 
-
Exemple: pour un lot de n=100 graines d'Arabidopsis ayant chacune une probabilité de germination p=0.8 (B(100,0.8)), je veux:
+
 
 +
'''Application à la forêt de séquoias'''
 +
 
 +
<span style='color: #990000;'>
 +
Utilisez la moyenne et l'écart-type des tailles des séquoias chargées précédemment pour calculer la probabilité qu'un arbre mesure 80m selon la loi normale paramétrée ainsi.
 +
</span>
 +
 
 +
<span style='color: #990000;'>
 +
Quelle est la probabilité qu'un séquoia mesure plus de 110m ?
 +
</span>
 +
 
 +
<!--span style='color: #990000;'>
 +
Quelle est la probabilité qu'un arbre pris au hasard dans cette forêt et qui mesure 70m soit un séquioa  ?
 +
</span-->
 +
 
 +
=== Calculs basés sur une loi Binomiale <math>B(n,p)</math> ===
 +
 
 +
Exemple: pour un lot de n=100 graines d'Arabidopsis ayant chacune une probabilité de germination p=0.8 donc selon <math>B(100,0.8)</math>
La probabilité que k=80 graines germent (k=nombre de succès)
La probabilité que k=80 graines germent (k=nombre de succès)
Line 153: Line 203:
</source>
</source>
-
La probabilité que plus de 80 graines germente
+
La probabilité que plus de 80 graines germent
<source lang='rsplus'>
<source lang='rsplus'>
pbinom(80,100,prob=0.8, lower.tail=F)
pbinom(80,100,prob=0.8, lower.tail=F)
</source>
</source>
-
<span style='color: #990000;'>Ajoutez ces parties à votre compte rendu.</span>
+
<span style='color: #990000'>
 +
En moyenne, combien de graines devraient germer ?
 +
</span>
 +
La distribution Binomiale <math>B(100,0.8)</math>:
 +
<source lang='rsplus'>
 +
k=seq(0,100,1)    # k=nombre de graines qui germent
 +
plot(dbinom(k,100,prob=0.8),xlab="k",ylab="Prob(X=k)")
 +
</source>
 +
 +
 +
 +
<span style='color: #990000;'>
 +
Ajoutez ces parties à votre compte rendu.
 +
</span>
= Estimer un intervalle de confiance (IC) d'un paramètre populationnel (µ, p,...) à partir d'un échantillon =
= Estimer un intervalle de confiance (IC) d'un paramètre populationnel (µ, p,...) à partir d'un échantillon =
-
== IC de la moyenne µ de la population, à partir de la moyenne m calculée sur un échantillon ==
+
== Intervalle de confiance de la moyenne théorique µ de la population, à partir de la moyenne m calculée sur un échantillon ==
-
<span style='color: #990000;'>Calculer l'IC du poids moyen d'un récolte de tomate cerises, à partir d'un échantillon.</span>
+
Calculer l'intervalle de confiance (IC) du poids moyen théorique d'une récolte de tomates cerises, à partir d'un échantillon.
 +
<span style='color: #990000;'>
Téléchargez et sauvegardez le fichier le fichier [[Media:tomates_cerises.txt|tomates_cerises.txt]] (click droit de la souris -- enregistrer la cible sous...).
Téléchargez et sauvegardez le fichier le fichier [[Media:tomates_cerises.txt|tomates_cerises.txt]] (click droit de la souris -- enregistrer la cible sous...).
 +
</span>
 +
<span style='color: #990000;'>
 +
Chargez le fichier avec la commande <tt>read.table</tt> dans une variable <tt>tomates</tt>, puis utilisez la commande <tt>attach</tt> pour accéder directement aux noms des colonnes. Utilisez également la commande <tt>names</tt> pour connaître le nom des colonnes dans le fichier. Affichez le contenu de la colonne <tt>poids_tomate</tt>.
 +
</span>
 +
<!--
Lecture du tableau de données:
Lecture du tableau de données:
<source lang='rsplus'>
<source lang='rsplus'>
Line 173: Line 242:
attach(tomates);names(tomates)
attach(tomates);names(tomates)
poids_tomate
poids_tomate
-
</source>
+
</source> -->
Script pour calculer l'IC:
Script pour calculer l'IC:
Line 181: Line 250:
n = length(poids_tomate)# taille de l'échantillon
n = length(poids_tomate)# taille de l'échantillon
</source>
</source>
-
<span style='color: #990000;'># si n<100, on utilise la distribution de Student pour calculer l'erreur.</span>
+
 
 +
Si n&le;100, on utilise la distribution de Student pour calculer la marge d'erreur. Ici, nous allons calculer l'intervalle ayant 95% de chance de contenir la vraie moyenne de la population (et pas celle de l'échantillon).
<source lang='rsplus'>
<source lang='rsplus'>
-
error <- qt(0.975,df=n-1)*s/sqrt(n) # qt(0.975) = quantile 97,5% de la distribution de student; s/sqrt(n) =erreur standard de la moyenne
+
error <- qt(0.975,df=n-1)*s/sqrt(n) # qt(0.975) = quantile 97,5% de la distribution de Student ; s/sqrt(n) = erreur standard de la moyenne
left  <- m-error
left  <- m-error
right <- m+error
right <- m+error
Line 192: Line 262:
</source>
</source>
-
<span style='color: #990000;'># si n>100 on utilise la distribution normale pour calculer l'erreur.</span>
+
Si n>100 on utilise la distribution normale pour calculer l'erreur.
<source lang='rsplus'>
<source lang='rsplus'>
error <- qnorm(0.975)*s/sqrt(n) # qnorm(0.975) = quantile 97,5% de la distribution normale
error <- qnorm(0.975)*s/sqrt(n) # qnorm(0.975) = quantile 97,5% de la distribution normale
Line 198: Line 268:
right <- m+error
right <- m+error
</source>
</source>
-
L'intervalle ayant 95% de chance de contenir la vraie moyenne µ de la population est alors
+
L'intervalle ayant 95% de chance de contenir la vraie moyenne µ de la population est alors:
<source lang='rsplus'>
<source lang='rsplus'>
left;right  
left;right  
</source>
</source>
-
<span style='color: #990000;'> REMARQUE IMPORTANTE: si l'on veut un IC à 99% ==> remplacer 0.975 par 0.995.</span>
+
'''REMARQUE IMPORTANTE :''' si l'on veut un IC à 99% ==> remplacer 0.975 par 0.995.
 +
== IC de la proportion théorique p d'un caractère dans la population, à partir de la fréquence f calculée sur un échantillon ==
-
== IC de la fréquence f d'un caractère dans la population, à partir de la proportion p calculée sur un échantillon ==
+
Calcul de l'IC de la proportion théorique de pommes rouges dans une récolte comportant des pommes rouges et vertes, à partir d'un échantillon de n=125 pommes contenant f=0.4 de pommes rouges (40%).
-
 
+
On utilise une approximation normale de la distribution d'une proportion :
-
<span style='color: #990000;'> Calculer l'IC de la fréquence de pommes rouges dans une récolte comportant des pommes rouges ET vertes, à partir d'un échantillon de n=125 pommes contenant p=40% de pommes rouges</span>
+
-
On utilise une approximation normale de la distribution d'une proportion:
+
<source lang='rsplus'>
<source lang='rsplus'>
-
p <- 0.4 # proportion de pommes rouges dans l'échantillon
+
f <- 0.4 # fréquence de pommes rouges dans l'échantillon
n <- 125 # taille de l'échantillon
n <- 125 # taille de l'échantillon
-
error <- qnorm(0.975)*sqrt((p*(1-p)/n))
+
error <- qnorm(0.975)*sqrt(f*(1-f)/n)
-
left  <- p-error
+
left  <- f-error
-
right <- p+error
+
right <- f+error
</source>
</source>
-
L'intervalle ayant 95% de chance de contenir la vraie proportion de la population est:
+
L'intervalle ayant 95% de chance de contenir la vraie proportion p de la population est:
<source lang='rsplus'>
<source lang='rsplus'>
left;right
left;right
</source>
</source>
-
Faites varier la taille de l'échantillon ainsi que la confiance ...
+
<span style='color: #990000;'>
 +
Faites varier la taille de l'échantillon ainsi que la confiance.
 +
</span>
-
 
+
<span style='color: #990000;'>Ajoutez tout cela au compte rendu de TP avant de l'envoyer à votre enseignant par mail ([mailto:bonhomme@lrsv.ups-tlse.fr bonhomme@lrsv.ups-tlse.fr] ou [mailto:roland.barriot@univ-tlse3.fr roland.barriot@univ-tlse3.fr]). Le compte rendu est à envoyer '''avant de commencer le TP3'''. Envoyez les 2 fichiers (.Rmd et .html). Envoyez-vous aussi le mail en copie pour pouvoir vérifier que tout est bien passé. Mettez un titre tel que "Compte rendu TP2 TDB de -et votre Nom et Prénom-".</span>
-
<span style='color: #990000;'>Ajoutez tout cela au compte rendu de TP avant de l'envoyer à votre enseignant par mail ([mailto:bonhomme@lrsv.ups-tlse.fr bonhomme@lrsv.ups-tlse.fr] ou [mailto:barriot@biotoul.fr barriot@biotoul.fr]). Le compte rendu est à envoyer '''avant de commencer le TP3'''. Envoyez les 2 fichiers (.Rmd et .html). Envoyez-vous aussi le mail en copie pour pouvoir vérifier que tout est bien passé. Mettez un titre tel que "Compte rendu TP2 TDB de -et votre Nom et Prénom-".</span>
+
= Liens =
= Liens =
Line 232: Line 302:
* Site de R : http://www.r-project.org et sites miroirs (dont ceux en France) pour télécharger le logiciel et les librairies : https://cran.r-project.org/mirrors.html
* Site de R : http://www.r-project.org et sites miroirs (dont ceux en France) pour télécharger le logiciel et les librairies : https://cran.r-project.org/mirrors.html
* RStudio : https://www.rstudio.com
* RStudio : https://www.rstudio.com
 +
* Utilisation de R depuis un navigateur : http://www.r-fiddle.org
 +
Chargement des données avec l'adresse des fichiers
 +
<source lang='rsplus'>
 +
genomes=read.table("http://silico.biotoul.fr/site/images/d/de/Bacterial_genomes.txt", sep="\t", header=TRUE)
 +
sequoia=read.table("http://silico.biotoul.fr/site/images/e/e1/Sequoia.txt", sep="\t", header=TRUE)
 +
tomates=read.table("http://silico.biotoul.fr/site/images/3/3d/Tomates_cerises.txt", sep="\t", header=TRUE)
 +
</source>

Current revision as of 08:10, 21 September 2021

Contents

Régression linéaire, probabilités et intervalles de confiance

Ce TP comporte trois parties:

  • Régression linéaire
  • Calculs simples de probabilité
  • Estimation d'un intervalle de confiance d'un paramètre (moyenne, proportion)

Créer un répertoire de travail sur le bureau (par exemple TDB-TP2) et commencez par télécharger le fichier source que vous allez utiliser et compléter pour générer le compte rendu de TP : M1.TDB.TP_reglin_proba_ic_R.Rmd (click droit de la souris -- enregistrer la cible sous...). Ouvrez le logiciel RStudio et chargez ce fichier puis lancez sa compilation pour voir le compte rendu. Pour cela cliquez sur le bouton Knit HTML ou bien utilisez la combinaison de touches Ctrl + shift + K.

Régression linéaire

On cherche à détecter une corrélation entre deux variables quantitatives mesurées sur les mêmes individus. Notre exemple: pour différentes espèces de bactéries, nous connaissons (i) la taille du génome grâce au séquençage, (ii) le nombre de séquences codantes prédites par un algorithme bioinformatique. Téléchargez et sauvegardez le fichier bacterial_genomes.txt dans le répertoire crée pour ce TP (click droit de la souris -- enregistrer la cible sous...).

Quelle(s) question(s) peut-on se poser ?

Lecture et exploration du tableau de données

Lecture:

genomes=read.table("bacterial_genomes.txt", sep="\t", header=TRUE)

Accéder directement aux colonnes du data.frame genomes simplement à partir de leurs noms (fonction attach):

attach(genomes)
names(genomes)
Genome_size
ORF_number

Quelle est l'unité de mesure pour Genome_size ?

Faire un graphique:

plot(Genome_size,ORF_number,pch=16)

Quantifier la relation entre ces 2 variables:

  • covariance
cov(Genome_size,ORF_number)
  • coefficient de corrélation r de Pearson
cor(Genome_size,ORF_number)

Régression

Calculer une fonction linéaire qui relie les 2 variables, avec la commande lm() :

lm(ORF_number ~ Genome_size) # remarque: c'est une régression de "y" sur "x", d'où lm(y~x)

Gardons en mémoire le résultat de la régression:

reglin=lm(ORF_number ~ Genome_size)

On peut vérifier la significativité des coefficients de la droite avec:

summary(reglin)

Quelle est l'équation de la droite de régression ?

Calculez le coefficient de détermination R2 (% de variance expliquée par le modèle linéaire ==> bien si > 70%):

cor(Genome_size,ORF_number)^2 #(stocké aussi dans summary(reglin))

Représenter le nuage de points avec la droite de régression:

plot(Genome_size,ORF_number,pch=16)			
abline(reglin,col="red",lwd=2)

Ajoutez ces parties à votre compte rendu.


Calculs de probabilité

Calculs basés sur des données

L'exemple porte sur la taille des arbres dans une forêt de séquoias. Téléchargez et sauvegardez le fichier sequoia.txt (click droit de la souris -- enregistrer la cible sous...).

Chargez le fichier avec la commande read.table dans une variable sequoia, puis utilisez la commande attach pour accéder directement aux noms des colonnes. Utilisez également la commande names pour connaître le nom des colonnes dans le fichier.


Représentation graphique des données:

hist(taille_arbre)

Quelle est la probabilité qu'un arbre mesure 80m ?

length(taille_arbre[taille_arbre==80])/length(taille_arbre)

ou mieux :

length( which(taille_arbre==80) ) / length(taille_arbre)


Quelle est la probabilité qu'un arbre mesure plus de 100 ?

length( which(taille_arbre>100) ) / length(taille_arbre)

Calculs basés sur des lois (ou distributions) de probabilité

Différentes lois sont programmées (binomiale, Poisson, normale, chi2 ou χ2, ...) Il y a plusieurs fonctions pour chaque loi. Par exemple, pour la loi normale:

  • dnorm() : fonction densité (density)
  • pnorm() : fonction de répartition (probability)
  • qnorm() : fonction quantile (quantile)
  • rnorm() : générateur aléatoire (random)

Calculs basés sur une loi normale centrée réduite N(0,1)

Rappels :

  • Une loi normale est définie par 2 paramètres : la moyenne et l'écart-type.
  • La loi normale centrée réduite est la loi normale dont la moyenne est 0 et l'écart-type 1. Elle est notée généralement N(0,1) (0 pour la moyenne et 1 pour l'écart-type de la variable X).
Fonction de densité de la loi normale centrée réduite N(0,1)

Probabilité de X = 0, notée P(X = 0)

dnorm(0,0,1)

P(X \le 0) :

pnorm(0,0,1)

P(X > 0) :

1-pnorm(0,0,1) # ou pnorm(0,0,1, lower.tail=F)

Histogramme pour 100 valeurs échantillonnées aléatoirement suivant une N(0,1) :

hist(rnorm(100,0,1))


Application à la forêt de séquoias

Utilisez la moyenne et l'écart-type des tailles des séquoias chargées précédemment pour calculer la probabilité qu'un arbre mesure 80m selon la loi normale paramétrée ainsi.

Quelle est la probabilité qu'un séquoia mesure plus de 110m ?


Calculs basés sur une loi Binomiale B(n,p)

Exemple: pour un lot de n=100 graines d'Arabidopsis ayant chacune une probabilité de germination p=0.8 donc selon B(100,0.8)

La probabilité que k=80 graines germent (k=nombre de succès)

dbinom(80,100,prob=0.8)

La probabilité que au maximum 80 graines germent

pbinom(80,100,prob=0.8)

La probabilité que plus de 80 graines germent

pbinom(80,100,prob=0.8, lower.tail=F)

En moyenne, combien de graines devraient germer ?

La distribution Binomiale B(100,0.8):

k=seq(0,100,1)    # k=nombre de graines qui germent
plot(dbinom(k,100,prob=0.8),xlab="k",ylab="Prob(X=k)")


Ajoutez ces parties à votre compte rendu.

Estimer un intervalle de confiance (IC) d'un paramètre populationnel (µ, p,...) à partir d'un échantillon

Intervalle de confiance de la moyenne théorique µ de la population, à partir de la moyenne m calculée sur un échantillon

Calculer l'intervalle de confiance (IC) du poids moyen théorique d'une récolte de tomates cerises, à partir d'un échantillon. Téléchargez et sauvegardez le fichier le fichier tomates_cerises.txt (click droit de la souris -- enregistrer la cible sous...).

Chargez le fichier avec la commande read.table dans une variable tomates, puis utilisez la commande attach pour accéder directement aux noms des colonnes. Utilisez également la commande names pour connaître le nom des colonnes dans le fichier. Affichez le contenu de la colonne poids_tomate.

Script pour calculer l'IC:

m = mean(poids_tomate)  # moyenne calculée sur l'échantillon
s = sd(poids_tomate)    # écart-type calculé sur l'échantillon
n = length(poids_tomate)# taille de l'échantillon

Si n≤100, on utilise la distribution de Student pour calculer la marge d'erreur. Ici, nous allons calculer l'intervalle ayant 95% de chance de contenir la vraie moyenne de la population (et pas celle de l'échantillon).

error <- qt(0.975,df=n-1)*s/sqrt(n) # qt(0.975) = quantile 97,5% de la distribution de Student ; s/sqrt(n) = erreur standard de la moyenne
left  <- m-error
right <- m+error

L'intervalle ayant 95% de chance de contenir la vraie moyenne µ de la population est:

left;right

Si n>100 on utilise la distribution normale pour calculer l'erreur.

error <- qnorm(0.975)*s/sqrt(n) # qnorm(0.975) = quantile 97,5% de la distribution normale
left  <- m-error
right <- m+error

L'intervalle ayant 95% de chance de contenir la vraie moyenne µ de la population est alors:

left;right

REMARQUE IMPORTANTE : si l'on veut un IC à 99% ==> remplacer 0.975 par 0.995.

IC de la proportion théorique p d'un caractère dans la population, à partir de la fréquence f calculée sur un échantillon

Calcul de l'IC de la proportion théorique de pommes rouges dans une récolte comportant des pommes rouges et vertes, à partir d'un échantillon de n=125 pommes contenant f=0.4 de pommes rouges (40%). On utilise une approximation normale de la distribution d'une proportion :

f <- 0.4	# fréquence de pommes rouges dans l'échantillon
n <- 125	# taille de l'échantillon
error <- qnorm(0.975)*sqrt(f*(1-f)/n)
left  <- f-error
right <- f+error

L'intervalle ayant 95% de chance de contenir la vraie proportion p de la population est:

left;right

Faites varier la taille de l'échantillon ainsi que la confiance.


Ajoutez tout cela au compte rendu de TP avant de l'envoyer à votre enseignant par mail (bonhomme@lrsv.ups-tlse.fr ou roland.barriot@univ-tlse3.fr). Le compte rendu est à envoyer avant de commencer le TP3. Envoyez les 2 fichiers (.Rmd et .html). Envoyez-vous aussi le mail en copie pour pouvoir vérifier que tout est bien passé. Mettez un titre tel que "Compte rendu TP2 TDB de -et votre Nom et Prénom-".

Liens

Chargement des données avec l'adresse des fichiers

genomes=read.table("http://silico.biotoul.fr/site/images/d/de/Bacterial_genomes.txt", sep="\t", header=TRUE)
sequoia=read.table("http://silico.biotoul.fr/site/images/e/e1/Sequoia.txt", sep="\t", header=TRUE)
tomates=read.table("http://silico.biotoul.fr/site/images/3/3d/Tomates_cerises.txt", sep="\t", header=TRUE)