ggplot2：如何优雅地绘制箱线图

2019-04-01 本文已影响6人 MC学公卫

数据可视化是将数据直观展示出来的一个必不可少的步骤。ggplot2是一个非常常见的绘图R包，在文献中，我们常常会看见用漂亮箱线图来展示原始数据的分布。这里我们就用一个基因表达水平的例子来展示用R包绘图神器ggplot2绘制的过程吧！

Step1：数据的预处理

这里我们利用GEO中单细胞RNA测序数据来练练手。
数据来源：
https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE128147

在R中读取这个数据：

counts = read.table("~/Downloads/GSE128147_raw_counts_w_ercc.txt", header = TRUE)

我们可以得到行是基因名，列是细胞名，中间的数值为表达量。也就是每个基因在每个细胞中的表达量，这种数据称为raw counts。如图所示：

数据展示

dim(counts)
[1] 23512   193

这里数据一共有193个细胞，通常我们有某种研究目的时，需要看某几类细胞的表达情况，那么根据细胞名或其他指标将它们区分开。

由于这里仅为了学习如何作图，那么为了方便，只取几种基因，并将这些细胞按照顺序粗暴地分为3组，每组64个，每组分别命名为CellType1, CellType2, CellType3，计算每个基因在各类细胞中的表达水平（求和），并整理成一个数据框，变量名为expres。整理代码如下：

## 将原始数据截取出左、中、右三段
type1 = counts[,1:65]
type2 = counts[,c(1,66:129)]
type3 = counts[,c(1,130:193)]

## 将每个基因在每种细胞内的表达量求和，并将结果添加到数据框中
type1 = within(type1,{sum.of.expression = rowSums(type1[-1])})
type2 = within(type2,{sum.of.expression = rowSums(type2[-1])})
type3 = within(type3,{sum.of.expression = rowSums(type3[-1])})

## 只取出第1列基因名与最后一列表达水平之和
type1 = type1[,c(1,66)]
type2 = type2[,c(1,66)]
type3 = type3[,c(1,66)]

## 将三个向量合并为一个数据框
expres = data.frame(Gene = type1[,1], 
                    CellType1 = type1[,2], 
                    CellType2 = type2[,2], 
                    CellType3 = type3[,2])

这么一来得到的数据就是这样：

数据展示2

> head(expres)
    Gene CellType1 CellType2 CellType3
1   Xkr4         0         0         0
2    Rp1         0      1249         4
3  Sox17         0         0         0
4 Mrpl15        12      2172      1483
5 Lypla1       268      1956      2235
6  Tcea1       206       508      1443

将数据变形，利于后续作图。

library(reshape2)
expres.melt = melt(expres, value.name = "Counts")
# Using Gene as id variables
# ggplot2常用melt型数据
expres.melt = subset(expres.melt, Counts > 300 & Counts < 10000) 
# 筛选掉一些表达量过少和过多的基因

expres.melt$Counts = log(expres.melt$Counts) 
#由于表达量绝对值差距太大，因此用对数将数据标准化

head(expres.melt)
      Gene  variable   Counts
8  Atp6v1h CellType1 6.967909
11  Rb1cc1 CellType1 7.502186
14  Pcmtd1 CellType1 5.910797
16    Rrs1 CellType1 7.092574
21  Vcpip1 CellType1 8.051978
30 Arfgef1 CellType1 7.686621

Step2：绘制箱线图

1. 画一个基础的箱线图：geom_boxplot()

library(ggplot2)
ggplot(expres.melt, aes(x=variable,y=Counts)) +     
  geom_boxplot(aes(fill=variable))

expres.melt即我们的数据；
aes(x = 横坐标向量, y = 纵坐标向量)；
geom_boxplot 箱线图；
fill=variable 按照variable类型来填充颜色。

这么一来就得到一个还不错的箱线图：

箱线图1

2. 更改坐标轴标题、图片标题与图例标题：labs()

在刚才的基础上我们想要将图的横坐标、纵坐标以及图的标题进行修改，那么可以使用labs()函数在之前的代码后面进行添加。

ggplot(expres.melt, aes(x=variable,y=Counts)) +     
  geom_boxplot(aes(fill=variable)) +
  labs(title="Expression Level of 3 Types of Cells",
       x="Cell Type", y = "Counts/(log)", fill = "Cell Type")  # fill为修改图例标题

箱线图2

3. 图例的设置：theme()

图例标题的修改比较特殊，不能再用labs了，而是用theme()，同时theme()这个函数还能设置图例的标题的字体、颜色、大小。也能修改总标题的位置。

ggplot(expres.melt, aes(x=variable,y=Counts)) +     
  geom_boxplot(aes(fill=variable)) +
  labs(title="Expression Level of 3 Types of Cells",
       x="Cell Type", y = "Counts/(log)", fill = "Cell Type") +
  theme(plot.title = element_text(hjust = 0.5), # 将图表标题居中
        legend.title=element_text(face="italic",  # 图例标题改为斜体
                                  family="Times",  # 图例标题字体调为Times
                                  colour="red"))  # 图例标题颜色改为红色

箱线图3.jpeg

4.将箱线图转置：coord_flip()

有时候如果数据类型比较多，为了排版方便，可能需要将箱线图转置，这种情况下怎么处理呢？

ggplot(expres.melt, aes(x=variable,y=Counts)) +     
  geom_boxplot(aes(fill=variable)) +
  labs(title="Expression Level of 3 Types of Cells",
       x="Cell Type", y = "Counts/(log)",
       fill = "Cell Type") +
  theme(plot.title = element_text(hjust = 0.5) +
  coord_flip()

箱线图4.jpeg

图形就横过来了。

5.将背景变为白色：theme_bw()

背景的颜色也可以进行修改，这样一来做PPT上组会汇报时白色背景显得更为学术和美观。

ggplot(expres.melt, aes(x=variable,y=Counts)) +     
  geom_boxplot(aes(fill=variable)) +
  labs(title="Expression Level of 3 Types of Cells",
       x="Cell Type", y = "Counts/(log)",
       fill = "Cell Type") +
  theme_bw()

箱线图5.jpeg

除了这些，其实还有许多的功能可以叠加，如增加数据点可叠加散点+geom_dotplot，标注异常值可利用outlier，代码示例：

ggplot(expres.melt, aes(x=variable,y=Counts)) +     
  geom_boxplot(aes(fill=variable),
               outlier.colour="red", 
               outlier.shape=8, 
               outlier.size=4)

但由于数据本身的原因不含异常值，因此显示不出来，代码贴出来供学习使用。

文章已发布到微信公众号：百味科研芝士，欢迎关注。