是否有一种巧妙的方法可以使用 geom_quantile() 中的方程和其他统计数据来标记 ggplot 图?

问题描述

我想以与处理 import subprocess for i in range(n): cmd = "cp methods.py driver.py subdir"+str(i) p = subprocess.Popen(cmd,shell=True) p.wait() #once every subdirectory has driver.py and methods.py,start running these codes for i in range(n): cmd = "cd subdir" + str(i) +" && python driver.py" p = subprocess.Popen(cmd,shell=True) p.wait() 拟合线性回归类似的方式包含来自 geom_quantile() 拟合线的相关统计数据(我之前使用过 ggpmisc awesome)。例如,这段代码:

geom_smooth(method="lm")

生成这个:

ggplot with linear ols equation

对于分位数回归,您可以将 # quantile regression example with ggpmisc equation # basic quantile code from here: # https://ggplot2.tidyverse.org/reference/geom_quantile.html library(tidyverse) library(ggpmisc) # see ggpmisc vignette for stat_poly_eq() code below: # https://cran.r-project.org/web/packages/ggpmisc/vignettes/user-guide.html#stat_poly_eq my_formula <- y ~ x #my_formula <- y ~ poly(x,3,raw = TRUE) # linear ols regression with equation labelled m <- ggplot(mpg,aes(displ,1 / hwy)) + geom_point() m + geom_smooth(method = "lm",formula = my_formula) + stat_poly_eq(aes(label = paste(stat(eq.label),"*\" with \"*",stat(rr.label),"*\",\"*",stat(f.value.label),and \"*",stat(p.value.label),"*\".\"",sep = "")),formula = my_formula,parse = TRUE,size = 3) 替换为 geom_smooth() 并绘制一条可爱的分位数回归线(在本例中为中位数):

geom_quantile()

geom_quantile plot

您将如何将汇总统计信息输出到标签中,或者在旅途中重新创建它们? (即除了在调用 ggplot 之前进行回归然后将其传入然后进行注释(例如类似于为线性回归所做的 herehere

解决方法

@mark-neal stat_fit_glance()quantreg::rq() 一起工作。然而,使用 stat_fit_glance() 涉及更多。此统计信息“不知道”对 glance() 的期望,因此必须手动组装 label

人们需要知道什么是可用的。可以在 ggplot 之外运行拟合模型并使用 glance() 找出它返回的列,或者可以在包 'gginnards' 的帮助下在 ggplot 中执行此操作。我将展示这个替代方案,从上面的代码示例继续。

library(gginnards)

m + 
  geom_quantile(quantiles = 0.5) +
  stat_fit_glance(method = "rq",method.args = list(formula = y ~ x),geom = "debug")

geom_debug() 默认仅将其输入打印到 R 控制台,其输入是统计信息返回的内容。

# A tibble: 1 x 11
   npcx  npcy   tau logLik    AIC    BIC df.residual     x      y PANEL group
  <dbl> <dbl> <dbl>  <dbl>  <dbl>  <dbl>       <int> <dbl>  <dbl> <fct> <int>
1    NA    NA   0.5   816. -1628. -1621.         232  1.87 0.0803 1        -1

我们可以使用 after_stat() 访问每一列(早期的化身是 stat() 并包含名称 ...。我们需要使用 {{1} 的编码符号进行格式化}}. 如果像这种情况我们组装一个需要解析成表达式的字符串,那么还需要sprintf()

parse = TRUE

这个例子产生了下图。 enter image description here

对于 m + geom_quantile(quantiles = 0.5) + stat_fit_glance(method = "rq",mapping = aes(label = sprintf('italic(tau)~"="~%.2f~~AIC~"="~%.3g~~BIC~"="~%.3g',after_stat(tau),after_stat(AIC),after_stat(BIC))),parse = TRUE) ,同样的方法应该有效。但是,在 'ggpmisc' (= 0.3.8) 中修复,现在在 CRAN 中。

以下示例仅适用于“ggpmisc”(>= 0.3.8)

剩下的问题是 stat_fit_tidy()tibble 返回的 glance() 是否包含想要添加到情节中的信息,{{1} 似乎并非如此}},至少在默认情况下。但是,tidy() 有一个参数 tidy.qr(),用于确定 tidy.rq() 中返回的值。修改后的 se.type 接受传递给 tibble 的命名参数,从而实现以下功能。

stat_fit_tidy()

这个例子产生了下图。

enter image description here

定义一个新的 stat tidy() 会更简单:

m + 
  geom_quantile(quantiles = 0.5) +
  stat_fit_tidy(method = "rq",tidy.args = list(se.type = "nid"),mapping = aes(label = sprintf('y~"="~%.3g~+~%.3g~x*",with "*italic(P)~"="~%.3f',after_stat(Intercept_estimate),after_stat(x_estimate),after_stat(x_p.value))),parse = TRUE)

答案变成:

stat_rq_eq()
,

请将此视为 Pedro 出色答案的附录,他在其中完成了大部分繁重工作 - 这增加了一些演示调整(颜色和线型)和代码以简化多个分位数,生成如下图:

quantile plot

library(tidyverse)
library(ggpmisc) #ensure version 0.3.8 or greater
library(quantreg)
library(generics)

my_formula <- y ~ x
#my_formula <- y ~ poly(x,3,raw = TRUE)

# base plot
m <- ggplot(mpg,aes(displ,1 / hwy)) +
  geom_point()

# function for labelling
# Doesn't neatly handle P values (e.g return "P<0.001 where appropriate)

stat_rq_eqn <- function(formula = y ~ x,tau = 0.5,colour = "red",label.y = 0.9,...) {
  stat_fit_tidy(method = "rq",method.args = list(formula = formula,tau = tau),mapping = aes(label = sprintf('italic(tau)~"="~%.3f~";"~y~"="~%.3g~+~%.3g~x*",with "~italic(P)~"="~%.3g',after_stat(x_tau),parse = TRUE,colour = colour,label.y = label.y,...)
}

# This works,though with double entry of plot specs
# custom colours and linetype
# https://stackoverflow.com/a/44383810/4927395
# https://stackoverflow.com/a/64518380/4927395


m + 
  geom_quantile(quantiles = c(0.1,0.5,0.9),aes(colour = as.factor(..quantile..),linetype = as.factor(..quantile..))
                )+
  scale_color_manual(values = c("red","purple","darkgreen"))+
  scale_linetype_manual(values = c("dotted","dashed","solid"))+
  stat_rq_eqn(tau = 0.1,label.y = 0.9)+
  stat_rq_eqn(tau = 0.5,colour = "purple",label.y = 0.95)+
  stat_rq_eqn(tau = 0.9,colour = "darkgreen",label.y = 1.0)+
  theme(legend.position = "none") # suppress legend


# not a good habit to have double entry above
# modified with reference to tibble for plot specs,# though still a stat_rq_eqn call for each quantile and manual vertical placement
# https://www.r-bloggers.com/2019/06/curly-curly-the-successor-of-bang-bang/

my_tau = c(0.1,0.9)
my_colours = c("red","darkgreen")
my_linetype = c("dotted","solid")

quantile_plot_specs <- tibble(my_tau,my_colours,my_linetype)

m + 
  geom_quantile(quantiles = {{quantile_plot_specs$my_tau}},linetype = as.factor(..quantile..))
  )+
  scale_color_manual(values = {{quantile_plot_specs$my_colours}})+
  scale_linetype_manual(values = {{quantile_plot_specs$my_linetype}})+
  stat_rq_eqn(tau = {{quantile_plot_specs$my_tau[1]}},colour = {{quantile_plot_specs$my_colours[1]}},label.y = 0.9)+
  stat_rq_eqn(tau = {{quantile_plot_specs$my_tau[2]}},colour = {{quantile_plot_specs$my_colours[2]}},label.y = 0.95)+
  stat_rq_eqn(tau = {{quantile_plot_specs$my_tau[3]}},colour = {{quantile_plot_specs$my_colours[3]}},label.y = 1.0)+
  theme(legend.position = "none")

相关问答

错误1:Request method ‘DELETE‘ not supported 错误还原:...
错误1:启动docker镜像时报错:Error response from daemon:...
错误1:private field ‘xxx‘ is never assigned 按Alt...
报错如下,通过源不能下载,最后警告pip需升级版本 Requirem...