简介
在查阅文献的过程中，看到了几幅非常不错的出版图，今天就跟着小编一起学习下，他们是怎么使用 R 绘制出来的。
今天主要介绍 第四幅图（D) —— 实现双 Y 轴，并且添加坐标轴的微小刻度线。这个图在科研绘图中较为常用，例如：将算法的收敛情况和计算所耗时间同时绘制。
前三幅图的详细代码介绍可见：基于 R 语言的科研论文绘图技巧详解（3）基于 R 语言的科研论文绘图技巧详解（2）基于 R 语言的科研论文绘图技巧详解（1）。后面几幅图会一一介绍，读者在学习过程中，可以将内部学到的知识点应用到自己的图形绘制中。推文已经将主要知识点进行罗列，更有利于读者学习和查阅。
那我们来看看，作者是怎么实现这个功能的吧，本文知识点较多，大家耐心学习，建议自己实践。对应代码、数据可在 GitHub - marco-meer/scifig_plot_examples_R: Scientific publication figure plotting examples with R[1] 中找到。
主要知识点
- 实现双 Y 轴；
- 学会修改坐标轴为对数尺度；
- 添加坐标轴的微小刻度线。
绘图
加载包
首先加载一些需要使用到的包。
代码语言：javascript复制library(ggplot2) # Grammar of graphics
设置主题
接下来，为了方便起见，作者在绘图前设置好了主题，并将该函数命名为 my_theme。
这一部分在第一篇推文 基于 R 语言的科研论文绘图技巧详解（1）给出，代码将在文末中完整代码给出。
手动修改大部分面板，具体可以参考本篇文章[2]。或者观看我在 B 站发布的《R 语言可视化教程》，里面也有一些简单主题设置介绍。
导入数据
首先使用 read.csv() 导入数据，其中一个数据前几行如下所示。
data_D1 = read.csv("./data_D1.csv")
data_D2 = read.csv("./data_D2.csv")
# format data to ggplot's liking
data_D = data.frame("width"=c(data_D1$width,data_D2$width),
"unit"=c(rep("shear_stress",nrow(data_D1)),
rep("velocity",nrow(data_D2))),
"value"=c(data_D1$shear_stress,data_D2$velocity)
)
head(data_D)
# width unit value
# 1 20.0 shear_stress 0.00174000
# 2 6.0 shear_stress 0.01622000
# 3 2.0 shear_stress 0.16065000
# 4 1.0 shear_stress 0.65696000
# 5 4.0 shear_stress 0.02312000
# 6 0.6 velocity 0.01726544
这里得到的数据，一共有三列，两个数据集的值在 value 中，width 放了两个数据集各自的width， unit 为离散数据。
绘图步骤详解
关键在于如何构建双 Y 轴，下面来看看作者是怎么设置的吧。
绘制单轴
首先，处理下第一个线性图所需要的数据，一共是两列。
代码语言：javascript复制curve_D1 = data.frame(width=data_D1$width,
shear_stress=33.28/(pi*18*data_D1$width^2))
curve_D1
# width shear_stress
# 1 20 0.001471299
# 2 6 0.016347767
# 3 2 0.147129903
# 4 1 0.588519612
# 5 4 0.036782476
这里绘制，小编带大家一步步解释，尤其注意作者的思想。
顺便提一下，很多公众号只是给了搬运的代码，但是并没有解释其中的奥秘。这对于初学者而言，就很难理解了。小编这个系列就是带大家一起学习作者画图思路。学会融会贯通，用到自己的科研绘图中。
先简单绘制出线性图，可以看到：在 x 轴附近， y 轴下降的很快。
代码语言：javascript复制ggplot(data=data_D1, aes(x=width,y=shear_stress))
geom_point(fill="red",size=3,pch=22)
geom_line(data=curve_D1)
我们可以对其做下对数变换，使得数据呈线性形式。使用 scale_x_log10() 和 scale_y_log10() 对刻度进行对数变换。内部参数这里不做解释，大家看着修改，就知道内部含义了。
scale_x_log10(expand=c(0,0),
breaks=c(0.5,1,2,5,10,20,50),
labels=c(0.5,1,2,5,10,20,50),
limits=c(0.5,50))
scale_y_log10(expand = c(0, 0),
labels = trans_format('log10', math_format(10^.x)),
breaks=c(0.001,0.01,0.1,1),
limits = c(0.001,1))
之后，我们对主题进行调整。尤其注意这里的 annotation_logticks(sides = "l")，这个代码可以增加对数坐标轴的标记（左边位置）。
annotation_logticks(sides = "l")
theme_void()
theme(
line = element_blank(),
# exclude everything outside axes bc it messes with positioning of grob in panel_D
text = element_blank(),
title = element_blank(),
axis.line.y = element_line(colour = "black")
)
ylab("shear stress (Pa)")
添加到另一张图形中
之后，将前面的图添加到另一张线性图中。首先把另一张图绘制出来：
代码语言：javascript复制ggplot(data=data_D2, aes(x=width,y=velocity))
geom_point(fill="blue",size=3,pch=21)
之后使用 annotation_custom(ggplotGrob(panel_D1)) 将前面那幅图添加到该图中。此时结果如下：
注意:
annotation_custom()是一个特殊的集合对象，用于静态注释。注释不会影响缩放。
这时，恭喜你两幅图已经合并啦！但是存在几个问题：
- 两幅图的 Y 轴重复了。这时候使用
scale_y_continuous()将原图的 Y 轴位置往右放置（position = "right"）。 - 但是变换完之后，左边标签没有，而左边的 Y 轴就是第一幅图得到的结果，我们需要添加缺失的标签。处理方式为：使用
sec.axis = sec_axis(~., name = "shear stress (Pa)",breaks=rescale(c(-3,-2,-1,0), to = c(0,1)),labels = c(expression("10"^"-3","10"^"-2", "10"^"-1", "10"^"0"))))手动添加坐标标签和数值。 - 两幅图的 x 轴不一致，使用
scale_x_log10()修改结果。 - 使用
annotation_logticks(sides = "b")添加 x 轴的 ticks。
代码语言：javascript复制所以，如果读者想使用这里的代码，还是需要一定能力看懂这些代码的～ 小编已经教会你们啦，以后要学会根据自己的结果更改代码噢～
scale_y_continuous(expand = c(0,0),
breaks = seq(0,1,0.1),
limits = c(0,1),
position = "right",
sec.axis = sec_axis(~., name = "shear stress (Pa)",breaks=rescale(c(-3,-2,-1,0), to = c(0,1)),labels = c(expression("10"^"-3","10"^"-2", "10"^"-1", "10"^"0"))))
scale_x_log10(expand=c(0,0),
breaks=c(0.5,1,2,5,10,20,50),
labels=c(0.5,1,2,5,10,20,50),
limits=c(0.5,50))
annotation_logticks(sides = "b")
后面就是添加一些线段，文字以及修改主题啦～前面几个推送都已经介绍过了，这里不做解释了。
代码语言：javascript复制 my_theme()
geom_line(color="blue")
geom_vline(xintercept = 1.1,
linetype="dashed")
geom_hline(yintercept = 0.9,
linetype="dashed")
theme(
plot.margin = unit(c(0.1,0,0,0.5), "cm"), # to match other panels
axis.title.y = element_text(margin = margin(r=1)),
axis.text.y = element_text(margin = margin(r=6)),
axis.text.y.right = element_text(margin = margin(l=7)),
axis.title.y.right = element_text(angle = 90)
)
xlab(expression(lumen~width~(mu*m)))
ylab("relative flow velocity")
annotate(geom = "text",x =6 ,y =0.85 ,label = "tau == 0.5~Pa",parse=T)
annotate(geom = "text",x =1.4 ,y =0.4 ,label = "b == 1.1*mu*m",parse=T,angle=90)
annotate(geom = "text",x =5.6 ,y =0.6 ,label = "tau == frac(4*mu*Q,pi*a*b^2)",parse=T)
完整代码
代码语言：javascript复制library(ggplot2)
base_size = 12
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_D1 = read.csv("./data_D1.csv")
data_D2 = read.csv("./data_D2.csv")
# format data to ggplot's liking
data_D = data.frame("width"=c(data_D1$width,data_D2$width),
"unit"=c(rep("shear_stress",nrow(data_D1)),
rep("velocity",nrow(data_D2))),
"value"=c(data_D1$shear_stress,data_D2$velocity)
)
head(data_D)
curve_D1 = data.frame(width=data_D1$width,
shear_stress=33.28/(pi*18*data_D1$width^2))
curve_D1
panel_D1 <- ggplot(data=data_D1,
aes(x=width,
y=shear_stress))
geom_point(fill="red",
size=3,
pch=22)
geom_line(data=curve_D1)
scale_x_log10(expand=c(0,0), # prevent gap between origin and first tick
breaks=c(0.5,1,2,5,10,20,50),
labels=c(0.5,1,2,5,10,20,50),
limits=c(0.5,50))
scale_y_log10( expand = c(0, 0),
# using trans_format from the scales package, but one can also use expressions
labels = trans_format('log10', math_format(10^.x)),
breaks=c(0.001,0.01,0.1,1),
limits = c(0.001,1)
)
annotation_logticks(sides = "l")
theme_void()
theme(
line = element_blank(),
# exclude everything outside axes bc it messes with positioning of grob in panel_D
text = element_blank(),
title = element_blank(),
axis.line.y = element_line(colour = "black")
)
ylab("shear stress (Pa)")
panel_D1
panel_D <-
ggplot(data=data_D2,
aes(x=width,
y=velocity))
# add plot of first dataset as grob as a trick to introduce two y-axes with different scalings
geom_point(fill="blue",
size=3,
pch=21)
annotation_custom(ggplotGrob(panel_D1))
scale_y_continuous(expand = c(0,0),
breaks = seq(0,1,0.1),
limits = c(0,1),
# putting the y axis of the second plot to the right
position = "right",
# now the secondary axis becomes the left axis
# we need the axis text title for panel_D1
# They were excluded in panel_D1 bc they were messing with the positioning
sec.axis = sec_axis(~.,
name = "shear stress (Pa)",
# rescale breaks bc sec_axis inherits scale from primary y axis
breaks=rescale(c(-3,-2,-1,0),
to = c(0,1)),
labels = c(expression("10"^"-3",
"10"^"-2",
"10"^"-1",
"10"^"0")))
)
scale_x_log10(expand=c(0,0),
breaks=c(0.5,1,2,5,10,20,50),
labels=c(0.5,1,2,5,10,20,50),
limits=c(0.5,50))
annotation_logticks(sides = "b")
my_theme()
geom_line(color="blue")
geom_vline(xintercept = 1.1,
linetype="dashed")
geom_hline(yintercept = 0.9,
linetype="dashed")
theme(
plot.margin = unit(c(0.1,0,0,0.5), "cm"), # to match other panels
axis.title.y = element_text(margin = margin(r=1)),
axis.text.y = element_text(margin = margin(r=6)),
axis.text.y.right = element_text(margin = margin(l=7)),
axis.title.y.right = element_text(angle = 90)
)
xlab(expression(lumen~width~(mu*m)))
ylab("relative flow velocity")
annotate(geom = "text",x =6 ,y =0.85 ,label = "tau == 0.5~Pa",parse=T)
annotate(geom = "text",x =1.4 ,y =0.4 ,label = "b == 1.1*mu*m",parse=T,angle=90)
annotate(geom = "text",x =5.6 ,y =0.6 ,label = "tau == frac(4*mu*Q,pi*a*b^2)",parse=T)
panel_D
小编有话说
本文主要学到的知识点如下：
- 使用
annotation_custom(ggplotGrob())图中添加其他图形； - 使用
scale_x_log10()和scale_y_log10()对刻度进行对数变换； - 使用
annotation_logticks(sides = "b")添加 x 轴的 ticks； - 使用
scale_y_continuous(position = "right")改变 Y 轴位置。
看完这篇文章，相信老板以后让你绘制双 Y 轴图，应该不在话下啦～ 如果觉得内容有用的话，小编写的有心的话。给小编来杯咖啡吧！动车上赶这篇推文也不容易额
