具有多个分组级别的 R Barplot答案

【问题标题】：R Barplot with multiple grouping levels具有多个分组级别的 R Barplot
【发布时间】：2018-09-17 13:20:47
【问题描述】：

我一直在尝试创建一个有点复杂的分组条形图。

这是我的数据文件的 sn-p

"location"|"region"|"treatment"|"GB"
"Georgia"|"Keys"|"pre"|354
"Georgia"|"Keys"|"pre"|183
"Georgia"|"Keys"|"pre"|182
"Georgia"|"Keys"|"pre"|133
"Georgia"|"North East"|"pre"|44
"Georgia"|"North East"|"pre"|19
"Georgia"|"North East"|"pre"|70
"Georgia"|"North East"|"pre"|66
"Georgia"|"North West"|"pre"|102
"Georgia"|"North West"|"pre"|33
"Georgia"|"North West"|"pre"|106
"Georgia"|"North West"|"pre"|61
"Georgia"|"North West"|"pre"|101
"Georgia"|"Texas"|"pre"|150
"Georgia"|"Texas"|"pre"|187
"Georgia"|"Texas"|"pre"|152
"Georgia"|"Texas"|"pre"|148
"Georgia"|"Texas"|"pre"|100
"Maryland"|"Keys"|"pre"|637
"Maryland"|"Keys"|"pre"|52
"Maryland"|"Keys"|"pre"|43
"Maryland"|"Keys"|"pre"|156
"Maryland"|"Keys"|"pre"|38
"Maryland"|"North East"|"pre"|166
"Maryland"|"North East"|"pre"|91
"Maryland"|"North East"|"pre"|167
"Maryland"|"North East"|"pre"|104
"Maryland"|"North East"|"pre"|113
"Maryland"|"North West"|"pre"|370
"Maryland"|"North West"|"pre"|895
"Maryland"|"North West"|"pre"|198
"Maryland"|"North West"|"pre"|137
"Maryland"|"North West"|"pre"|168
"Maryland"|"Texas"|"pre"|95
"Maryland"|"Texas"|"pre"|331
"Maryland"|"Texas"|"pre"|163
"Maryland"|"Texas"|"pre"|90
"North Carolina"|"Keys"|"pre"|217
"North Carolina"|"Keys"|"pre"|91
"North Carolina"|"Keys"|"pre"|148
"North Carolina"|"Keys"|"pre"|208
"North Carolina"|"Keys"|"pre"|18
"North Carolina"|"North East"|"pre"|49
"North Carolina"|"North East"|"pre"|60
"North Carolina"|"North East"|"pre"|167
"North Carolina"|"North East"|"pre"|82
"North Carolina"|"North East"|"pre"|31
"North Carolina"|"North West"|"pre"|47
"North Carolina"|"North West"|"pre"|10
"North Carolina"|"North West"|"pre"|207
"North Carolina"|"North West"|"pre"|70
"North Carolina"|"North West"|"pre"|214
"North Carolina"|"Texas"|"pre"|183
"North Carolina"|"Texas"|"pre"|162
"North Carolina"|"Texas"|"pre"|94
"North Carolina"|"Texas"|"pre"|102
"South Carolina"|"Keys"|"pre"|101
"South Carolina"|"Keys"|"pre"|155
"South Carolina"|"Keys"|"pre"|85
"South Carolina"|"Keys"|"pre"|67
"South Carolina"|"Keys"|"pre"|60
"South Carolina"|"North East"|"pre"|173
"South Carolina"|"North East"|"pre"|148
"South Carolina"|"North East"|"pre"|575
"South Carolina"|"North East"|"pre"|96
"South Carolina"|"North West"|"pre"|51
"South Carolina"|"North West"|"pre"|86
"South Carolina"|"North West"|"pre"|34
"South Carolina"|"North West"|"pre"|67
"South Carolina"|"Texas"|"pre"|124
"South Carolina"|"Texas"|"pre"|155
"South Carolina"|"Texas"|"pre"|183
"South Carolina"|"Texas"|"pre"|101
"Georgia"|"Keys"|"post"|344
"Georgia"|"Keys"|"post"|241
"Georgia"|"Keys"|"post"|486
"Georgia"|"Keys"|"post"|191
"Georgia"|"North East"|"post"|128
"Georgia"|"North East"|"post"|14
"Georgia"|"North East"|"post"|192
"Georgia"|"North East"|"post"|298
"Georgia"|"North West"|"post"|540
"Georgia"|"North West"|"post"|236
"Georgia"|"North West"|"post"|172
"Georgia"|"North West"|"post"|87
"Georgia"|"Texas"|"post"|357
"Georgia"|"Texas"|"post"|221
"Georgia"|"Texas"|"post"|131
"Georgia"|"Texas"|"post"|55
"Maryland"|"Keys"|"post"|196
"Maryland"|"Keys"|"post"|85
"Maryland"|"Keys"|"post"|90
"Maryland"|"Keys"|"post"|530
"Maryland"|"North East"|"post"|477
"Maryland"|"North East"|"post"|447.253
"Maryland"|"North East"|"post"|509
"Maryland"|"North East"|"post"|64
"Maryland"|"North West"|"post"|1204
"Maryland"|"North West"|"post"|756
"Maryland"|"North West"|"post"|635
"Maryland"|"North West"|"post"|948
"Maryland"|"Texas"|"post"|740
"Maryland"|"Texas"|"post"|567
"Maryland"|"Texas"|"post"|549
"Maryland"|"Texas"|"post"|271
"North Carolina"|"Keys"|"post"|173
"North Carolina"|"Keys"|"post"|114
"North Carolina"|"Keys"|"post"|1159
"North Carolina"|"Keys"|"post"|113
"North Carolina"|"North East"|"post"|176
"North Carolina"|"North East"|"post"|187
"North Carolina"|"North East"|"post"|279
"North Carolina"|"North East"|"post"|182
"North Carolina"|"North West"|"post"|103
"North Carolina"|"North West"|"post"|230
"North Carolina"|"North West"|"post"|117
"North Carolina"|"North West"|"post"|143
"North Carolina"|"Texas"|"post"|358
"North Carolina"|"Texas"|"post"|458
"North Carolina"|"Texas"|"post"|102

我的水平是“治疗”（前后）、“地区”[收集样本的地方]（钥匙、西北、东北和德克萨斯州）和“地点”[进行实验的地方]（马里兰州，乔治亚州、北卡罗来纳州和南卡罗来纳州）。我正在测量每个样本的“GB”。

我想得到一个图，显示每个实验位置内每个样本区域的前后 GB

我可以绘制一个图，显示区域内治疗前和治疗后的平均 GB，没有问题。

当我按如下方式向脚本添加位置时

data <- tapply(dat$GB, list(treat,regi,loci), mean)

我可以让 R 计算每个位置的每个区域的前后的平均 GB。快乐的日子！

但是当我尝试将其绘制成条形图时，我在 R 中得到一个错误，说“高度”必须是向量或矩阵'

这是我写的脚本。您会看到我对每个级别的数据进行了排序，然后使用这个新的“排序”数据为每个位置的每个区域内的每个前后样本创建均值和 SE。

然后我尝试将这个“平均”数据用于绘图。

谢谢

马特

#treatments pre and post between regions across locations
loci = factor (dat$location, levels=c("Georgia","Maryland","North Carolina","South Carolina"))
regi = factor(dat$region, levels=c("Keys","North West","North East", "Texas"))
treat = factor(dat$treatment, levels=c("Pre","Post"),ordered=TRUE)

data <- tapply(dat$GB, list(treat,regi,loci), mean)
ses<- tapply(dat$GB, list(treat,regi,loci),function (x) sd(x)/sqrt(length(x)))
lower<-data-ses
upper<-data+ses

my.plot<-barplot(data,beside=TRUE, legend= F, main="",ylim=c(0,400),xlab="",ylab="",cex=1.3, 
                 cex.lab=1.3, cex.axis=1.3)
arrows(my.plot, data, my.plot, lower, angle=90,length=.1)
arrows(my.plot, data, my.plot, upper, angle=90, length=.1)
mtext(expression(ug~GB~(ml^-1)), side=2, line=3, cex=1.5)

数据。

dat1 <-
structure(list(location = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L), .Label = c("Georgia", "Maryland", "North Carolina", 
"South Carolina"), class = "factor"), region = structure(c(1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 
4L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 
3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 4L, 4L, 4L), .Label = c("Keys", "North East", 
"North West", "Texas"), class = "factor"), treatment = structure(c(2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L), .Label = c("post", "pre"), class = "factor"), 
    GB = c(354, 183, 182, 133, 44, 19, 70, 66, 102, 33, 106, 
    61, 101, 150, 187, 152, 148, 100, 637, 52, 43, 156, 38, 166, 
    91, 167, 104, 113, 370, 895, 198, 137, 168, 95, 331, 163, 
    90, 217, 91, 148, 208, 18, 49, 60, 167, 82, 31, 47, 10, 207, 
    70, 214, 183, 162, 94, 102, 101, 155, 85, 67, 60, 173, 148, 
    575, 96, 51, 86, 34, 67, 124, 155, 183, 101, 344, 241, 486, 
    191, 128, 14, 192, 298, 540, 236, 172, 87, 357, 221, 131, 
    55, 196, 85, 90, 530, 477, 447.253, 509, 64, 1204, 756, 635, 
    948, 740, 567, 549, 271, 173, 114, 1159, 113, 176, 187, 279, 
    182, 103, 230, 117, 143, 358, 458, 102)), .Names = c("location", 
"region", "treatment", "GB"), class = "data.frame", row.names = c(NA, 
-120L))

【问题讨论】：

嗨@matt，请提供一些数据（例如dput(dat) 的输出），以便其他人可以重现您的问题。您还应该展示您如何“尝试添加位置”，这将有助于人们调试您的代码。
嗨@dash2，谢谢你的建议。我在帖子中添加了一些数据，并解释了我是如何按照您的建议尝试将位置添加到脚本中的
现在数据在这里，很明显这是缺少分隔符。从dput 发布输出或发布原始文本文件。当应该有分隔符时，我拒绝修复数据输入问题。
可能没有那么多数据......你能做一个minimal reproducible example吗？

标签： r statistics bar-chart

【解决方案1】：

为了找到解决您问题的方法，我们先来看看您的代码。

当我按如下方式向脚本添加位置时
data <- tapply(dat$GB, list(treat,regi,loci), mean)
我可以让 R 计算每个位置的每个区域的前后的平均 GB。

你是怎么做到的？我的 R 给了我一个空的三维数组，其中包含 2 行、4 列和 4 个三维元素。如果我跑

all(is.na(data))

我得到了结果TRUE。所以我在这里没有快乐的日子。但是，tapply (?tapply) 的帮助页面帮助我更进一步：表示索引的列表应该包含一个或多个因子长度为 X，但 lapply 会将合适的向量转换为本身的因素。所以如果我跑

data <- tapply(dat$GB, list(dat$treatment,dat$region,dat$location), mean)

上面提到的数组填充了平均值。置信区间边界也是如此。

但现在开始真正的麻烦。 barplot() 仍然抛出错误“'height' must be a vector or a matrix”。意味着它不会接受一个三维数组作为条形高度的输入，即作为数据。我想这是因为没有显示 3-D 数据的“标准”，即条形图中每个元素有 4 个信息（它是值，它是处理，它是采样区域，它是采样子区域）。即使使用适当的 3-D 条形图，如 this post 所示，您仍然需要在一个轴上压缩两个信息（这里仅通过结合“位置”和“区域”来初步完成，请参阅后面的图表以获取详细解决方案以及如何准确区域和位置已合并为“ExactReg”）：

library(latticeExtra)
cloud(GB ~ exactReg + treatment, data = dat, panel.3d.cloud=panel.3dbars, xbase=0.5, ybase=0.2)

产生：（以后可以添加标准偏差甚至置信区间，如here 所示（我很确定这也是 R 上的 SO 帖子的一部分，但我找不到它）。）

所以我们必须想出一个替代解决方案。在这里，这取决于您要强调数据的哪个方面。恕我直言，Hmisc::bpplot 提供了有关数据分布的最佳信息，因此很可能有关形成此数据的过程的信息。但它侧重于稳健的统计数据，因此不显示平均值和标准差（但显示中位数和四分位数）。但是，这些可以稍后添加，但同时添加它们最终会导致相当拥挤的情节，在我看来这已经超载了。因此，我只添加了手段，展示了这个概念。 bpplots 也确实比箱线图（以下称为箱线图）更不清晰地显示异常值。另一个缺点是，到目前为止，它们还没有“at”和“add”参数，因此很难将它们很好地分组在一行中，“pre”和“post”治疗测量彼此相邻，这是最好的用于比较“前”和“后”。我现在会这样描述它们：

library(TeachingDemos)
set.seed(7)
ms.sun <- rbind(ms.polygon(n = 50, r = 1),
                c(NA, NA),
                ms.polygon(n = 50, r = 0.03),
                c(0, 0),
                c(0, 0.0075),
                c(0.0075, 0.0075),
                c(0.0075, 0),
                c(0.0075, -0.0075),
                c(0, -0.0075),
                c(-0.0075, -0.0075),
                c(-0.0075, 0),
                c(-0.0075, 0.0075),
                matrix(data = runif(80, min = -0.029, max = 0.029), ncol = 2),
                c(NA, NA))


library(Hmisc)
library(pBrackets)
x11(10,8)
par(mfrow = c(2,1), mar = c(0,4,2,1), oma = c(8,0,0,0))
      #2 plots in the same columns (one above the other) -- no margin at bottom -- more space for x-axis below lower plot
xpos <- bpplot(split(dat[which(dat$treatment == "pre"), "GB"], dat[which(dat$treatment == "pre"), "exactReg"]), xlab = "", name = FALSE, ylab = "pre treatment", main = "GB", plotopts = list("ylim" = c(0, 1300)))
      #we need the positions of the barplots on the x-axis later, that's what is saved in "xpos"
      #I want to give pre treatment and post treatment the same value range on the y-axis, for comparability. Therefore the "ylim" argument in plotopts (in both calls to bpplot)
  abline(v = c(seq(from = 0, by = 1.2, length.out = length(levels(dat$exactReg)))), lty = "dashed", col = gray(0.7))
      #add vertical lines in both plots to make assignment easy
  my.symbols(xpos, by(dat[which(dat$treatment == "pre"), "GB"], dat[which(dat$treatment == "pre"), "exactReg"], mean),
             symb = ms.sun, xsize = 0.3)
      #add mean values
  legend("topright", pch = c(NA), legend = c("mean value"), bg = "white")
  my.symbols(16.85487, 1238.545, symb = ms.sun, xsize = 0.3)
      #add a legend to clarify what this extra circle means

par(mar = c(1,4,0,1))
      #no margin on top of plot, same margins left and right!
bpplot(split(dat[which(dat$treatment == "post"), "GB"], dat[which(dat$treatment == "post"), "exactReg"]), main = "", ylab = "post treatment", plotopts = list("ylim" = c(0, 1300)))
  abline(v = c(seq(from = 0, by = 1.2, length.out = length(levels(dat$exactReg)))), lty = "dashed", col = gray(0.7))
  my.symbols(xpos, by(dat[which(dat$treatment == "post"), "GB"], dat[which(dat$treatment == "post"), "exactReg"], mean),
             symb = ms.sun, xsize = 0.3)
  #The Error occurs, because not every group has got values to calculate a mean from
  #I think, in this case there is no extra legend needed
  mtext("|", side = 1, line = -0.15, at = xpos)   #axis ticks
  text(x = xpos, y = par("usr")[3] - 1/2.5* diff(c(par("usr")[3], par("usr")[4])),
       labels = rep(levels(dat$region), length(levels(dat$location))), srt = 90, adj = c(0, NA), xpd = NA)
  brackets(xpos[4] + 0.5*diff(xpos[4:5]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           xpos[1] - 0.5*diff(xpos[1:2]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  brackets(xpos[8] + 0.5*diff(xpos[8:9]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           xpos[5] - 0.5*diff(xpos[4:5]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  brackets(xpos[12] + 0.5*diff(xpos[12:13]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           xpos[9] - 0.5*diff(xpos[8:9]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  brackets(xpos[16] + 0.5*diff(xpos[15:16]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           xpos[13] - 0.5*diff(xpos[12:13]), par("usr")[3] - 0.45*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  #mtext(dat$location[1], side = 1, line = 5, at = mean(xpos[2:3]), outer = T)   #I can't tell for now, why this doesn't work
  text(x = seq(from = mean(xpos[2:3]), by = length(levels(dat$region))*diff(xpos[1:2]), length.out = length(levels(dat$location))),
       y = par("usr")[3] - (1/1.8)* diff(c(par("usr")[3], par("usr")[4])), labels = levels(dat$location), xpd = NA)
  rm(xpos)

给予：

注意：我的解决方案只能像这样工作，如果您在“帖子”处理组的每个位置的每个区域也至少获得一个值。现在“帖子”中的框治疗”情节与所述位置不一致。可能还有很多工作要做（比如让太阳符号中的点更大，在大括号之间添加一点空间......），但概念应该很清楚......

另一方面，箱线图提供的数据分布信息较少，但可以更清楚地显示异常值，并且可以理想地分组以比较治疗前后的数据。它们也是稳健统计的工具，不显示均值和标准差，而是显示中位数和四分位数。我再次添加了仅表示：

x11(10,5)
par(mfrow = c(1,1), mar = c(8,4,3,1))
boxplot(dat[which(dat$treatment == "pre"), "GB"] ~ dat[which(dat$treatment == "pre"), "exactReg"],
        at = seq(from = 1, by = 3, length.out = length(levels(dat$exactReg))), col = "brown1",
        xlim = c(1, 3*length(levels(dat$exactReg))-1), ylim = c(0, max(dat$GB)),
        names = character(length(levels(dat$exactReg))), xaxt = "n", ylab = "GB", main = "Header")
  my.symbols(seq(from = 1, by = 3, length.out = length(levels(dat$exactReg))),
             by(dat[which(dat$treatment == "pre"), "GB"], dat[which(dat$treatment == "pre"), "exactReg"], mean),
             symb = ms.sun, xsize = 0.7)
  boxplot(dat[which(dat$treatment == "post"), "GB"] ~ dat[which(dat$treatment == "post"), "exactReg"],
        at = seq(from = 2, by = 3, length.out = length(levels(dat$exactReg))), col = "lightblue",
        names = character(length(levels(dat$exactReg))), xaxt = "n", ylab = "", add = TRUE)
  my.symbols(seq(from = 2, by = 3, length.out = length(levels(dat$exactReg))),
             by(dat[which(dat$treatment == "post"), "GB"], dat[which(dat$treatment == "post"), "exactReg"], mean),
             symb = ms.sun, xsize = 0.7)
        #The Error occurs, because not every group has got values to calculate a mean from
  legend("topright", fill = c("brown1", "lightblue", NA), border = c("black", "black", NA), legend = c("prior to treatment", "after treatment", "mean value"))
  my.symbols(40.5, 964.8124, symb = ms.sun, xsize = 0.7)
  mtext("|", side = 1, line = -0.15, at = seq(from = 1.5, by = 3, length.out = length(levels(dat$exactReg))))   #axis ticks
  text(x = seq(from = 1.5, by = 3, length.out = length(levels(dat$exactReg))),
       y = par("usr")[3] - 1/2.5* diff(c(par("usr")[3], par("usr")[4])),
       labels = rep(levels(dat$region), length(levels(dat$location))), srt = 90, adj = c(0, NA), xpd = TRUE)
  brackets(11.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           0.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  brackets(23.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           12.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  brackets(35.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           24.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  brackets(47.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           36.5, par("usr")[3] - 0.41*diff(c(par("usr")[3], par("usr")[4])),
           h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
  mtext(levels(dat$location), side = 1, line = 6.5, at = seq(6, by = 12, length.out = length(levels(dat$location))) , xpd = TRUE)

这个解决方案在缺失值方面与其他解决方案存在相同的问题：这些框不是在它们应该在的位置绘制，而是一个接一个地绘制，即一旦一个框完全缺失，下面的框就是一个插槽太左等等。 Unfurtanately 我不能用 bpplots 做到这一点。

最后，如果您只想显示您的处理确实有所作为（而不是如何），并且如果异常值在您的数据中不是太大的问题，那么条形图确实最适合这个任务，因为您可以最好地使用经典统计数据。
我用你的方法和tapply 来计算数据。它清楚地向我表明，确实缺少最后四组的数据，因此其他图中的顺序没有混乱。您可以实现类似于箱线图的条形图：

data <- tapply(dat$GB, list(dat$exactReg, dat$treatment), mean)
ses <- tapply(dat$GB, list(dat$exactReg, dat$treatment), function(x){sd(x)/sqrt(length(x))})
lower <- data-ses
upper <- data+ses

x11(10,8)
par(mar = c(12,5,3,1), las = 3)
my.plot <- barplot(data[, "pre"], space = c(0, rep(2, 15)), main = "Header", ylim = c(0, 1100), xlab = "", ylab = "", names.arg = "", col = "lightblue")
arrows(my.plot, data[, "pre"], my.plot, lower[, "pre"], angle = 90, length = 0.08)
arrows(my.plot, data[, "pre"], my.plot, upper[, "pre"], angle = 90, length = 0.08)
my.plot <- barplot(data[, "post"], space = c(1, rep(2, 15)), main = "", ylim = c(0, 1100), xlab = "", ylab = "", names.arg = "", col = "brown1", add = TRUE)
arrows(my.plot, data[, "post"], my.plot, lower[, "post"], angle = 90, length = 0.08)
arrows(my.plot, data[, "post"], my.plot, upper[, "post"], angle = 90, length = 0.08)
mtext(expression(ug~GB~(ml^-1)), side=2, line=3, cex=1.5)
legend("topright", fill = c("lightblue", "brown1"), legend = c("prior to treatment", "after treatment"))
axis(1, at = seq(from = 1.5, by = 3, length.out = length(levels(dat$exactReg))), labels = FALSE, lwd.ticks = 2)
text(x = seq(from = 1.5, by = 3, length.out = length(levels(dat$exactReg))),
     y = par("usr")[3] - 1/5.1* diff(c(par("usr")[3], par("usr")[4])),
     labels = rep(levels(dat$region), length(levels(dat$location))), srt = 90, adj = c(0, NA), xpd = TRUE)
brackets(11.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         0.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
brackets(23.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         12.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
brackets(35.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         24.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
brackets(47.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         36.5, par("usr")[3] - 0.2*diff(c(par("usr")[3], par("usr")[4])),
         h = 0.04*diff(c(par("usr")[3], par("usr")[4])), xpd = NA)
mtext(levels(dat$location), side = 1, line = 6.2, at = seq(6, by = 12, length.out = length(levels(dat$location))) , xpd = TRUE)

【讨论】：

感谢您的有用建议@mattu！这是一种展示数据的有趣方式，我想我会进一步研究。
非常酷！！这种呈现数据的方法看起来棒极了！我有每个地区的数据，所以我认为如果我使用完整的数据文件，它可能会起作用。我会试一试。非常感谢您的帮助！