Adicione a equação da linha de regressão e R ^ 2 no gráfico

228

Gostaria de saber como adicionar a equação da linha de regressão e R ^ 2 no ggplot. Meu código é:

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

Qualquer ajuda será muito apreciada.

MYaseen208
fonte
1
Para gráficos de treliça , consulte latticeExtra::lmlineq().
Josh O'Brien

Respostas:

234

Aqui está uma solução

# 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)

EDITAR. Eu descobri a fonte de onde eu escolhi esse código. Aqui está o link para a postagem original nos grupos do Google ggplot2

Resultado

Ramnath
fonte
1
O comentário de @ JonasRaedle sobre como obter textos com melhor aparência annotateestava correto na minha máquina.
IRTFM
2
Isso não se parece nada com a saída publicada em minha máquina, onde o rótulo é substituído tantas vezes quanto os dados são chamados, resultando em um texto grosso e embaçado. Passando os rótulos para um data.frame primeiros trabalhos (ver a minha sugestão em um comentário abaixo.
PatrickT
@PatrickT: remova o aes(e o correspondente ). aesé para mapear variáveis ​​do quadro de dados para variáveis ​​visuais - isso não é necessário aqui, pois há apenas uma instância, para que você possa colocar tudo na geom_textchamada principal . Vou editar isso na resposta.
naught101
O problema com esta solução parece ser que, se o conjunto de dados for maior (o meu era 370000 observações), a função parece falhar. Eu recomendaria a solução da @kdauria que faz o mesmo, mas muito, muito mais rápido.
Benjamin Benjamin
3
para quem deseja valores de rp em vez de R2 e equação: eq <- substitute (itálico (r) ~ "=" ~ rvalue * "," ~ itálico (p) ~ "=" ~ pvalor, lista (rvalue = sprintf ("% .2f", sinal (coef (m) [2]) * sqrt (resumo (m) $ r.squared)), pvalor = formato (resumo (m) $ coeficientes [2,4], dígitos = 2 )))
Jerry T
135

Eu incluí uma estatística stat_poly_eq()no meu pacote ggpmiscque permite esta resposta:

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

insira a descrição da imagem aqui

Essa estatística funciona com qualquer polinômio sem termos ausentes e, esperançosamente, possui flexibilidade suficiente para ser geralmente útil. As etiquetas R ^ 2 ou R ^ 2 ajustadas podem ser usadas com qualquer fórmula de modelo equipada com lm ​​(). Sendo uma estatística ggplot, ela se comporta conforme o esperado, tanto em grupos quanto em facetas.

O pacote 'ggpmisc' está disponível no CRAN.

A versão 0.2.6 foi aceita no CRAN.

Ele aborda os comentários de @shabbychef e @ MYaseen208.

@ MYaseen208 mostra como adicionar um chapéu .

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

insira a descrição da imagem aqui

@shabbychef Agora é possível combinar as variáveis ​​na equação com as utilizadas para os rótulos dos eixos. Para substituir x por dizer z e y por h , usaria:

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

insira a descrição da imagem aqui

Sendo essas expressões R analisadas normais, as letras gregas também podem agora ser usadas nos lhs e rhs da equação.

[08-03-2017] @elarry Edite para abordar com mais precisão a pergunta original, mostrando como adicionar uma vírgula entre os rótulos de equação e 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

insira a descrição da imagem aqui

[2019-10-20] @ helen.h A seguir, exemplos de uso de stat_poly_eq()com agrupamento.

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

insira a descrição da imagem aqui

insira a descrição da imagem aqui

[2020-01-21] @Herman Pode ser um pouco contra-intuitivo à primeira vista, mas para obter uma única equação ao usar o agrupamento, é necessário seguir a gramática dos gráficos. Restrinja o mapeamento que cria o agrupamento para camadas individuais (mostradas abaixo) ou mantenha o mapeamento padrão e substitua-o por um valor constante na camada em que você não deseja o agrupamento (por exemplo colour = "black").

Continuando do exemplo anterior.

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

insira a descrição da imagem aqui

[2020-01-22] Para completar, um exemplo com facetas, demonstrando que também nesse caso as expectativas da gramática dos gráficos são cumpridas.

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

insira a descrição da imagem aqui

Pedro Aphalo
fonte
1
Deve-se notar que o xe yna fórmula se referem aos dados xe ynas camadas do gráfico, e não necessariamente aos que estão no escopo no momento my.formula. Assim, a fórmula deve sempre usar variáveis ​​x e y?
shabbychef
É bem verdade que xe se yreferem às variáveis ​​que estão mapeadas para essa estética. Essa é a expectativa também para geom_smooth () e como a gramática dos gráficos funciona. Poderia ter sido mais claro usar nomes diferentes dentro do quadro de dados, mas apenas os mantive como na pergunta original.
Pedro Aphalo 6/02/16
Será possível na próxima versão do ggpmisc. Obrigado pela sugestão!
Pedro Aphalo 25/02
3
Bom ponto @elarry! Isso está relacionado a como a função parse () de R funciona. Por tentativa e erro, descobri que aes(label = paste(..eq.label.., ..rr.label.., sep = "*plain(\",\")~"))faz o trabalho.
Pedro Aphalo
1
@HermanToothrot Normalmente, o R2 é preferido para uma regressão, portanto, não há uma etiqueta r predefinida nos dados retornados por stat_poly_eq(). Você também pode usar o stat_fit_glance()pacote 'ggpmisc', que retorna R2 como um valor numérico. Veja exemplos na página de ajuda e substitua stat(r.squared)por sqrt(stat(r.squared)).
Pedro Aphalo 14/03
99

Alterei algumas linhas da fonte stat_smoothe das funções relacionadas para criar uma nova função que adiciona a equação de ajuste e o valor de R ao quadrado. Isso também funcionará em gráficos de faceta!

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)

insira a descrição da imagem aqui

Eu usei o código na resposta de @ Ramnath para formatar a equação. A stat_smooth_funcfunção não é muito robusta, mas não deve ser difícil brincar com ela.

https://gist.github.com/kdauria/524eade46135f6348140 . Tente atualizar ggplot2se você receber um erro.

kdauria
fonte
2
Muito Obrigado. Este não funciona apenas para facetas, mas também para grupos. Eu acho muito útil para regressões por partes, por exemplo stat_smooth_func(mapping=aes(group=cut(x.val,c(-70,-20,0,20,50,130))),geom="text",method="lm",hjust=0,parse=TRUE), em combinação com EvaluateSmooths de stackoverflow.com/questions/19735149/…
Julian
1
@aelwan, altere estas linhas: gist.github.com/kdauria/… como desejar. Em seguida, sourceo arquivo inteiro em seu script.
Kdauria
1
@kdauria E se eu tiver várias equações em cada facet_wraps e diferentes y_values ​​em cada facet_wrap. Alguma sugestão de como fixar as posições das equações? Eu tentei várias opções de hjust, vjust eo ângulo usando este exemplo dropbox.com/s/9lk9lug2nwgno2l/R2_facet_wrap.docx?dl=0 mas eu não poderia trazer todas as equações no mesmo nível em cada um dos facet_wrap
brilhante
3
@aelwan, a posição da equação é determinada por estas linhas: gist.github.com/kdauria/… . Eu fiz xpose yposargumentos da função no Gist. Então, se você quiser que todas as equações se sobreponham, basta definir xpose ypos. Caso contrário, xpose ypossão calculados a partir dos dados. Se você deseja algo mais sofisticado, não deve ser muito difícil adicionar alguma lógica dentro da função. Por exemplo, talvez você possa escrever uma função para determinar qual parte do gráfico tem mais espaço vazio e colocar a função lá.
Kdauria # 10/16
6
Eu encontrei um erro com o source_gist: Erro nos arquivos r_ [[qual]]: tipo de subscrito inválido 'encerramento'. Veja este post para a solução: stackoverflow.com/questions/38345894/r-source-gist-not-working
Matifou
73

Modifiquei o post de Ramnath para: a) tornar mais genérico para que ele aceite um modelo linear como parâmetro, em vez do quadro de dados eb) exiba os negativos de maneira mais apropriada.

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));                 
}

O uso mudaria para:

p1 = p + geom_text(aes(x = 25, y = 300, label = lm_eqn(lm(y ~ x, df))), parse = TRUE)
Jayden
fonte
17
Isso parece ótimo! Mas estou plotando geom_points em várias facetas, onde o df difere com base na variável da faceta. Como faço isso?
12123
24
A solução de Jayden funciona muito bem, mas o tipo de letra parece muito feio. Eu recomendaria alterar o uso para: p1 = p + annotate("text", x = 25, y = 300, label = lm_eqn(lm(y ~ x, df)), colour="black", size = 5, parse=TRUE)edit: isso também resolve qualquer problema que você possa ter com as letras que aparecem na sua legenda.
Jonas Raedle
1
@ Jonas, por algum motivo eu estou recebendo "cannot coerce class "lm" to a data.frame". Esta alternativa funciona: df.labs <- data.frame(x = 25, y = 300, label = lm_eqn(df))e p <- p + geom_text(data = df.labs, aes(x = x, y = y, label = label), parse = TRUE)
PatrickT
1
@PatrickT - Essa é a mensagem de erro que você receberia se lm_eqn(lm(...))ligasse com a solução de Ramnath. Você provavelmente já tentou este depois de tentar que um, mas esqueceu-se de garantir que você redefiniulm_eqn
Hamy
@PatrickT: você poderia fazer da sua resposta uma resposta separada? Eu ficaria feliz em votar!
JelenaČuklina 2/11
11

realmente amo a solução @Ramnath. Para permitir o uso para personalizar a fórmula de regressão (em vez de fixos como yex como nomes de variáveis ​​literais) e também adicionar o valor p na impressão (como @Jerry T comentou), aqui está o 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')

insira a descrição da imagem aqui Infelizmente, isso não funciona com facet_wrap ou facet_grid.

XX
fonte
Muito arrumado, eu referenciei aqui . Um esclarecimento - o seu código está faltando ggplot(mtcars, aes(x = wt, y = mpg, group=cyl))+antes do geom_point ()? Uma pergunta semi-relacionada - se nos referirmos a hp e wt no aes()for ggplot, podemos pegá- los para usar na chamada para lm_eqn, então precisamos codificar apenas em um lugar? Eu sei que poderíamos configurar xvar = "hp"antes da chamada ggplot () e usar xvar nos dois locais para substituir o hp , mas isso parece que deve ser desnecessário.
Mark Neal
9

Usando 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)

insira a descrição da imagem aqui

# 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

insira a descrição da imagem aqui

zx8754
fonte
Você viu uma maneira programática clara de especificar um número label.y?
Mark Neal
@ MarkNeal talvez obtenha o máximo de y e multiplique por 0,8. label.y = max(df$y) * 0.8
zx8754 19/03
1
Os pontos positivos do @NarkNeal, talvez enviem um problema como solicitação de recurso no GitHub ggpubr.
zx8754 19/03
1
Problema na localização automática enviado aqui
Mark Neal
1
@ zx8754, no seu gráfico é mostrado rho e não R², alguma maneira fácil de mostrar R²?
matmar 27/04
6

Aqui está o código mais simples para todos

Nota: Mostrando Rho de Pearson e não 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

Um exemplo com meu próprio conjunto de dados

Sork-kal
fonte
Mesmo problema que o anterior, no seu gráfico é mostrado rho e não R²!
matmar 27/04
3

Inspirado no estilo de equação fornecido nesta resposta , uma abordagem mais genérica (mais de um preditor + saída de látex como opção) pode ser:

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)
}

O modelargumento espera um lmobjeto, o latexargumento é um booleano para solicitar um caractere simples ou uma equação em formato de látex, e o ...argumento passa seus valores para o parâmetroformat função.

Também adicionei uma opção para produzi-lo como látex, para que você possa usar esta função em um rmarkdown como este:


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

Agora usando:

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)

Este código produz: y = 11.3382963933174 + 2.5893419 * x + 0.1002227 * z

E se pedirmos uma equação do látex, arredondar os parâmetros para 3 dígitos:

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

Isso produz: equação do látex

rvezy
fonte
0

Tenho uma dúvida, como colocar uma estatística significativa de t.test para bheta na equação, usando ggpmisc::stat_poly_eq()?

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

Jean Karlos
fonte