gap.boxplot()添加刻度及解决线段、轴尺寸异常问题

原数据与代码 #

数据集 #

data_N3=structure(c(2.7240367196496, 4.36841055897346, 2.71177214885543, 
4.36849175181861, 17.7213767231357, 17.626304994383, 17.6468190136497, 
17.6321641647487, 1121.5558117271, 1116.19429181322, 1116.20063832111, 
1116.23541006948, NA, NA, NA, NA, 15.5774104850416, 5.54498649819063, 
15.4691309692694, 5.54438655223586, 11.0853210461256, 7.34658260987168, 
25.2590667564114, 12.1012005413784, 17.5713360077048, 8.44869502735733, 
17.4844139957775, 8.44912288685442), .Dim = c(4L, 7L))

原绘图代码 #

library("plotrix")
data_N3[,4] <- 0 #gap.boxplot不能包含NA值
par(mgp=c(4,2.3, 0), mar=c(10, 8, 2, 2), family="Times")

gap.boxplot(data_N3, axes=FALSE,
            col='white', family="Times", 
            whisklty = 0, staplelty = 0, medcol="white", xlab="", ylab="",
            border = 'white',
            gap = list(top=c(30,1110), bottom=c(NA,NA)) 
)

label=c("Control", 1, 2, 3, 1, 2, 3)
axis(1, c(1,2,3,4,5,6,7), labels=label, cex=2.8)
axis(2, c(10, 20, 30, 1110, 1120), cex=2.8)

par(family="Times")

#x轴线条和叉号
abline(v=4.5, lwd=4)
abline(v=1.5, lwd=4)
points(x=4, y=3.3, pch=4, cex=8, lwd=15, col='#238b45')

#x轴文本
mtext(c("Plants with aphids", "Floral resources"), side=1, at=c(3.1, 6), cex=3, line=5)
mtext("Resource removed", side=1, cex=3, line=7.8, at=3.5)
mtext(expression(paste("Population density " (number/m^{2}))), side=2, cex=3, line=4.2)


#可视化年度平均值
n=ncol(data_N3)
x0s=1:n-0.38  #适配箱线图宽度的线段起始x坐标
x1s=1:n+0.38

#添加数据为线段，替代箱线图
data_N3[,4] <- NA
y0s1=data_N3[1,]
y0s2=data_N3[2,]
y0s3=data_N3[3,]
y0s4=data_N3[4,]
segments(x0=x0s, x1=x1s, y0=y0s1, lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=y0s2, lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=y0s3, lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=y0s4, lwd=8, col='#238b45')

#无生境移除时的平均密度红线
abline(h=3.543178, lwd=6, col="red")

问题解决方案 #

问题1：y轴间隙上方的刻度及对应数值未显示 #

gap.boxplot会对间隙上方的y轴坐标做压缩转换，直接使用axis(2)无法正确映射刻度位置。需要借助plotrix包的gap.axis函数，并与boxplot的gap参数保持一致：

# 接收gap.boxplot返回值，包含坐标转换信息
bp <- gap.boxplot(data_N3, axes=FALSE,
            col='white', family="Times", 
            whisklty = 0, staplelty = 0, medcol="white", xlab="", ylab="",
            border = 'white',
            gap = list(top=c(30,1110), bottom=c(NA,NA)) 
)

# 绘制下方y轴刻度
axis(2, at=c(10,20,30), cex.axis=2.8)
# 绘制上方间隙处的y轴刻度
gap.axis(2, at=c(1110,1120), gap=list(top=c(30,1110)), cex.axis=2.8)

问题2：轴标签的尺寸参数（cex=2.8）未生效 #

axis函数中控制刻度标签大小的参数是cex.axis，而非cex。修改后代码：

label=c("Control", 1, 2, 3, 1, 2, 3)
axis(1, c(1,2,3,4,5,6,7), labels=label, cex.axis=2.8)

问题3：设置间隙后，data_N3[,3]的线段数据未在间隙上方显示 #

gap.boxplot对y轴坐标做了转换，直接使用原始y值绘制会超出可见区域。需要用返回的ytrans函数转换原始y值：

# 添加数据为线段，替代箱线图
data_N3[,4] <- NA
y0s1=data_N3[1,]
y0s2=data_N3[2,]
y0s3=data_N3[3,]
y0s4=data_N3[4,]

# 使用转换后的y坐标绘制线段
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s1), lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s2), lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s3), lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s4), lwd=8, col='#238b45')

完整修正后的代码 #

library("plotrix")
data_N3=structure(c(2.7240367196496, 4.36841055897346, 2.71177214885543, 
4.36849175181861, 17.7213767231357, 17.626304994383, 17.6468190136497, 
17.6321641647487, 1121.5558117271, 1116.19429181322, 1116.20063832111, 
1116.23541006948, NA, NA, NA, NA, 15.5774104850416, 5.54498649819063, 
15.4691309692694, 5.54438655223586, 11.0853210461256, 7.34658260987168, 
25.2590667564114, 12.1012005413784, 17.5713360077048, 8.44869502735733, 
17.4844139957775, 8.44912288685442), .Dim = c(4L, 7L))

data_N3[,4] <- 0 #gap.boxplot不能包含NA值
par(mgp=c(4,2.3, 0), mar=c(10, 8, 2, 2), family="Times")

# 接收gap.boxplot返回值，获取坐标转换函数
bp <- gap.boxplot(data_N3, axes=FALSE,
            col='white', family="Times", 
            whisklty = 0, staplelty = 0, medcol="white", xlab="", ylab="",
            border = 'white',
            gap = list(top=c(30,1110), bottom=c(NA,NA)) 
)

label=c("Control", 1, 2, 3, 1, 2, 3)
axis(1, c(1,2,3,4,5,6,7), labels=label, cex.axis=2.8)
# 分两部分绘制y轴刻度
axis(2, at=c(10, 20, 30), cex.axis=2.8)
gap.axis(2, at=c(1110, 1120), gap=list(top=c(30,1110)), cex.axis=2.8)

par(family="Times")

#x轴线条和叉号
abline(v=4.5, lwd=4)
abline(v=1.5, lwd=4)
# 叉号的y值也需要转换
points(x=4, y=bp$ytrans(3.3), pch=4, cex=8, lwd=15, col='#238b45')

#x轴文本
mtext(c("Plants with aphids", "Floral resources"), side=1, at=c(3.1, 6), cex=3, line=5)
mtext("Resource removed", side=1, cex=3, line=7.8, at=3.5)
mtext(expression(paste("Population density " (number/m^{2}))), side=2, cex=3, line=4.2)


#可视化年度平均值
n=ncol(data_N3)
x0s=1:n-0.38  #适配箱线图宽度的线段起始x坐标
x1s=1:n+0.38

#添加数据为线段，替代箱线图
data_N3[,4] <- NA
y0s1=data_N3[1,]
y0s2=data_N3[2,]
y0s3=data_N3[3,]
y0s4=data_N3[4,]

# 使用转换后的y坐标绘制线段
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s1), lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s2), lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s3), lwd=8, col='#238b45')
segments(x0=x0s, x1=x1s, y0=bp$ytrans(y0s4), lwd=8, col='#238b45')

#无生境移除时的平均密度红线，y值也需要转换
abline(h=bp$ytrans(3.543178), lwd=6, col="red")

内容的提问来源于stack exchange，提问作者Laura