Ajouter l'équation de la ligne de régression et R ^ 2 sur le graphique

228

Je me demande comment ajouter l'équation de la ligne de régression et R ^ 2 sur le ggplot. Mon code est:

library(ggplot2)

df <- data.frame(x = c(1:100))
df$y <- 2 + 3 * df$x + rnorm(100, sd = 40)
p <- ggplot(data = df, aes(x = x, y = y)) +
            geom_smooth(method = "lm", se=FALSE, color="black", formula = y ~ x) +
            geom_point()
p

Toute aide sera grandement appréciée.

MYaseen208
la source
1
Pour les graphiques en treillis , voir latticeExtra::lmlineq().
Josh O'Brien

Réponses:

234

Voici une solution

# GET EQUATION AND R-SQUARED AS STRING
# SOURCE: https://groups.google.com/forum/#!topic/ggplot2/1TgH-kG5XMA

lm_eqn <- function(df){
    m <- lm(y ~ x, df);
    eq <- substitute(italic(y) == a + b %.% italic(x)*","~~italic(r)^2~"="~r2, 
         list(a = format(unname(coef(m)[1]), digits = 2),
              b = format(unname(coef(m)[2]), digits = 2),
             r2 = format(summary(m)$r.squared, digits = 3)))
    as.character(as.expression(eq));
}

p1 <- p + geom_text(x = 25, y = 300, label = lm_eqn(df), parse = TRUE)

ÉDITER. J'ai trouvé la source d'où j'ai choisi ce code. Voici le lien vers le message d'origine dans les groupes Google ggplot2

Production

Ramnath
la source
1
@ Le commentaire de JonasRaedle à propos de l'amélioration des textes avec annotateétait correct sur ma machine.
IRTFM
2
Cela ne ressemble en rien à la sortie publiée sur ma machine, où l'étiquette est remplacée autant de fois que les données sont appelées, ce qui donne un texte d'étiquette épais et flou. Passer les étiquettes à un data.frame fonctionne d'abord (voir ma suggestion dans un commentaire ci-dessous.
PatrickT
@PatrickT: supprimez le aes(et le correspondant ). aessert à mapper des variables de trame de données à des variables visuelles - ce n'est pas nécessaire ici, car il n'y a qu'une seule instance, vous pouvez donc tout mettre dans l' geom_textappel principal . Je vais modifier cela dans la réponse.
naught101
Le problème avec cette solution semble être que si le jeu de données est plus grand (le mien était 370000 observations), la fonction semble échouer. Je recommanderais la solution de @kdauria qui fait de même, mais beaucoup plus rapidement.
Benjamin
3
pour ceux qui veulent des valeurs r et p au lieu de R2 et de l'équation: eq <- substitute (italic (r) ~ "=" ~ rvalue * "," ~ italic (p) ~ "=" ~ pvalue, list (rvalue = sprintf ("% .2f", signe (coef (m) [2]) * sqrt (résumé (m) $ r.squared)), pvalue = format (résumé (m) $ coefficients [2,4], chiffres = 2 )))
Jerry T
135

J'ai inclus une statistique stat_poly_eq()dans mon package ggpmiscqui permet cette réponse:

library(ggplot2)
library(ggpmisc)
df <- data.frame(x = c(1:100))
df$y <- 2 + 3 * df$x + rnorm(100, sd = 40)
my.formula <- y ~ x
p <- ggplot(data = df, aes(x = x, y = y)) +
   geom_smooth(method = "lm", se=FALSE, color="black", formula = my.formula) +
   stat_poly_eq(formula = my.formula, 
                aes(label = paste(..eq.label.., ..rr.label.., sep = "~~~")), 
                parse = TRUE) +         
   geom_point()
p

entrez la description de l'image ici

Cette statistique fonctionne avec tout polynôme sans termes manquants et, espérons-le, a suffisamment de flexibilité pour être généralement utile. Les étiquettes R ^ 2 ou R ^ 2 ajustées peuvent être utilisées avec n'importe quelle formule de modèle équipée de lm (). Étant une statistique ggplot, elle se comporte comme prévu à la fois avec les groupes et les facettes.

Le package 'ggpmisc' est disponible via CRAN.

La version 0.2.6 vient d'être acceptée au CRAN.

Il répond aux commentaires de @shabbychef et @ MYaseen208.

@ MYaseen208 cela montre comment ajouter un chapeau .

library(ggplot2)
library(ggpmisc)
df <- data.frame(x = c(1:100))
df$y <- 2 + 3 * df$x + rnorm(100, sd = 40)
my.formula <- y ~ x
p <- ggplot(data = df, aes(x = x, y = y)) +
   geom_smooth(method = "lm", se=FALSE, color="black", formula = my.formula) +
   stat_poly_eq(formula = my.formula,
                eq.with.lhs = "italic(hat(y))~`=`~",
                aes(label = paste(..eq.label.., ..rr.label.., sep = "~~~")), 
                parse = TRUE) +         
   geom_point()
p

entrez la description de l'image ici

@shabbychef Il est maintenant possible de faire correspondre les variables de l'équation à celles utilisées pour les étiquettes d'axe. Pour remplacer le x par disons z et y par h, on utiliserait:

p <- ggplot(data = df, aes(x = x, y = y)) +
   geom_smooth(method = "lm", se=FALSE, color="black", formula = my.formula) +
   stat_poly_eq(formula = my.formula,
                eq.with.lhs = "italic(h)~`=`~",
                eq.x.rhs = "~italic(z)",
                aes(label = ..eq.label..), 
                parse = TRUE) + 
   labs(x = expression(italic(z)), y = expression(italic(h))) +          
   geom_point()
p

entrez la description de l'image ici

Étant ces expressions R normales, les lettres grecques peuvent désormais être utilisées à la fois dans les lhs et les rhs de l'équation.

[08/03/2017] @elarry Edit pour répondre plus précisément à la question d'origine, montrant comment ajouter une virgule entre l'équation et les étiquettes R2.

p <- ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "lm", se=FALSE, color="black", formula = my.formula) +
  stat_poly_eq(formula = my.formula,
               eq.with.lhs = "italic(hat(y))~`=`~",
               aes(label = paste(..eq.label.., ..rr.label.., sep = "*plain(\",\")~")), 
               parse = TRUE) +         
  geom_point()
p

entrez la description de l'image ici

[2019-10-20] @ helen.h Je donne ci-dessous des exemples d'utilisation de stat_poly_eq()avec groupement.

library(ggpmisc)
df <- data.frame(x = c(1:100))
df$y <- 20 * c(0, 1) + 3 * df$x + rnorm(100, sd = 40)
df$group <- factor(rep(c("A", "B"), 50))
my.formula <- y ~ x
p <- ggplot(data = df, aes(x = x, y = y, colour = group)) +
  geom_smooth(method = "lm", se=FALSE, formula = my.formula) +
  stat_poly_eq(formula = my.formula, 
               aes(label = paste(..eq.label.., ..rr.label.., sep = "~~~")), 
               parse = TRUE) +         
  geom_point()
p

p <- ggplot(data = df, aes(x = x, y = y, linetype = group)) +
  geom_smooth(method = "lm", se=FALSE, formula = my.formula) +
  stat_poly_eq(formula = my.formula, 
               aes(label = paste(..eq.label.., ..rr.label.., sep = "~~~")), 
               parse = TRUE) +         
  geom_point()
p

entrez la description de l'image ici

entrez la description de l'image ici

[2020-01-21] @Herman Cela peut être un peu contre-intuitif à première vue, mais pour obtenir une seule équation lors de l'utilisation du regroupement, il faut suivre la grammaire des graphiques. Limitez le mappage qui crée le regroupement à des calques individuels (illustré ci-dessous) ou conservez le mappage par défaut et remplacez-le par une valeur constante dans le calque où vous ne souhaitez pas le regroupement (par exemple colour = "black").

Suite de l'exemple précédent.

p <- ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "lm", se=FALSE, formula = my.formula) +
  stat_poly_eq(formula = my.formula, 
               aes(label = paste(..eq.label.., ..rr.label.., sep = "~~~")), 
               parse = TRUE) +         
  geom_point(aes(colour = group))
p

entrez la description de l'image ici

[2020-01-22] Par souci d'exhaustivité un exemple à facettes, démontrant que dans ce cas également les attentes de la grammaire des graphismes sont remplies.

library(ggpmisc)
df <- data.frame(x = c(1:100))
df$y <- 20 * c(0, 1) + 3 * df$x + rnorm(100, sd = 40)
df$group <- factor(rep(c("A", "B"), 50))
my.formula <- y ~ x

p <- ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "lm", se=FALSE, formula = my.formula) +
  stat_poly_eq(formula = my.formula, 
               aes(label = paste(..eq.label.., ..rr.label.., sep = "~~~")), 
               parse = TRUE) +         
  geom_point() +
  facet_wrap(~group)
p

entrez la description de l'image ici

Pedro Aphalo
la source
1
Il convient de noter que le xet ydans la formule se réfèrent aux données xet ydans les couches de l'intrigue, et pas nécessairement à celles de la portée au moment de la my.formulaconstruction. Ainsi, la formule devrait toujours utiliser des variables x et y?
shabbychef
Il est très vrai que xet se yréfèrent aux variables qui sont mappées à ces esthétiques. C'est l'attente également pour geom_smooth () et comment fonctionne la grammaire des graphiques. Il aurait pu être plus clair d'utiliser des noms différents dans le bloc de données, mais je les ai simplement conservés comme dans la question d'origine.
Pedro Aphalo
Sera possible dans la prochaine version de ggpmisc. Merci pour la suggestion!
Pedro Aphalo
3
Bon point @elarry! Ceci est lié au fonctionnement de la fonction parse () de R. Par essais et erreurs, j'ai trouvé que cela aes(label = paste(..eq.label.., ..rr.label.., sep = "*plain(\",\")~"))faisait l'affaire.
Pedro Aphalo
1
@HermanToothrot Habituellement, R2 est préféré pour une régression, il n'y a donc pas de r.label prédéfini dans les données renvoyées par stat_poly_eq(). Vous pouvez utiliser stat_fit_glance(), également à partir du package 'ggpmisc', qui retourne R2 comme valeur numérique. Consultez les exemples dans la page d'aide et remplacez-les stat(r.squared)par sqrt(stat(r.squared)).
Pedro Aphalo
99

J'ai changé quelques lignes de la source des stat_smoothfonctions associées pour créer une nouvelle fonction qui ajoute l'équation d'ajustement et la valeur R au carré. Cela fonctionnera également sur les tracés à facettes!

library(devtools)
source_gist("524eade46135f6348140")
df = data.frame(x = c(1:100))
df$y = 2 + 5 * df$x + rnorm(100, sd = 40)
df$class = rep(1:2,50)
ggplot(data = df, aes(x = x, y = y, label=y)) +
  stat_smooth_func(geom="text",method="lm",hjust=0,parse=TRUE) +
  geom_smooth(method="lm",se=FALSE) +
  geom_point() + facet_wrap(~class)

entrez la description de l'image ici

J'ai utilisé le code dans la réponse de @ Ramnath pour formater l'équation. La stat_smooth_funcfonction n'est pas très robuste, mais il ne devrait pas être difficile de jouer avec.

https://gist.github.com/kdauria/524eade46135f6348140 . Essayez de mettre à jour ggplot2si vous obtenez une erreur.

kdauria
la source
2
Merci beaucoup. Celui-ci ne fonctionne pas seulement pour les facettes, mais même pour les groupes. Je le trouve très utile pour les régressions par morceaux, par exemple stat_smooth_func(mapping=aes(group=cut(x.val,c(-70,-20,0,20,50,130))),geom="text",method="lm",hjust=0,parse=TRUE), en combinaison avec EvaluateSmooths de stackoverflow.com/questions/19735149/…
Julian
1
@aelwan, changez ces lignes: gist.github.com/kdauria/… comme vous le souhaitez. Ensuite, sourcele fichier entier dans votre script.
kdauria
1
@kdauria Et si j'ai plusieurs équations dans chacun des facet_wraps et j'ai différentes valeurs y dans facet_wrap. Des suggestions sur la façon de fixer les positions des équations? J'ai essayé plusieurs options de hjust, vjust et d'angle en utilisant cet exemple dropbox.com/s/9lk9lug2nwgno2l/R2_facet_wrap.docx?dl=0 mais je n'ai pas pu apporter toutes les équations au même niveau dans chacune des facet_wrap
brillant
3
@aelwan, la position de l'équation est déterminée par ces lignes: gist.github.com/kdauria/… . J'ai fait xposet des yposarguments de la fonction dans le Gist. Donc, si vous voulez que toutes les équations se chevauchent, définissez simplement xposet ypos. Sinon, xposet ypossont calculés à partir des données. Si vous voulez quelque chose de plus sophistiqué, il ne devrait pas être trop difficile d'ajouter de la logique à l'intérieur de la fonction. Par exemple, vous pourriez peut-être écrire une fonction pour déterminer quelle partie du graphique a l'espace le plus vide et y placer la fonction.
kdauria
6
J'ai rencontré une erreur avec source_gist: erreur dans r_files [[qui]]: type d'indice non valide «fermeture». Voir cet article pour la solution: stackoverflow.com/questions/38345894/r-source-gist-not-working
Matifou
73

J'ai modifié le post de Ramnath pour a) rendre plus générique afin qu'il accepte un modèle linéaire comme paramètre plutôt que le bloc de données et b) affiche les négatifs de manière plus appropriée.

lm_eqn = function(m) {

  l <- list(a = format(coef(m)[1], digits = 2),
      b = format(abs(coef(m)[2]), digits = 2),
      r2 = format(summary(m)$r.squared, digits = 3));

  if (coef(m)[2] >= 0)  {
    eq <- substitute(italic(y) == a + b %.% italic(x)*","~~italic(r)^2~"="~r2,l)
  } else {
    eq <- substitute(italic(y) == a - b %.% italic(x)*","~~italic(r)^2~"="~r2,l)    
  }

  as.character(as.expression(eq));                 
}

L'utilisation changerait en:

p1 = p + geom_text(aes(x = 25, y = 300, label = lm_eqn(lm(y ~ x, df))), parse = TRUE)
Jayden
la source
17
Ça a l'air super! Mais je trace des geom_points sur plusieurs facettes, où le df diffère en fonction de la variable de facette. Comment je fais ça?
bshor
24
La solution de Jayden fonctionne assez bien, mais la police de caractère semble très moche. Je recommanderais de changer l'utilisation de ceci: p1 = p + annotate("text", x = 25, y = 300, label = lm_eqn(lm(y ~ x, df)), colour="black", size = 5, parse=TRUE)edit: cela résout également tous les problèmes que vous pourriez avoir avec les lettres apparaissant dans votre légende.
Jonas Raedle
1
@ Jonas, pour une raison que j'obtiens "cannot coerce class "lm" to a data.frame". Cette alternative fonctionne: df.labs <- data.frame(x = 25, y = 300, label = lm_eqn(df))et p <- p + geom_text(data = df.labs, aes(x = x, y = y, label = label), parse = TRUE)
PatrickT
1
@PatrickT - C'est le message d'erreur que vous obtiendriez si vous lm_eqn(lm(...))appeliez avec la solution de Ramnath. Vous avez probablement essayé celui-ci après avoir essayé celui-ci, mais vous avez oublié de vous assurer que vous aviez redéfinilm_eqn
Hamy
@PatrickT: pourriez-vous faire de votre réponse une réponse distincte? Je serais ravi de voter!
JelenaČuklina
11

J'adore vraiment la solution @Ramnath. Pour autoriser l'utilisation pour personnaliser la formule de régression (au lieu de la fixer comme y et x comme noms de variables littérales), et ajouter la valeur p dans l'impression également (comme l'a commenté @Jerry T), voici le mod:

lm_eqn <- function(df, y, x){
    formula = as.formula(sprintf('%s ~ %s', y, x))
    m <- lm(formula, data=df);
    # formating the values into a summary string to print out
    # ~ give some space, but equal size and comma need to be quoted
    eq <- substitute(italic(target) == a + b %.% italic(input)*","~~italic(r)^2~"="~r2*","~~p~"="~italic(pvalue), 
         list(target = y,
              input = x,
              a = format(as.vector(coef(m)[1]), digits = 2), 
              b = format(as.vector(coef(m)[2]), digits = 2), 
             r2 = format(summary(m)$r.squared, digits = 3),
             # getting the pvalue is painful
             pvalue = format(summary(m)$coefficients[2,'Pr(>|t|)'], digits=1)
            )
          )
    as.character(as.expression(eq));                 
}

geom_point() +
  ggrepel::geom_text_repel(label=rownames(mtcars)) +
  geom_text(x=3,y=300,label=lm_eqn(mtcars, 'hp','wt'),color='red',parse=T) +
  geom_smooth(method='lm')

entrez la description de l'image ici Malheureusement, cela ne fonctionne pas avec facet_wrap ou facet_grid.

XX
la source
Très soigné, je l'ai référencé ici . Une clarification - votre code manque-t-il ggplot(mtcars, aes(x = wt, y = mpg, group=cyl))+avant le geom_point ()? Une question semi-liée - si nous nous référons à hp et wt dans le aes()for ggplot, pouvons-nous alors les saisir pour les utiliser dans l'appel à lm_eqn, donc nous n'avons qu'à coder en un seul endroit? Je sais que nous pourrions configurer xvar = "hp"avant l'appel à ggplot () et utiliser xvar dans les deux emplacements pour remplacer hp , mais cela semble être inutile.
Mark Neal
9

Utilisation de ggpubr :

library(ggpubr)

# reproducible data
set.seed(1)
df <- data.frame(x = c(1:100))
df$y <- 2 + 3 * df$x + rnorm(100, sd = 40)

# By default showing Pearson R
ggscatter(df, x = "x", y = "y", add = "reg.line") +
  stat_cor(label.y = 300) +
  stat_regline_equation(label.y = 280)

entrez la description de l'image ici

# Use R2 instead of R
ggscatter(df, x = "x", y = "y", add = "reg.line") +
  stat_cor(label.y = 300, 
           aes(label = paste(..rr.label.., ..p.label.., sep = "~`,`~"))) +
  stat_regline_equation(label.y = 280)

## compare R2 with accepted answer
# m <- lm(y ~ x, df)
# round(summary(m)$r.squared, 2)
# [1] 0.85

entrez la description de l'image ici

zx8754
la source
Avez-vous vu une manière programmatique soignée de spécifier un nombre label.y?
Mark Neal
@MarkNeal obtient peut-être le maximum de y puis multiplie par 0,8. label.y = max(df$y) * 0.8
zx8754
1
@MarkNeal bons points, peut-être soumettre le problème comme demande de fonctionnalité à GitHub ggpubr.
zx8754
1
Problème sur la localisation automatique soumis ici
Mark Neal
1
@ zx8754, dans votre intrigue, il est affiché rho et non R², un moyen facile de montrer R²?
matmar
6

Voici le code le plus simple pour tout le monde

Remarque: montrant Rho de Pearson et non R ^ 2.

library(ggplot2)
library(ggpubr)

df <- data.frame(x = c(1:100)
df$y <- 2 + 3 * df$x + rnorm(100, sd = 40)
p <- ggplot(data = df, aes(x = x, y = y)) +
        geom_smooth(method = "lm", se=FALSE, color="black", formula = y ~ x) +
        geom_point()+
        stat_cor(label.y = 35)+ #this means at 35th unit in the y axis, the r squared and p value will be shown
        stat_regline_equation(label.y = 30) #this means at 30th unit regresion line equation will be shown

p

Un exemple avec mon propre ensemble de données

Sork-kal
la source
Même problème que ci-dessus, dans votre tracé, il est indiqué rho et non R²!
matmar
3

Inspirée par le style d'équation fourni dans cette réponse , une approche plus générique (plus d'un prédicteur + sortie latex en option) peut être:

print_equation= function(model, latex= FALSE, ...){
    dots <- list(...)
    cc= model$coefficients
    var_sign= as.character(sign(cc[-1]))%>%gsub("1","",.)%>%gsub("-"," - ",.)
    var_sign[var_sign==""]= ' + '

    f_args_abs= f_args= dots
    f_args$x= cc
    f_args_abs$x= abs(cc)
    cc_= do.call(format, args= f_args)
    cc_abs= do.call(format, args= f_args_abs)
    pred_vars=
        cc_abs%>%
        paste(., x_vars, sep= star)%>%
        paste(var_sign,.)%>%paste(., collapse= "")

    if(latex){
        star= " \\cdot "
        y_var= strsplit(as.character(model$call$formula), "~")[[2]]%>%
            paste0("\\hat{",.,"_{i}}")
        x_vars= names(cc_)[-1]%>%paste0(.,"_{i}")
    }else{
        star= " * "
        y_var= strsplit(as.character(model$call$formula), "~")[[2]]        
        x_vars= names(cc_)[-1]
    }

    equ= paste(y_var,"=",cc_[1],pred_vars)
    if(latex){
        equ= paste0(equ," + \\hat{\\varepsilon_{i}} \\quad where \\quad \\varepsilon \\sim \\mathcal{N}(0,",
                    summary(MetamodelKdifEryth)$sigma,")")%>%paste0("$",.,"$")
    }
    cat(equ)
}

L' modelargument attend un lmobjet, l' latexargument est un booléen pour demander un caractère simple ou une équation au format latex, et l' ...argument transmet ses valeurs à laformat fonction.

J'ai également ajouté une option pour le produire sous forme de latex afin que vous puissiez utiliser cette fonction dans un rmarkdown comme celui-ci:


```{r echo=FALSE, results='asis'}
print_equation(model = lm_mod, latex = TRUE)
```

Maintenant, je l'utilise:

df <- data.frame(x = c(1:100))
df$y <- 2 + 3 * df$x + rnorm(100, sd = 40)
df$z <- 8 + 3 * df$x + rnorm(100, sd = 40)
lm_mod= lm(y~x+z, data = df)

print_equation(model = lm_mod, latex = FALSE)

Ce code donne: y = 11.3382963933174 + 2.5893419 * x + 0.1002227 * z

Et si nous demandons une équation en latex, en arrondissant les paramètres à 3 chiffres:

print_equation(model = lm_mod, latex = TRUE, digits= 3)

Cela donne: équation de latex

rvezy
la source
0

J'ai un doute, comment mettre une statistique significative de t.test pour bheta dans l'équation, en utilisant ggpmisc::stat_poly_eq()?

ex: expression(hat(Y)== 0000*"**"+0000*"x"*"*"-0000*"x"^2*"**"~~~~"R"^2*":"~~0.000)

Jean Karlos
la source