Denis Puthier
'2015-10-07'
Premiers pas
x <- 15
print(x)
[1] 15
show(x)
[1] 15
x
[1] 15
x = 23
# x <- 32
x
[1] 23
Les instructions peuvent être séparées par un retour à la ligne ou par le caractère ;.
x <- 12; y <- 13
cat(x) # print(x)
12
cat(y) # print(y)
13
ls()
[1] "x" "y"
rm(list=ls())
# Arg1 et arg2 (...) correspondent au noms des arguments
# a et b correspondent aux objets que l’on souhaite passer à la fonction
nomdelafonction(arg1= a, arg2= b,...)
y <- c(-50, 1:8, 70) ; print(y)
[1] -50 1 2 3 4 5 6 7 8 70
length(y)
[1] 10
is(y)
[1] "numeric" "vector"
mean(y)
[1] 5.6
mean(y, trim = 0.2) # 10% deleted at each ends
[1] 4.5
mean(y, tr = 0.2)
[1] 4.5
mean(y, 0.2)
[1] 4.5
noms <- c("celine", "alain", "robert")
class(noms) # is(noms)
[1] "character"
logic <- c(TRUE, FALSE, TRUE) # logic <- c(T, F, T)
class(logic) # is(logic)
[1] "logical"
help(mean) # ?mean
x <- 1:10 # une suite d'entiers
y <- c(2, 5, 7) # combiner 2 5 et 7 dans un vecteur
z <- rep(x, times = 5) # répéter x cinq fois.
On peut utiliser par exemple:
par(mfrow=c(2, 2))
x <- rnorm(n=1000, mean=10, sd=2)
hist(x, main="Distribution normale (m=10, sd=2)")
x <- runif(n=1000, min = 10, max=20)
hist(x, main="Distribution uniforme (min=10, max=20)")
x <- rpois(1000, lambda=10)
hist(x, main="distribution de Poisson (lambda=10)")
x <- rhyper(nn=1000, m=100, n=100, k = 20)
hist(x, main="Distribution hypergéométrique (m=100, n=100, k = 20)")
par(mfrow=c(1, 1))
x <- 1:10
sort(x, decreasing = TRUE)
[1] 10 9 8 7 6 5 4 3 2 1
sample(x)
[1] 10 6 8 3 2 4 5 7 1 9
sample(x, replace = TRUE)
[1] 6 8 4 3 10 2 8 2 5 3
Pour interroger les éléments d’un vecteur on utilisera:
1. Un vecteurs d’indices (ie; de positions)
* les positions sont notées de 1 à n (taille du vecteur).
2. Un vecteur logique
3. les noms des éléments du vecteur.
x <- c(-2, 4, -5, 1, 7)
x[3] # position 3
[1] -5
x[3:5] # position 3 à 5
[1] -5 1 7
x[c(2,5)] # position 2 et 5
[1] 4 7
x[-3] # x sans la position 4 (-5)
[1] -2 4 1 7
x <- c(-2, 4, -5, 1, 7)
x < 4
[1] TRUE FALSE TRUE TRUE FALSE
x[x < 4] # renvoie les positions pour lesquelles le test est vrai
[1] -2 -5 1
x[x < 4 & x > 0] # & correspond au et logique
[1] 1
x[x < 4 | x == 7] # | correspond au ou logique, == à l'opérateur d'égalité
[1] -2 -5 1 7
x <- c(-2, 4, -5, 1, 7)
names(x) <- letters[1:length(x)]
x
a b c d e
-2 4 -5 1 7
x["a"]
a
-2
x[c("a", "d")]
a d
-2 1
Fonction read.table
Fonction read.table, arguments importants
Exercice:
base.dir <- "http://dputhier.github.io/STAT1_R/data"
file.name <- "yeast_up500_6nt-1str-ovlp_MET_vs_all_sub.tsv"
yeast.df <- read.table(file.path(base.dir,file.name),
sep="\t",
header=T,
row=1)
class(yeast.df)
[1] "data.frame"
mon.data.frame[ligne, colonne]
yeast.df[1,] # première ligne
yeast.df[1:10, 1:2] # 10 premières lignes des colonnes 1 et 2
yeast.df[c(1,3), 1] # lignes 1 et 3 de la colonne 1
yeast.df[ ,"MET_occ"] # la colonne MET_occ
yeast.df$MET_occ # la colonne MET_occ
yeast.df$MET_occ[1:3] # les trois premiers éléments de la colonne MET_occ
Exemple: le k-mère W = TTCACG
\( N_{W,all}=424 \) Nombre d'occurrences d'un k-mère dans l'ensemble des séquences promotrices (colonne “all_occ”).
\( N_{W,MET}=7 \) Nombre d'occurrences dans les gènes de la méthionine (colonne “MET_occ”“)
\( N_{tous,all}=2421832 \) Nombre total de mots comptés dans tous les promoteurs
\( N_{tous,MET}=9350 \) Nombre total de mots comptés dans les promoteurs de méthionine
sum_met <- sum(yeast.df$MET_occ)
print(sum_met)
[1] 9350
sum_all <- sum(yeast.df$all_occ)
print(sum_all)
[1] 2421832
for( i in vec){
faireQuelqueChose(i)
...
}
x <- 1:10
y <- -x
z <- vector()
for (i in 1:length(x)){
z[i] <- x[i] + y[i]
}
z <- x + y
z
[1] 0 0 0 0 0 0 0 0 0 0
# Autres exemples
x/2
[1] 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0
x^2
[1] 1 4 9 16 25 36 49 64 81 100
Définissons par \( P(W|all) \), ou \( p_W \), la probabilité d'observer un mot (\( W \)) donné à une position donnée, dans l'ensemble des promoteurs (\( all \)). On appelle cette probabilité probabilité a priori (prior).
Pour le mot “TTCACG”, on obtient
\( P(W|all) = p_W = \frac{N_{W,all}}{N_{tous,all}} = \frac{424}{2421832} = 0.000175 \)
En vous basant sur ce que nous venons de voir, calculez la fréquence de chaque mots dans la classe MET et la class ALL (variables MET_freq et all_freq).
Ajouter ces variables comme nouvelles colonnes au data.frame
Testez que le résultat obtenu est conforme pour le k-mer “TTCACG”
MET_freq <- yeast.df$MET_occ / sum_met
all_freq <- yeast.df$all_occ / sum_all
yeast.df$MET_freq <- MET_freq
yeast.df$all_freq <- all_freq
head(yeast.df, n=2)
MET_occ MET_mseq all_occ all_mseq MET_freq all_freq
CACGTG 16 14 317 285 0.0017112299 1.308926e-04
CGTGAC 6 6 197 192 0.0006417112 8.134338e-05
yeast.df["TTCACG","all_freq"]
[1] 0.0001750741
round(yeast.df["TTCACG","all_freq"], 6) == 0.000175
[1] TRUE
On peut calculer le nombre d'occurrences attendues de la manière suivante:
yeast.df$MET_exp_occ <- yeast.df$all_freq * sum_met
plot(yeast.df$all_occ, yeast.df$MET_occ)
plot(yeast.df$all_occ, yeast.df$MET_occ,
pch=16, col="blue", panel.first=grid())
col.palette <- densCols(yeast.df$all_occ, yeast.df$MET_occ)
plot(x=yeast.df$all_occ,
y=yeast.df$MET_occ,
pch=16,
col=col.palette,
panel.first=grid(),
xlab="Nombre d'occurence dans l'ensemble des promoters",
ylab="Nombre d'occurence dans les promoters MET",
main="Nombre d'occurrences des k-mers"
)
lm <- lm(yeast.df$MET_occ ~ yeast.df$all_occ)
abline(lm, col="red", lwd=2)