文章/答案/技术大牛

发布

社区首页 >专栏 >基于 R 语言的科研论文绘图技巧详解（2）

基于 R 语言的科研论文绘图技巧详解（2）

庄闪闪

发布于 2022-05-24 07:44:05

1K00

代码可运行

文章被收录于专栏：庄闪闪的R语言手册庄闪闪的R语言手册

运行总次数：0

代码可运行

简介

在查阅文献的过程中，看到了几幅非常不错的出版图，今天就跟着小编一起学习下，他们是怎么使用 R 绘制出来的。

今天主要介绍 第二幅图（B) ，直观来看是由两幅图所构成的。

绘制带误差项的柱状图并添加密度函数线。
简单的曲线图并添加公式。

之后，将两幅图合并。值得注意的是：x 轴数值使用不同图形进行描绘（小编不是很懂，作者想表达什么，不过这种技巧是第一次见，可以学习学习）。

第一幅图的介绍可见：基于 R 语言的科研论文绘图技巧详解（1）。后面几幅图会一一介绍，读者在学习过程中，可以将内部学到的知识点应用到自己的图形绘制中。推文已经将主要知识点进行罗列，更有利于读者学习和查阅。

那我们来看看，作者是怎么实现这个功能的吧，本文知识点较多，大家耐心学习，建议自己实践。对应代码、数据可在 GitHub - marco-meer/scifig_plot_examples_R: Scientific publication figure plotting examples with R[1] 中找到。

主要知识点

学会定义密度函数，并在图形中将其添加；
学会设置自定义主题，简化代码，统一主题，方便绘制其他图形使用；
学会添加子图、给坐标轴添加修饰图形；
学会添加带有特殊符号的公式。

绘图

加载包

首先加载一些需要使用到的包。

library(ggplot2) # Grammar of graphics
library(cowplot) # Arranging multiple plots into a grid
library(png)     # Load JPEG, PNG and TIFF format 
library(scales)  # Generic plot scaling methods
library(viridis) # Default color maps from 'matplotlib'
library(grid)    # A rewrite of the graphics layout capabilities
library(magick)  # graphics and image processing
library(rsvg)    # Render svg image into a high quality bitmap
library(ggforce) # Collection of additional ggplot stats + geoms

设置主题

接下来，为了方便起见，作者在绘图前设置好了主题，并将该函数命名为 my_theme。

这一部分在第一篇推文给出，代码将在文末中完整代码给出。需要数据的朋友可以从 GitHub 上下载（或者文末链接直达）。

手动修改大部分面板，具体可以参考本篇文章[2]。或者观看我在 B 站发布的《R 语言可视化教程》，里面也有一些简单主题设置介绍。

导入数据

首先使用 read.csv() 导入数据，使用data.frame()转化数据成 ggplot 所需要的格式。

data_B = read.csv("./data_B.csv")
data_B = data.frame("n"=c(data_B$n,data_B$n),
                    "sample"=c(rep("apical side",nrow(data_B)),
                               rep("basal side",nrow(data_B))),
                    "fraction"=c(data_B$fraction1,data_B$fraction2),
                    "err"=c(data_B$err1,data_B$err2)
)
head(data_B)
#  n      sample fraction   err
# 1 3 apical side    0.010 0.002
# 2 4 apical side    0.050 0.010
# 3 5 apical side    0.240 0.040
# 4 6 apical side    0.450 0.060
# 5 7 apical side    0.190 0.034
# 6 8 apical side    0.045 0.009

定义密度函数

这里作者直接定义了对数正态密度函数，并确定了对应参数的值（mu=log(6),sigma= 0.14)。

sigma = 0.14
n = seq(3,10,1)
logn_dist <- function(n) exp(-(log(n)-log(6))^2/(2*sigma^2))/(sqrt(2*pi)*sigma*n) #函数里面就是密度函数的公式，可以在数理统计相关书籍中找到。

小编表示疑问：参数应该得通过参数估计的方式得到吧，但是这里作者直接给定了。并且函数设定来看，均值直接写在函数内部了，并没有赋值为 mu=log(6)。虽然结果相同，但是可读性不强。

绘图步骤详解

由于代码复杂，知识点较多，为了读者更好理解代码逻辑和含义，小编将其分布讲解。最后再将完整代码放到本节末。

绘制子图一（条形、密度函数、误差图）

使用 geom_bar() 绘制条形图，注意这里使用了 position=position_dodge() 将小类并列放置，具体细节可以参考：《R语言教程》[3]；使用 geom_errorbar() 添加误差项；使用 stat_function() 将对数正态的密度函数加入图中（当然也可以使用 geom_line()）。

 ggplot(data=data_B,(aes(x=n, # aes: aesthetics
                          y=fraction,
                          fill=sample)))+
  geom_bar(stat = "identity",
           position=position_dodge()) + # dodge overlapping objects side-to-side
  stat_function(fun=logn_dist, 
                geom="line",
                linetype="solid",
                aes(x=n,
                    y=logn_dist(n))) +
  geom_errorbar(aes(ymin=fraction-err, 
                    ymax=fraction+err), 
                width=.2, 
                position=position_dodge(.9))

微调主题

之后添加主题，使用先前设定好的主题函数 my_theme() 与其他细节调整。

  my_theme() +
  # some extra theme tweaking 
  theme(legend.position = c(0.18,0.95),
        legend.title = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.x = element_text(vjust=-0.5),
        axis.title.x = element_text(vjust=-0.5),
        legend.key.size = unit(0.4, "lines")
  ) +
  xlab(expression(paste("number of neighbors ",italic("n")))) +
  ylab("fraction of cells (%)")

这时候图中的第一个子图就得到了，接下来绘制第二个子图。

绘制子图二（曲线图）

方式一：定义曲线函数，然后通过 geom_line() 进行绘制。
方式二：加载图形（公式），通过 annotation_custom() 结合 rasterGrob()，将该图导入到子图中。

方式二看着比较麻烦，不过思路可以学习下，如果以后公式太难/复杂，可以使用这种方式展现。下面展示第一种结果，第二种结果可在官网代码中找到。

其他代码流程和上一个子图类似，这里不做过多介绍。

注意：关于特殊符号的输入，可以使用 expression() 函数，可以参考这个教程[4]。

定义曲线函数
inset_curve <- function(n) 180*(n-2)/n 

inset <- ggplot() + 
  geom_line(aes(x=n,y=inset_curve(n))) +
  annotate(geom="text",label="180*degree*(n-2)/n",x = 5.5,y=140,
           parse=TRUE,size=5) +
  scale_x_continuous(expand = c(0, 0),
                     breaks=c(seq(3,10,1)),
                     limits = c(3,10)) +
  scale_y_continuous(expand = c(0, 0),
                     breaks=c(seq(60,150,20)),
                     limits = c(60,150)
  ) +
  my_theme() +
  xlab(expression(italic(n)))+ 
  ylab(expression(interior~angle~italic(theta[n])(degree))) +
  # some extra theme tweaking 
  theme( 
    panel.background = element_blank(),
    plot.background = element_blank(),
    axis.ticks.length = unit(-0.1, "lines"), # make ticks a bit shorter
    axis.title.x = element_text(size = 0.5*base_size, 
                                color = "black", 
                                margin = margin(t = -4)), 
    # top (t), right (r), bottom (b), left (l) 
    axis.title.y = element_text(size = 0.5*base_size, 
                                color = "black", angle = 90,
                                margin = margin(r = -4)),  
    axis.text.x = element_text(size = base_size*0.5, 
                               color = "black", 
                               lineheight = 0.9,
                               margin=unit(c(0.1,0.1,0.1,0.1), "cm")),
    axis.text.y = element_text(size = base_size*0.5, color = "black", 
                               lineheight = 0.9,
                               margin=unit(c(0.1,0.1,0.1,0.1), "cm"))
  ) 

inset

修改坐标轴值样式

使用 ggforce 包中的 geom_regon() 函数，不同图形主要是在内部参数 size 进行设置。

polygons <- 
  ggplot() +
  ggforce::geom_regon(aes(x0=c(seq(3,10,1)),
                          y0=rep(-0.3,8),
                          sides = c(seq(3,10,1)), 
                          angle=0, r=0.3),
                      fill=NA,
                      color="black") +
  theme_void() + coord_fixed()
polygons

完整版

最后将前面的内容整合在一起。还是使用 annotation_custom()，感觉这个函数无敌！

panel_B <- 
  panel_B + 
  annotation_custom(grob = panel_B_inset, 
                    xmin=6.6, xmax=10.1, ymin=0.1, ymax=0.6) +
  annotation_custom(grob = ggplotGrob(polygons), 
                    xmin = 2.36, xmax = 10.67, ymin = -0.172, ymax = 0.115) 
panel_B

完整版

完整代码

library(ggplot2) # Grammar of graphics
library(cowplot) # Arranging multiple plots into a grid
library(png)     # Load JPEG, PNG and TIFF format 
library(scales)  # Generic plot scaling methods
library(viridis) # Default color maps from 'matplotlib'
library(grid)    # A rewrite of the graphics layout capabilities
library(magick)  # graphics and image processing
library(rsvg)    # Render svg image into a high quality bitmap
library(ggforce) # Collection of additional ggplot stats + geoms

base_size = 12 

# Panel B ----
my_theme <-  function() {
  theme(
    aspect.ratio = 1,
    axis.line =element_line(colour = "black"),  
    
    # shift axis text closer to axis bc ticks are facing inwards
    axis.text.x = element_text(size = base_size*0.8, color = "black", 
                               lineheight = 0.9,
                               margin=unit(c(0.3,0.3,0.3,0.3), "cm")), 
    axis.text.y = element_text(size = base_size*0.8, color = "black", 
                               lineheight = 0.9,
                               margin=unit(c(0.3,0.3,0.3,0.3), "cm")),  
    
    axis.ticks = element_line(color = "black", size  =  0.2),  
    axis.title.x = element_text(size = base_size, 
                                color = "black", 
                                margin = margin(t = -5)), 
    # t (top), r (right), b (bottom), l (left)
    axis.title.y = element_text(size = base_size, 
                                color = "black", angle = 90,
                                margin = margin(r = -5)),  
    axis.ticks.length = unit(-0.3, "lines"),  
    legend.background = element_rect(color = NA, 
                                     fill = NA), 
    legend.key = element_rect(color = "black",  
                              fill = "white"),  
    legend.key.size = unit(0.5, "lines"),  
    legend.key.height =NULL,  
    legend.key.width = NULL,      
    legend.text = element_text(size = 0.6*base_size, 
                               color = "black"),  
    legend.title = element_text(size = 0.6*base_size, 
                                face = "bold", 
                                hjust = 0, 
                                color = "black"),  
    legend.text.align = NULL,  
    legend.title.align = NULL,  
    legend.direction = "vertical",  
    legend.box = NULL, 
    panel.background = element_rect(fill = "white", 
                                    color  =  NA),  
    panel.grid.major = element_blank(),
    panel.grid.minor = element_blank(),
    plot.title = element_text(size = base_size, 
                              color = "black"), 
  ) 
  
  
}


data_B = read.csv("./data_B.csv")

# format data to ggplot's liking
data_B = data.frame("n"=c(data_B$n,data_B$n),
                    "sample"=c(rep("apical side",nrow(data_B)),
                               rep("basal side",nrow(data_B))),
                    "fraction"=c(data_B$fraction1,data_B$fraction2),
                    "err"=c(data_B$err1,data_B$err2)
)
head(data_B)

# define lognormal distribution to be called via ggplot2::stat_function
sigma = 0.14
n = seq(3,10,1)
logn_dist <- function(n) exp(-(log(n)-log(6))^2/(2*sigma^2))/(sqrt(2*pi)*sigma*n)



panel_B <- 
  ggplot(data=data_B,(aes(x=n, # aes: aesthetics
                          y=fraction,
                          fill=sample)))+
  geom_bar(stat = "identity",
           position=position_dodge()) + # dodge overlapping objects side-to-side
  stat_function(fun=logn_dist, 
                geom="line",
                linetype="solid",
                aes(x=n,
                    y=logn_dist(n))) +
  geom_errorbar(aes(ymin=fraction-err, 
                    ymax=fraction+err), 
                width=.2, 
                position=position_dodge(.9)) +
  scale_fill_manual(values=c('#f2a340','#998fc2'))+
  scale_x_continuous(expand = c(0, 0), # prevent gap between origin and first tick
                     breaks=c(seq(from =  min(data_B$n), 
                                  to = max(data_B$n),
                                  by = 1)),
                     limits = c(min(data_B$n),max(data_B$n))) +
  scale_y_continuous(expand = c(0, 0),
                     breaks=c(seq(0,0.6,0.1)),
                     limits = c(0,0.6)
  )+
  my_theme() +
  # some extra theme tweaking 
  theme(legend.position = c(0.18,0.95),
        legend.title = element_blank(),
        axis.ticks.x=element_blank(),
        axis.text.x = element_text(vjust=-0.5),
        axis.title.x = element_text(vjust=-0.5),
        legend.key.size = unit(0.4, "lines")
  ) +
  xlab(expression(paste("number of neighbors ",italic("n")))) +
  ylab("fraction of cells (%)") 

panel_B


inset_curve <- function(n) 180*(n-2)/n 
# The equation is added as an image bc the annotate() function output looks too ugly
# inset_equation <- image_read_svg("./panel_b_inset_equation.svg",width = 583,height = 240)

# now comes the inset plot

inset <- ggplot() + 
  geom_line(aes(x=n,y=inset_curve(n))) +
  annotate(geom="text",label="180*degree*(n-2)/n",x = 5.5,y=140,
           parse=TRUE,size=3) +
  # annotation_custom(rasterGrob(image = inset_equation, 
  #                              x=0.5,
  #                              y=0.8,
  #                              width = unit(0.3125,"npc"),
  #                              height = unit(0.126,"npc")), 
  #                   -Inf, Inf, -Inf, Inf) +
  scale_x_continuous(expand = c(0, 0),
                     breaks=c(seq(3,10,1)),
                     limits = c(3,10)) +
  scale_y_continuous(expand = c(0, 0),
                     breaks=c(seq(60,150,20)),
                     limits = c(60,150)
  ) +
  my_theme() +
  xlab(expression(italic(n)))+ 
  # expression doc: https://stat.ethz.ch/R-manual/R-devel/library/grDevices/html/plotmath.html
  ylab(expression(interior~angle~italic(theta[n])(degree))) +
  # some extra theme tweaking 
  theme( 
    panel.background = element_blank(),
    plot.background = element_blank(),
    axis.ticks.length = unit(-0.1, "lines"), # make ticks a bit shorter
    axis.title.x = element_text(size = 0.5*base_size, 
                                color = "black", 
                                margin = margin(t = -4)), 
    # top (t), right (r), bottom (b), left (l) 
    axis.title.y = element_text(size = 0.5*base_size, 
                                color = "black", angle = 90,
                                margin = margin(r = -4)),  
    axis.text.x = element_text(size = base_size*0.5, 
                               color = "black", 
                               lineheight = 0.9,
                               margin=unit(c(0.1,0.1,0.1,0.1), "cm")),
    axis.text.y = element_text(size = base_size*0.5, color = "black", 
                               lineheight = 0.9,
                               margin=unit(c(0.1,0.1,0.1,0.1), "cm"))
  ) 

inset
# create grob (grid graphical object) from the inset plot
panel_B_inset <- ggplotGrob(x = inset) 
panel_B_inset

# now for some regular polygon drawing
polygons <- 
  ggplot() +
  ggforce::geom_regon(aes(x0=c(seq(3,10,1)),
                          y0=rep(-0.3,8),
                          sides = c(seq(3,10,1)), 
                          angle=0, r=0.3),
                      fill=NA,
                      color="black") +
  theme_void() + coord_fixed()
polygons


# finally, we can combine everything to plot panel B

panel_B <- 
  panel_B + 
  annotation_custom(grob = panel_B_inset, 
                    xmin=6.6, xmax=10.1, ymin=0.1, ymax=0.6) +
  annotation_custom(grob = ggplotGrob(polygons), 
                    xmin = 2.36, xmax = 10.67, ymin = -0.172, ymax = 0.115) 
panel_B

小编有话说

本文主要学到的知识点如下：

自定义密度函数，并使用 stat_function() 在图形中将其添加；
设置自定义主题(my_theme)，简化代码，统一主题，方便绘制其他图形使用；
使用 annotation_custom() 添加子图；
使用 ggforce 包中的 geom_regon() 函数绘制修饰图形；
使用 expression() 函数添加带有特殊符号的公式。
如果觉得文章对你有一丝丝用处，帮忙点赞，在看，分享吧～

参考资料

[1]

GitHub - marco-meer/scifig_plot_examples_R: Scientific publication figure plotting examples with R: https://github.com/marco-meer/scifig_plot_examples_R

[2]

文章: https://ggplot2.tidyverse.org/reference/theme.html

[3]

《R语言教程》: https://www.math.pku.edu.cn/teachers/lidf/docs/Rbook/html/_Rbook/ggplotvis.html

[4]

教程: https://stat.ethz.ch/R-manual/R-devel/library/grDevices/html/plotmath.html

本文参与腾讯云自媒体同步曝光计划，分享自微信公众号。

原始发表：2022-04-26，如有侵权请联系 cloudcommunity@tencent.com 删除

https