【问题标题】:How do I restrict the standard deviation bars in my barplot to a maximum value?如何将条形图中的标准偏差条限制为最大值?
【发布时间】:2017-03-18 11:22:06
【问题描述】:

我正在使用ggplot2 创建带有标准偏差条的条形图。我的数据框很大,但这里是一个截断的版本,例如:

SampleName  Target.ID   Maj.Allele.Freq SD  AVG.MAF
W15-P2-1    rs1005533   99.74811083 24.98883743 93.70753223
W15-P2-2    rs1005533   100 24.98883743 93.70753223
W15-P2-3    rs1005533   100 24.98883743 93.70753223
W15-P2-4    rs1005533   100 24.98883743 93.70753223
W15-P2-1    rs1005533   99.94819995 24.98883743 93.70753223
W15-P2-2    rs1005533   100 24.98883743 93.70753223
W15-P2-3    rs1005533   100 24.98883743 93.70753223
W15-P2-4    rs1005533   100 24.98883743 93.70753223
W21-P2-1    rs1005533   100 24.98883743 93.70753223
W21-P2-2    rs1005533   100 24.98883743 93.70753223
W21-P2-3    rs1005533   99.90044798 24.98883743 93.70753223
W21-P2-4    rs1005533   99.72375691 24.98883743 93.70753223
W21-P2-1    rs1005533   100 24.98883743 93.70753223
W21-P2-2    rs1005533   100 24.98883743 93.70753223
W21-P2-3    rs1005533   100 24.98883743 93.70753223
W21-P2-4    rs1005533   0   24.98883743 93.70753223
W15-P2-1    rs10092491  52.40641711 1.340954343 51.8604281
W15-P2-2    rs10092491  53.69923603 1.340954343 51.8604281
W15-P2-3    rs10092491  52.56689284 1.340954343 51.8604281
W15-P2-4    rs10092491  50.11764706 1.340954343 51.8604281
W15-P2-1    rs10092491  50.30094583 1.340954343 51.8604281
W15-P2-2    rs10092491  50.96277279 1.340954343 51.8604281
W15-P2-3    rs10092491  50.94102886 1.340954343 51.8604281
W15-P2-4    rs10092491  51.2849162  1.340954343 51.8604281
W21-P2-1    rs10092491  53.56976202 1.340954343 51.8604281
W21-P2-2    rs10092491  50.27861123 1.340954343 51.8604281
W21-P2-3    rs10092491  52.8358209  1.340954343 51.8604281
W21-P2-4    rs10092491  51.42585551 1.340954343 51.8604281
W21-P2-1    rs10092491  52.77890467 1.340954343 51.8604281
W21-P2-2    rs10092491  52.89017341 1.340954343 51.8604281
W21-P2-3    rs10092491  53.70786517 1.340954343 51.8604281
W21-P2-4    rs10092491  50  1.340954343 51.8604281

由于最后一列 (AVG.MAF) 中的平均值可能会产生超过最大值 100 的标准差条形图,因此该图显示的条形图超出了 y 轴 100 的限制。

这是创建上述图的代码:

pe1 = ggplot(half1, aes(x=Target.ID, y=AVG.MAF))+
 geom_bar(stat = "identity", position = "dodge", colour = "black", 
 width = 0.5, fill = "yellowgreen")+xlab("")+
 ylab("Average Major Allele Frequency")+
 labs(title="Allele Balance AmpliSeq Identity Sample P2")+
 geom_errorbar(aes(ymin = AVG.MAF-SD, ymax = AVG.MAF+SD), 
 width = 0.4, position = position_dodge(0.9), 
   size = 0.6)+
 theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = .5))

我尝试使用coord_cartesian 截断情节,但这会使情节看起来像是在隐藏一些数据:

以下是创建标准差条截断图的代码:

pe1 = ggplot(half1, aes(x=Target.ID, y=AVG.MAF))+geom_bar(stat = "identity", position = "dodge", colour = "black", width = 0.5, fill = "yellowgreen")+xlab("")+ylab("Average Major Allele Frequency")+labs(title="Allele Balance AmpliSeq Identity Sample P2")+geom_errorbar(aes(ymin = AVG.MAF-SD, ymax = AVG.MAF+SD), width = 0.4, position = position_dodge(0.9), size = 0.6)+theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = .5))+coord_cartesian(ylim=c(0,100))

似乎必须有一种方法可以将标准偏差条限制为我预期的 ymax 100,并且仍然保持顶部水平条在图中可见。有人知道怎么做吗?

【问题讨论】:

  • 为什么要通过截断 std 开发栏的顶部来歪曲标准差?
  • ...geom_errorbar(aes(ymin = AVG.MAF-SD, ymax = pmin(AVG.MAF+SD,100)... 会做你想做的事吗?几乎可以肯定的是,您现在低估了不确定性,可能是因为使用的基础错误模型不合适。
  • @NathanDay 和 Miff 你都给了我一些思考。感谢你们的 cmets 和可能的解决方案。
  • 要连接到 Nathan Day 的评论,也许标准差并不是你真正应该追求的,如果你可以获得自举置信区间呢?

标签: r ggplot2 standard-deviation


【解决方案1】:

除了人们在 cmets 中提出的问题外,还有其他一些注意事项:

  1. 您无需为数据的每一行添加重复平均值的列。相反,您可以使用Maj.Allele.Freq 中的实际数据值计算并绘制 ggplot 中的平均值。 (事实上​​,通过为每个Target.ID 重复使用平均值的 y 值列,您实际上是在绘制平均值条的多个副本,一个在另一个之上。)

    您还可以在 ggplot 之外汇总数据(即计算均值和标准差),然后使用汇总的数据框进行绘图。这在更复杂的情况下有时是必要的,但您可以在此处的 ggplot 中完成所有操作。

  2. 在我看来,积分比这里的条更有效。

下面的代码提供了点和条形版本,还显示了如何添加数据的标准差或数据平均值的 95% 置信区间。蓝线代表标准差,红线代表 95% 置信区间。

我提供了自举置信区间。要提供经典的正态置信区间,请从 mean_cl_boot 切换到 mean_cl_normal

如果您希望 y 轴下降到零,请添加 coord_cartesian(ylim=c(0,150)) 或您希望的任何最大 y 值(正如 cmets 所讨论的,为避免误导性图表,它应该在误差条顶部上方,无论条形代表 SD 还是 CI)。

ggplot(half1, aes(x=Target.ID, y=Maj.Allele.Freq)) +
  stat_summary(fun.data=mean_sdl, geom="errorbar", width=0.1, colour="blue") +
  stat_summary(fun.data=mean_sdl, geom="point", colour="blue", size=3) +
  stat_summary(fun.data = mean_cl_boot, colour="red", geom="errorbar", width=0.1) +
  stat_summary(fun.data = mean_cl_boot, colour="red", geom="point") +
  labs(x="", y="Average Major Allele Frequency", 
       title="Allele Balance AmpliSeq\nIdentity Sample P2") +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = .5)) 

ggplot(half1, aes(x=Target.ID, y=Maj.Allele.Freq)) +
  stat_summary(fun.y=mean, geom="bar", fill="yellowgreen", colour="black") +
  stat_summary(fun.data=mean_sdl, geom="errorbar", width=0.1, size=1, colour="blue") +
  stat_summary(fun.data = mean_cl_boot, colour="red", geom="errorbar", width=0.1, size=0.7) +
  labs(x="", y="Average Major Allele Frequency", 
       title="Allele Balance AmpliSeq\nIdentity Sample P2") +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = .5))   

您也可以将 SD 和 95% CI 放在同一个图上:

pnp = position_nudge(x=0.1)
pnm = position_nudge(x=-0.1)

ggplot(half1, aes(x=Target.ID, y=Maj.Allele.Freq)) +
  stat_summary(fun.data=mean_sdl, geom="errorbar", width=0.1, position=pnp, aes(colour="SD")) +
  stat_summary(fun.data=mean_sdl, geom="point", position=pnp, aes(colour="SD")) +
  stat_summary(fun.data = mean_cl_boot, geom="errorbar", width=0.1, 
               position=pnm, aes(colour="95% CI")) +
  stat_summary(fun.data = mean_cl_boot, geom="point", position=pnm, aes(colour="95% CI")) +
  labs(x="", y="Average Major Allele Frequency", colour="",
       title="Allele Balance AmpliSeq\nIdentity Sample P2") +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = .5))

【讨论】:

  • 我总是在这里学到新东西。我不知道我可以在 ggplot 中进行计算。谢谢你的详细回答。我可以使用它并根据需要进行一些调整。评论者是对的,切断 SD 条没有意义,但现在我知道该怎么做了。
  • 您可以创建一个汇总数据框,然后绘制它(在更复杂的情况下有时需要,但这里不是),或者您可以在 ggplot 中进行计算。无论哪种方式,您都不需要添加一个一遍又一遍地重复平均值的列,实际上这样做会导致每个条的多个副本被绘制在另一个之上。我已经更新了我的答案来讨论这个问题。
猜你喜欢
  • 1970-01-01
  • 2020-08-13
  • 2021-02-24
  • 2015-09-11
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2019-04-16
  • 2018-03-08
相关资源
最近更新 更多