ggplot2绘图基础：

今天我们学习第三章的最后几节，其中的“绘制地图”部分，因为我木有顺利安装maps package，而且在我们的工作中也不常用，暂时跳过。下面继续~

6. 添加误差线和误差范围

数据中的不确定信息的展示也很重要。

ggplot2中，四类几何对象可以用于这项工作，这取决于x的值是离散型还是连续型，以及我们是否想展示区间中的中心值：

• 离散型变量+区间： geom_errorbar(), geom_linerange()

• 离散型变量+区间+中间值： geom_crossbar(), geom_pointrange()

• 连续型变量+区间： geom_ribbon()

• 连续型变量+区间+中间值： geom_smooth(stat="identity")

以上函数默认，我们对给定x时y的值域和分布情况感兴趣，所以使用了图形属性 ymin和 ymax来确定y的值域。

如下例：

1. y <- c(18, 11, 16)

2. df <- data.frame(x = 1:3, y = y, se = c(1.2, 0.5, 1.0))

3. base <- ggplot(df, aes(x, y, ymin = y - se, ymax = y + se))

4. ## 箱线图的形式，离散型变量的区间和中间值

5. base + geom_crossbar()  

6. ## 点线图的形式，离散型变量的区间和中间值

7. base + geom_pointrange()

8. ## 连续型变量，展示区间和中间值

   base + geom_smooth(stat="identity")

9. ## 箱线图的形式，离散型变量只显示区间

10. base + geom_errorbar()

11. ## 点线图的形式，离散型变量只显示区间

12. base + geom_linerange()

13. ## 连续型变量，只显示区间

14. base + geom_ribbon()

以上六个命令依次生成如下六张图片：

实际上，添加标准误差线的方法有很多，除此之外还有其他可用的计算方法（PS. 这个方法我理解的前提是你已经知道误差范围，并具体到数值，才能通过这种代码运行）。在后面的第十一章还有更详细的从复杂数据模型中计算置信区间（ confidence intervals）的方法。

7. 加权数据（Weighted Data）

加权数与加权平均数的概念相同，不同于普通的平均数，在加权数的计算中需要考虑每个数据不同的比例权重。即将各数值乘以相应的权数，然后加总求和得到总体值，再除以总的单位数。加权平均值的大小不仅取决于总体中各单位的数值（变量值）的大小，而且取决于各数值出现的次数（频数），由于各数值出现的次数对其在平均数中的影响起着权衡轻重的作用，因此叫做权数。 （from. 百度百科）

在处理整合数据（aggregated data）的时候，数据的每一行可能代表了多种观测值，这是我们需要通过某种标准/方式把权重变量考虑到其中。

以2000年美国人口普查，东西部各州的统计数据集（midwest）为例。此数据主要包括的是比例型数据（eg. 白种人比例，贫困线下人口比例，大学以上学历人口比例）以及每个地区的信息（面积、人口总数、人口密度等）。

1. library(ggplot2)

2. midwest

3. # A tibble: 437 x 28

4.     PID county    state   area poptotal popdensity popwhite popblack popamerindian popasian

5.   <int> <chr>     <chr>  <dbl>    <int>      <dbl>    <int>    <int>         <int>    <int>

6. 1   561 ADAMS     IL    0.0520    66090       1271    63917     1702            98      249

7. 2   562 ALEXANDER IL    0.0140    10626        759     7054     3496            19       48

8. 3   563 BOND      IL    0.0220    14991        681    14477      429            35       16

9. 4   564 BOONE     IL    0.0170    30806       1812    29344      127            46      150

10. 5   565 BROWN     IL    0.0180     5836        324     5264      547            14        5

11. 6   566 BUREAU    IL    0.0500    35688        714    35157       50            65      195

12. 7   567 CALHOUN   IL    0.0170     5322        313     5298        1             8       15

13. 8   568 CARROLL   IL    0.0270    16805        622    16519      111            30       61

14. 9   569 CASS      IL    0.0240    13437        560    13384       16             8       23

15. 10   570 CHAMPAIGN IL    0.0580   173025       2983   146506    16559           331     8033

16. # ... with 427 more rows, and 18 more variables: popother <int>, percwhite <dbl>,

17. #   percblack <dbl>, percamerindan <dbl>, percasian <dbl>, percother <dbl>, popadults <int>,

18. #   perchsd <dbl>, percollege <dbl>, percprof <dbl>, poppovertyknown <int>,

19. #   percpovertyknown <dbl>, percbelowpoverty <dbl>, percchildbelowpovert <dbl>,

20. #   percadultpoverty <dbl>, percelderlypoverty <dbl>, inmetro <int>, category <chr>

权重变量能在很大程度上影响图形内容和观察结论。有两种参数可以调整权重属性。

• 首先，对于点、线这些简单的几何对象，可以根据 size改变点的大小：

例，白人百分比和贫困线下人口数量的散点图：

1. # Unweighted

2. ggplot(midwest, aes(percwhite, percbelowpoverty)) +

3.  geom_point()

按照人口总数的绝对数字设置权重 size=poptotal/1e6，可以看到图片中也自动生成了对应的图例：

1. # Weight by population

2. ggplot(midwest, aes(percwhite, percbelowpoverty)) +

3.  geom_point(aes(size = poptotal / 1e6)) +

4.  scale_size_area("Population\n(millions)", breaks = c(0.5, 1, 2, 4))

• 对于更加复杂、涉及到统计变换的情况，我们通过修改weight图形属性来表现权重。这些权重将被传递给统计汇总计算函数。在权重有意义的情况下，各种元素基本都支持权重的设定，例如各类smooth平滑器、箱线图、分位回归、直方图以及密度图等等。

它没有对应的图例，我们无法看到这个权重变量，但是它却会改变统计汇总的结果。

下图可以看到，作为权重的人口密度，会影响白种人比例和贫困线下人口比例这个散点图的分布情况：

1. # Unweighted

2. ggplot(midwest, aes(percwhite, percbelowpoverty)) +

3.  geom_point() +

4.  geom_smooth(method = lm, size = 1)

加权后（可以看到散点的位置没有变，但是大小变了，蓝色平滑线的角度也根据人口总量权重发生了改变）：

1. # Weighted by population

2. ggplot(midwest, aes(percwhite, percbelowpoverty)) +

3.  geom_point(aes(size = poptotal / 1e6)) +

4.  geom_smooth(aes(weight = poptotal), method = lm, size = 1) +

5.  scale_size_area(guide = "none")

我们使用总人口数作为权重修改直方图和密度图的时候，加入 weight 权重会大大改变分析结果：

首先简单地对郡的数量分布进行观察：

1. # Unweighted

2. ggplot(midwest, aes(percbelowpoverty)) +

3.  geom_histogram(binwidth = 1) +

4.  ylab("Counties")

加“人口数量”作为权重后，转向在人口数量分布的基础上对郡县数量的观察：

1. # Weighted by population

2. ggplot(midwest, aes(percbelowpoverty)) +

3.  geom_histogram(aes(weight = poptotal), binwidth = 1) +

4.  ylab("Population (1000s)")

8. 展示数据分布

首先我们引入一个大数据集——“钻石数据集”（diamonds）

1. library(ggplot2)

2. diamonds

3. # A tibble: 53,940 x 10

4.   carat cut       color clarity depth table price     x     y     z

5.   <dbl> <ord>     <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>

6. 1 0.230 Ideal     E     SI2      61.5  55.0   326  3.95  3.98  2.43

7. 2 0.210 Premium   E     SI1      59.8  61.0   326  3.89  3.84  2.31

8. 3 0.230 Good      E     VS1      56.9  65.0   327  4.05  4.07  2.31

9. 4 0.290 Premium   I     VS2      62.4  58.0   334  4.20  4.23  2.63

10. 5 0.310 Good      J     SI2      63.3  58.0   335  4.34  4.35  2.75

11. 6 0.240 Very Good J     VVS2     62.8  57.0   336  3.94  3.96  2.48

12. 7 0.240 Very Good I     VVS1     62.3  57.0   336  3.95  3.98  2.47

13. 8 0.260 Very Good H     SI1      61.9  55.0   337  4.07  4.11  2.53

14. 9 0.220 Fair      E     VS2      65.1  61.0   337  3.87  3.78  2.49

15. 10 0.230 Very Good H     VS1      59.4  61.0   338  4.00  4.05  2.39

16. # ... with 53,930 more rows

这个数据集中包含近54000种钻石的不同属性，如，4个C（克拉、切割、颜色、纯度）、以及各种长度属性和价格。我们借此数据集来展示数据分布。

数据分布可以分为条件分布和联合分布。

条件分布：对于二维随机变量(X，Y)，可以考虑在其中一个随机变量取得(可能的)固定值的条件下，另一随机变量的概率分布，这样得到的X或Y的概率分布叫做条件概率分布，简称条件分布。

联合分布：联合分布函数亦称多维分布函数，随机向量的分布函数，以二维情形为例，若（X，Y）是二维随机向量，x、y是任意两个实数，则称二元函数。

一些几何对象可以用于展示数据的分布，具体使用哪种，取决于分布的维度、分布是连续型还是离散型，以及我们使用条件分布还是联合分布

对于一维连续型分布，最重要的几何对象是直方图 geom_histogram()。

例如，展示diomands数据集中的depth变量。

1. ## 改变组距宽度 binwidth()

2. ## 限制x轴坐标的范围

3. ggplot(diamonds, aes(depth)) +

4.  geom_histogram(binwidth = 0.1) +

5.  xlim(55, 70)

6. ## Warning message: Removed 45 rows containing non-finite values (stat_bin).

如果你想展示更多组间信息，那你可以使用下面的几种参数：

• facet_wrap(~var)：分面

• geom_freqpoly()：颜色和频率多边形

• geom_histogram(position="fill")：条件密度图

按照cut变量的不同取值进行线条的上色：

1. ggplot(diamonds, aes(depth)) +

2.  geom_freqpoly(aes(colour = cut), binwidth = 0.1, na.rm = TRUE) +

3.  xlim(58, 68) +

4.  theme(legend.position = "none")

按照cut变量的不同取值进行填充柱状图：

1. ggplot(diamonds, aes(depth)) +

2.  geom_histogram(aes(fill = cut), binwidth = 0.1, position = "fill", na.rm = TRUE) +

3.  xlim(58, 68) +

4.  theme(legend.position = "none")

geom_density 一维密度曲线图：

如果不设置其他参数，实际上就是直方图的平滑曲线版本，无法回溯数据本身，意义不大。

1. ggplot(diamonds, aes(depth)) +

2.  geom_density(na.rm = TRUE) +

3.  xlim(58, 68) +

4.  theme(legend.position = "none")

如果用cut变量的不同取值进行上色，就会得到一个信息量更大的版本：

1. ggplot(diamonds, aes(depth, fill = cut, colour = cut)) +

2.  geom_density(alpha = 0.2, na.rm = TRUE) +

3.  xlim(58, 68) +

4.  theme(legend.position = "none")

如果想比较两个（以上）变量之间的关系，那就需要使用其他函数，诸如：

• geom_boxplot()：箱型图，也称箱须图(Box-Whisker-Plot）

横轴取值是离散型变量时：

1. ggplot(diamonds, aes(clarity, depth)) +

2.  geom_boxplot()

横轴取值是连续型变量时，用 cut_width设置离散区间，用 xlim设置横轴取值范围：

1. ggplot(diamonds, aes(carat, depth)) +

2.  geom_boxplot(aes(group = cut_width(carat, 0.1))) +

3.  xlim(NA, 2.05)

4. #> Warning: Removed 997 rows containing non-finite values

5. #> (stat_boxplot).

• geom_violin()：小提琴图

横轴取值是离散型变量时：

1. ggplot(diamonds, aes(clarity, depth)) +

2.  geom_violin()

横轴取值为连续型变量时：

1. ggplot(diamonds, aes(carat, depth)) +

2.  geom_violin(aes(group = cut_width(carat, 0.1))) +

3.  xlim(NA, 2.05)

---

参考资料：

1. Hadley Wickham(2016). ggplot2. Springer International Publishing. doi: 10.1007/978-3-319-24277-4

2. 《R语言应用系列丛书·ggplot2:数据分析与图形艺术》

快快和我一起上车，请关注：生信小白学习记

猜你喜欢

• 热文：1高分文章 2不可或缺的人 3图表规范

• 一文读懂：1微生物组 2寄生虫益处 3进化树

• 必备技能：1提问 2搜索  3Endnote

• 文献阅读 1热心肠 2SemanticScholar 3geenmedical

• 扩增子分析：1图表解读 2分析流程 3统计绘图  4功能预测

• 科研经验：1云笔记  2云协作 3公众号  

• 系列教程：1Biostar 2微生物组  3宏基因组

• 生物科普 1肠道细菌 2人体上的生命 3生命大跃进  4细胞的暗战 5人体奥秘

写在后面

为鼓励读者交流、快速解决科研困难，我们建立了“宏基因组”专业讨论群，目前己有国内外150+ PI，1300+ 一线科研人员加入。参与讨论，获得专业解答，欢迎分享此文至朋友圈，并扫码加主编好友带你入群，务必备注“姓名-单位-研究方向-职称/年级”。技术问题寻求帮助，首先阅读《如何优雅的提问》学习解决问题思路，仍末解决群内讨论，问题不私聊，帮助同行。

学习16S扩增子、宏基因组科研思路和分析实战，关注“宏基因组”