1. Grâce à ce formulaire vous pouvez proposer une astuce ou un script sur R.
  2. Votre script doit pouvoir être lancé en l'état. Veuillez penser à :
    • inclure le chargement des "library" nécessaires
    • construire un petit jeu de données si besoin est.
    • commenter les lignes de codes pour en faciliter la compréhension.
  3. Le titre de votre script ou astuce doit être clair et explicite.
  4. Pensez à mettre votre code entre les balises [R] et [/R]; Pour cela, vous pouvez utiliser le bouton

Vous pouvez utiliser vos comptes Facebook, twitter ou google pour vous identifer (google est compatible yahoo, openID...)
L'ideal étant de vous connecter si vous avez un compte utilisateur, ou faire une demande de compte utilisateur si vous n'en avez pas encore.
Créer un compte va vous permettre de pouvoir éditer vos codes et de mettre en avant votre site internet.
Sinon vous pouvez soumettre anonymement en remplissant les champs ci-après.







Choisissez les catégories correspondantes à votre Code:

  • algorithmique
  • Analyse de survie
  • base indispensable
  • bayésien
  • configuration de R
  • exportation de données
  • fonctions utiles
  • graphique
  • importation de données
  • inférence
  • manipulation de données
  • message d'erreur
  • modélisation
  • Non classé
  • optimisation
  • planification
  • programmer avec R
  • regression linéaire
  • Test
  • Transformation de données


Les scripts R de Melen

Comment afficher la palette de couleur ? palette

1
dans graphique
- ça ne sert à rien -- c\'est interessant - (score de +3 sur 3 votes)
Loading ... Loading ...
La palette correspond à l’ensemble des couleurs actives.

Vous pouvez l’afficher et la modifier avec la fonction palette()

 
palette()
 


Il existe de nombreuses palettes prédéfinies : hsv(), gray(), rainbow(), heat.colors(), terrain.colors(), topo.colors(), cm.colors().
Vous pouvez ainsi modifier la palette :

 
palette(gray(1:10/10))
palette()
 


Enfin il existe des librairies comme RColorBrewer qui permettent d’utiliser différentes palettes adaptées à divers graphiques.

Proposé par Melen.

1 commentaire. Cliquez ici pour réagir.

Comment afficher la liste des couleurs prédéfinies dans R ? colors

0
dans base indispensable, graphique
- ça ne sert à rien -- c\'est interessant - (score de +1 sur 1 votes)
Loading ... Loading ...
Vous en avez assez d’utiliser toujours les mêmes couleurs pour vos graphiques?
R a 657 couleurs prédéfinies…

 
 
colors()
 



Proposé par Melen.

Pas encore de commentaire, cliquez ici pour réagir.

Comment ajuster des fonctions avec les moindres carrés ? nls

7
dans inférence, modélisation
- ça ne sert à rien -- c\'est interessant - (score de +11 sur 13 votes)
Loading ... Loading ...
La manière la plus simple d’ajuster une fonction à des données est la méthode « géométrique » des moindres carrés (minimiser la somme des carrés des écarts correspond à maximiser la vraisemblance avec une loi normale).

La fonction nls de R permet de réaliser ceci de manière simple.

Voyons deux exemples :

 
#exemple modèle de croissance exponentiel
 
#on crée une fonction qui correspond à un modèle de Malthus
 
malthus<-function(t,N0,r){N0*exp(r*t)}
 
t<-seq(0:100)
NO<-2
r<-0.05
 
plot(malthus(t,NO,r)~t,type='l',col='green',lwd=2)
 
#on crée des données en ajoutant du bruit
sim<-malthus(t,NO,r)+rnorm(t,sd=0.3*malthus(t,NO,r))
plot(sim~t,pch=20)
 
#on ajuste la fonction sur les données simulées en utilisant les moindres carrés
fitmalthus<-nls(sim~malthus(t,a,b),start=list(a=1,b=0.01))
fitmalthus
summary(fitmalthus)
 
#on vérifie
plot(sim~t,pch=20)
lines(malthus(t,NO,r)~t,type='l',col='green',lwd=2)
lines(malthus(t,coef(fitmalthus)[1],coef(fitmalthus)[2])~t,type='l',col='red',lwd=2)
 
#exemple 2 avec une fonction logistique
 
logist<-function(t,v,w,z){v/(1+w*exp(-z*t))}
 
v<-10;w<-8;z<-0.09;x<-seq(from=-100,to=100)
 
plot(logist(x,v,w,z)~x,type='l',col='blue',lwd=2)
 
#on crée des données avec du bruit
noisy<-rnorm(x,mean=logist(x,v,w,z),sd=0.9)
 
plot(noisy,pch=20)#notez qu'il y a des valeurs négatives pas réalistes...
 
fitlog<-nls(noisy~logist(x,a,b,c),start=list(a=5,b=5,c=0.02))
fitlog
summary(fitlog)
 
#on vérifie
plot(noisy~x,pch=20,col="grey")
lines(logist(x,v,w,z)~x,type='l',col='blue',lwd=2)
lines(logist(x,coef(fitlog)[1],coef(fitlog)[2],coef(fitlog)[3])~x,type='l',col='red',lwd=2) 
 


Ces deux exemples cachent en fait de nombreuses difficultés :

Pour que l’algorithme converge, vous devez donner des valeurs de départ des paramètres assez proches de la réalité –> essayez d’ajuster vos courbes à l’œil afin d’obtenir les meilleurs paramètres de départ possibles. Sinon vous pouvez utiliser la fonction nls2 et l’algorithme brute force, ou passer à la vraisemblance voire au merveilleux monde des mcmc…

La fonction nls présente, comme beaucoup de fonctions R, de nombreux arguments qui vous permettront peut être d’ajuster vos fonctions comme vous le souhaitez : control, trace, weights etc.

Enfin vous pouvez récupérer de nombreuses informations à partir de l’objet créé lors de l’ajustement : résidus, valeur des paramètres estimés etc…

Proposé par Melen.

Déjà 7 commentaires. Cliquez ici pour réagir.

Comment obtenir le critère d’Akaike avec R ? : AIC

1
dans modélisation
- ça ne sert à rien -- c\'est interessant - (score de +7 sur 7 votes)
Loading ... Loading ...
Le critère d’Akaike (AIC) est un critère utilisé pour la sélection de modèles. Ce critère représente un compromis entre le biais diminuant avec le nombre de paramètres libres et la parcimonie, volonté de décrire les données avec le plus petit nombre de paramètres possible. Il se calcule de la façon suivante -2*log-likelihood + k*npar. Par défaut on a souvent k=2. Le meilleur modèle est celui qui possède l’AIC le plus faible.

On obtient ce critère en utilisant la fonction AIC(objet,k=?), k=2 par défaut. Prenons un exemple

 
 
 
 library(MASS) #pour la fonction fitdistr
 
#  z est un vecteur contenant les données, on essaie de modéliser ces données par une distribution de weibull ou gamma
 
z<-c(14,14,14,17,17,17,17,17,17,17,17,17,17,17,17,17,17,17,17,17,17,17,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,20,23)
 
# on utilise la fonction fitdistr pour une loi Weibull, regarder ?fitdistr
 
paraw <- fitdistr(z,densfun="weibull") 
loglikw<-logLik(paraw)  # on peut avoir le loglikelihood
loglikw
 
#on fait la même chose pour une loi gamma
 
parag<-fitdistr(z,densfun="gamma")
loglikg<-logLik(parag)
loglikg
 
#AIC avec k=2 par défaut
 
AIC(paraw)#environ 204
AIC(parag)#environ 209
 
#d'après ce critère on choisirait la loi de weibull
 
#on vérifie en recalculant les AIC directement
 
akaike<-function(npar,loglik,k){-2*loglik+k*npar}
 
akaike(2,loglikw[1],2)
akaike(2,loglikg[1],2)
 
 
 


L’AIC est un critère comme les autres, faites en bon usage et n’en n’abusez pas trop!


Proposé par Melen.

1 commentaire. Cliquez ici pour réagir.

Comment représenter des fonctions 3D avec persp ?

0
dans graphique
- ça ne sert à rien -- c\'est interessant - (score de +2 sur 2 votes)
Loading ... Loading ...
Il existe plusieurs façons de représenter des surfaces en 3D. Je vous en propose une qui utilise la fonction de base persp().

 
 
#on construit les fonctions à représenter f1 et f2
 
 f1<-function(x,y){
    15*sin(sqrt(x^2+y^2))
}
 
 
################
 
mu1<-0 
mu2<-0 
s11<-15 
s12<-20 
s22<-10 
rho<-0.5 
 
f2<-function(x,y){
 
term1<-1/(2*pi*sqrt(s11*s22*(1-rho^2)))
term2<--1/(2*(1-rho^2))
term3<-(x-mu1)^2/s11
term4<-(y-mu2)^2/s22
term5<--2*rho*((x-mu1)*(y-mu2))/(sqrt(s11)*sqrt(s22))
term1*exp(term2*(term3+term4-term5))}
 
# on définit deux vecteurs correspondant aux axes x et y
 
x<-seq(-15,15,length=50)
y<-x
 
# on calcule la valeur de z=f(x,y) pour tous les couples x[i],y[i] avec la fonction outer
 
z1<-outer(x,y,f1)
z2<-outer(x,y,f2)
 
#on utilise la fonction persp
 
x11()
persp(x,y,z1,theta=30,phi=40,expand=0.5,col="lightblue",ticktype="detailed")
x11()
persp(x,y,z2,theta=40,phi=30,expand=0.5,col="lightgreen",ticktype="detailed")
 
 
 


Il est possible d’ajouter des points au graphique en utilisant points(trans3d()).
Regardez bien tous les arguments de la fonction persp() : ?persp. Ils sont nombreux et vous permettront d’obtenir votre surface sous l’angle que vous voulez.

Enfin il existe comme souvent d’autres fonctions pour représenter des surfaces : plot3d du package rgl, plot.surface du package fields…à vous de choisir!

bonne 3D!

Proposé par Melen.

Pas encore de commentaire, cliquez ici pour réagir.

Comment ajouter des courbes, des points et des droites sur un graphique ? lines, points,abline

0
dans graphique
- ça ne sert à rien -- c\'est interessant - (score de +7 sur 11 votes)
Loading ... Loading ...
R permet de faire beaucoup de chose avec les graphiques mais il faut coder ce qu’on veut faire. Nous allons voir ici comment ajouter des courbes, des points et des droites sur un plot.

 
 
#on crée des données pour l'exemple
 
 x<-seq(0:100)
a<-2
b<-5
 
y1<-a*x+b
y2<-a*x^0.5+b
y3<-a*x^0.3+b*x
 
ynoisy1<-y1+rnorm(length(y1),sd=0.2*y1)
ynoisy2<-y2+rnorm(length(y2),sd=0.2*y2)
ynoisy3<-y3+rnorm(length(y3),sd=0.2*y3)
 
#étape 1 on trace la fonction y1 avec un plot
 
plot(y1~x,type='l',col="green",lwd=2,ylim=c(0,300))
 
#étape 2 on veut ajouter les fonctions y2 et y3 sur le même graphiques : on utilise lines
 
lines(y2~x,type='l',col="blue",lwd=2)
lines(y3~x,type='l',col="purple",lwd=2)
 
#étape 3 on veut ajouter les nuages de point ynoisy1 2 et 3 sur le graphique: on utilise points
 
points(ynoisy1~x,pch=20,col="grey")
points(ynoisy2~x,pch=20,col="black")
points(ynoisy3~x,pch=22,col="black")
 
#étape 4 on veut ajouter une droite verticale pour x=50 : on utilise abline(v=) v pour vertical
 
 abline(v=50,col="red")
 
#étape 5 on veut une droite horizontale à y=150 : idem avec h =
 
abline(h=150,col="red")
 
#étape 6 on veut une droite d'équation y= -2*x+300 : abline avec c(ordonnée,pente)
abline(c(300,-2),col="black") 
 
 


Notez que ce graphique est très moche…

Proposé par Melen.

Pas encore de commentaire, cliquez ici pour réagir.

Comment enlever les étoiles mises par défaut dans les résultats de p value?

0
dans configuration de R, regression linéaire
- ça ne sert à rien -- c\'est interessant - (pas encore de vote)
Loading ... Loading ...
Dans les résultats d’analyse statistique, R affiche souvent des étoiles à côté des p values avec le code suivant :

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Il est possible de les enlever:

 
 
#on utilise les données iris
 
iris
summary(iris)
 
#pour l'exemple on fait une analyse de variance : la longueur des sépales est-elle expliquée par l'espèce?
 
reg<-lm(iris$Sepal.Length~iris$Species)
anova(reg)
summary(reg)
 
 
 


R affiche les résultats avec les étoiles pour les p-values :

Analysis of Variance Table

Response: iris$Sepal.Length
              Df Sum Sq Mean Sq F value    Pr(>F)    
iris$Species   2 63.212  31.606  119.26 < 2.2e-16 ***
Residuals    147 38.956   0.265      

On décide de les enlever en utilisant options :

 
 
options(show.signif.stars=FALSE)
 
 
summary(reg)
anova(reg)
 
 


R affiche les résultats sans les étoiles :

Response: iris$Sepal.Length
              Df Sum Sq Mean Sq F value    Pr(>F)
iris$Species   2 63.212  31.606  119.26 < 2.2e-16
Residuals    147 38.956   0.265     


Proposé par Melen.

Pas encore de commentaire, cliquez ici pour réagir.

Comment enlever les outliers d’un boxplot?

6
dans base indispensable
- ça ne sert à rien -- c\'est interessant - (score de +4 sur 4 votes)
Loading ... Loading ...
Les boxplots mettent parfois en évidence des individus qu’on peut qualifier d’atypiques ou outliers.

Un fois mis en évidence graphiquement on peut les repérer et si nécessaire les enlever.

 
 
#on crée un jeu de donnée 
 
b1<-c(0.1, 0.2,6,5,5,6,7,8,8,9,9,9,10,10,25)
 
#on trace le boxplot
 
boxplot(b1) #il y a 3 outliers 
 
#on met le boxplot dans un objet box
 
box<-boxplot(b1)
boxplot(b1)
 
#box$out donne les outliers
 
#on crée des nouvelles données sans les outliers
 
b2<-b1[-which(b1%in%box$out)]
 
#on vérifie
 
boxplot(b2)
 
 


Proposé par Melen.

Déjà 6 commentaires. Cliquez ici pour réagir.

Comment choisir la graine pour générer des nombres aléatoires dans R ? : set.seed

0
dans base indispensable
- ça ne sert à rien -- c\'est interessant - (score de +4 sur 4 votes)
Loading ... Loading ...
Lorsqu’on génère des nombres aléatoires il est souvent utile de pouvoir choisir la graine (pour resimuler exactement un modèle stochastique par exemple).

Sous R ceci se fait via la fonction set.seed(« nombre entier »).

 
 
set.seed(3)
runif(3)
 
set.seed(4567)
runif(3)
 
set.seed(3)
runif(3)
 
 

Proposé par Melen.

Pas encore de commentaire, cliquez ici pour réagir.

Comment tester la conformité à une loi avec le test de Kolmogorov-Smirnov ?

2
dans Test
- ça ne sert à rien -- c\'est interessant - (score de +3 sur 5 votes)
Loading ... Loading ...
Le test de Kolmogorov-Smirnov est un test d’hypothèse utilisé pour décider si un échantillon suit une loi de probabilité donnée ou si deux échantillons suivent la même loi.

Sous R on peut réaliser ce test avec la fonction ks.test()

 
 
#on crée des échantillons
 
a<-rnorm(100,mean=0,sd=1)
 
b<-rgamma(100,shape=1,rate=0.8)
 
c<-rnorm(50,mean=0,sd=1)
 
#a et b proviennent-ils de la même loi?
 
ks.test(a,b)#p=7.5e-11 on rejette l'hypothèse nulle
 
#a et c?
 
ks.test(a,c)#p=0.35 on accepte l'hypothèse nulle
 
#a provient-il d'une loi gamma avec 3 comme paramètre de forme et 2 pour le taux?
 
ks.test(a,"pgamma",3,2)#p value très faible on rejette l'hypothèse
 
#a provient-il d'une loi normale?
 
ks.test(a,"pnorm")#p=0.13 on accepte l'hypothèse
 
 
 


Comme pour tous les tests, faites bien attention à ce que veut réellement dire ce test, prenez du recul sur la notion de p-value et ne basez pas vos analyses sur ce test seul !!!

Proposé par Melen.

Déjà 2 commentaires. Cliquez ici pour réagir.
Formation logiciel R