ggplot2绘制火山图

作者: 路人里的路人 | 来源:发表于2023-08-26 11:01 被阅读0次

可能需要写在前面的东西

gene_exp <- read.table("D:/bioinfotools/R/R/R-4.3.1/R-4.3.1/wkdir/rna-seq/genes.TMM.EXPR.matrix",
                       header = T, row.names = 1)
#加载表达矩阵

library(readr)
library(tidyverse)
gene_info <- read_delim("D:/share/R_data/data/rnaseq-apple/query_seqs.fa.emapper.annotations", 
                                    delim = "\t", escape_double = FALSE, 
                                    col_names = FALSE, comment = "#", trim_ws = TRUE) %>%
  select(Gene_Id = X1,
         Gene_Symbol = X6,
         GO = X7,
         Ko = X9,
         Pathway = X10,
         COG = X21,
         Gene_Name = X22,)
#加载基因信息表

de_result <- read.table(file = 'D:/bioinfotools/R/R/R-4.3.1/R-4.3.1/wkdir/rna-seq/genes.counts.matrix.KID_S1_vs_BLO_S1.DESeq2.DE_results',
                 header = T)
#加载差异表达数据

de_result <- mutate(de_result, direction = if_else(padj > 0.05, 'ns', 
            if_else(abs(log2FoldChange) < 1, 'ns', 
                    if_else(log2FoldChange >= 1, 'up', 'down')))) %>%
left_join(gene_info, by = c('id' = 'Gene_Id')) %>%
left_join(rownames_to_column(gene_exp, var = 'id'), by = 'id') %>%
dplyr::select(-c(2:6, 8:9)) %>%
arrange(desc(abs(log2FoldChange)))
#合并数据变量

library(tidyverse)
de_result <- mutate(de_result,
                    direction = factor(direction, levels = c('up', 'ns', 'down')))
#对差异表达数据的direction列进行一个离散型变量生成的过程

de_top <- filter(de_result, abs(log2FoldChange) > 2 & padj < 1e-5)
#de_top <- filter(de_result, id %in% c('HF04388', 'HF01838'))
#对想要着重标记的基因进行一个筛选

快速入门

导入差异表达数据

de_result <- read.table(file = 'D:/bioinfotools/R/R/R-4.3.1/R-4.3.1/wkdir/rna-seq/genes.counts.matrix.KID_S1_vs_BLO_S1.DESeq2.DE_results',
                 header = T)

创建一个背景画布ggplot()，并将x轴（log2FoldChange）和y轴（-log10(padj)）映射到画布上

library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj)))

添加散点几何对象geom_point()并添加背景主题theme_void()，几何对象可以是多种多样的，如散点几何对象geom_point()、棒状几何对象geom_bar(),根据自己的需要可以选择不同的几何对象，敲出geom后选择即可。背景主题也是一样的道理。

library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point() +
  theme_bw()

将direction这一列映射给点的颜色，点的颜色不属于背景画布ggplot()而属于散点几何对象图层geom_point()

library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction)) +
  theme_bw()

标度控制着映射的方式,这里添加的标度就是调色板scale_color_frontiers()，scale_color_frontiers()是别人已经设置好的对应的期刊的配色特征。

library(ggplot2)
library(ggsci)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction)) +
  scale_color_frontiers() +
  theme_bw()

ggplot2所涉及的思想
1.图层：ggplot2是由一层一层的图层构建起来的，首先建立一背景图图层，再建立一个几何对象，再建立主题颜色
2.映射：aes()可以将数值映射到x，y轴，也可以将颜色映射到点上
3.标度：控制映射的规律
4.主题：图片好不好看

深入了解

修改点的颜色

mutate()函数更改direction列的信息,factor(direction, levels = c('up', 'ns', 'down')))表示对direction列重新进行定义，新的列有因子，按照levels中的顺序进行排列使direction列变成了有顺序的离散型变量，数据映射的顺序与因子的顺序是一致的

library(tidyverse)
de_result <- mutate(de_result,
                    direction = factor(direction, levels = c('up', 'ns', 'down')))

在这里，我们使用的颜色标度进行了更改，不再是系统提供的，而是自己设定的。首先自定义一个调色板my_palette，在调色板上放上自己喜欢的颜色，再将标度控制映射的方式进行更改（我这里改成了scale_color_manual(),即手动调控）

my_palette <- c('#7FFFD4', '#DCDCDC', '#FA8072')
library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction)) +
  scale_color_manual(values = my_palette) +
  theme_bw()

修改点的透明度

点的透明度位于geom_point()图层中，通过映射alpha来控制点的透明度，透明度的高低仍由abs(log2FoldChange)来控制

my_palette <- c('#7FFFD4', '#DCDCDC', '#FA8072')
library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, 
                 alpha = abs(log2FoldChange))) +
  scale_color_manual(values = my_palette) +
  theme_bw()

修改点的大小

点的大小仍然位于geom_point（）这个图层中，为了保证上调和下调的点大小一致，所以取的是log2FoldChange的绝对值并通过scale_size（）这个标度控制了点的大小范围

my_palette <- c('#7FFFD4', '#DCDCDC', '#FA8072')
library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, 
                 size = abs(log2FoldChange))) +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  theme_bw()

添加阈值线

在ggplot画图的过程中，将所有的几何对象geom()放在一起，所有的标度scale()放一起，所有的主题theme()放一起。
首先添加的是水平线，阈值为-log10(0.05),水平线的几何对象是geom_hline(),yintercept = -log10(0.05)表示阈值为-log10(0.05)，linetype = 'dashed'表示填充的线条类型是虚线，线条类型包括 ‘solid‘：实线 ’dotted‘: 点线 ’dotdash‘: 点划线 ’longdash‘: 长虚线 ’twodash‘: 双虚线
yintercept表示y轴，xintercept表示x轴。
再添加的是竖直线，阈值为（-1, 1）竖直线的几何对象是geom_vline(),xintercept = c(-1, 1)表示虚线的取值范围是-1到1。
另外还需要明白一件事，就是几何对象中的两种操作模式：映射与直接设置，直接设置是对所有对象执行相同的操作，映射就是对每个点进行单独的设置

my_palette <- c('#7FFFD4', '#DCDCDC', '#FA8072')
library(ggplot2)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, 
                 size = abs(log2FoldChange))) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  theme_bw()

添加标签

对数据进行筛选，因为最终的图片我只需要显示我想要的几个基因或筛选后的基因，那么就要设置条件进行筛选，filter(de_result, abs(log2FoldChange) > 2 & padj < 1e-10)为设置条件进行筛选，将log2FoldChange绝对值大于2并且padj < 1e-10的值筛选出来。
filter(de_result, id %in% c('HF04388', 'HF01838'))为指定想要显示的几个基因，这里制定了两个想要显示的基因：'HF04388', 'HF01838'。如果只要显示一个基因则为filter(de_result, id == 'HF04388')

library(tidyverse)
de_top <- filter(de_result, abs(log2FoldChange) > 2 & padj < 1e-5)

de_top <- filter(de_result, id %in% c('HF04388', 'HF01838'))

对图中某些点打标签是通过映射完成的，需要使用的包是ggrepel，geom_label_repel(data = de_top, aes(label = id))实现了对图中特定的点打标签，data = de_top表示加载的数据是前面筛选后的数据，aes(label = id)表示映射到图中的标签是de_top中的id列

my_palette <- c('#7FFFD4', '#DCDCDC', '#FA8072')
library(ggplot2)
library(ggrepel)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, size = abs(log2FoldChange))) +
  geom_label_repel(data = de_top, aes(label = id)) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  theme_bw()

设置坐标轴的名称

坐标轴名称主要是通过labs()来实现更改的，x,y,title,size分别表示x轴名称,y轴名称，图片标题，图例标题

my_palette <- c('#FA8072', '#DCDCDC', '#7FFFD4')
library(ggplot2)
library(ggrepel)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, size = abs(log2FoldChange))) +
  geom_label_repel(data = de_top, aes(label = id)) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  labs(x = 'log2 fold change', 
       y = '-log10(pvalue)',
       title = 'volcano plot',
       size = 'log2 fold change') +
  theme_bw()

设置坐标轴范围

y轴的范围是由ylim()进行控制的，同理也可以通过xlim()控制x轴的坐标轴范围

my_palette <- c('#FA8072', '#DCDCDC', '#7FFFD4')
library(ggplot2)
library(ggrepel)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, size = abs(log2FoldChange))) +
  geom_label_repel(data = de_top, aes(label = id)) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  labs(x = 'log2 fold change', 
       y = '-log10(pvalue)',
       title = 'volcano plot',
       size = 'log2 fold change') +
  ylim(c(0, 15)) +
  theme_bw()

标题居中

由于使用的主题是theme_bw()，该主题默认情况下就是标题靠左，需要对其进行一些参数调整才可以达到目的
通过theme()函数实现目标，plot.title = element_text表示标题是一串文字，element_text(size = 18, hjust = 0.5)表示标题的字体大小为18，对齐方式未居中，0为左对齐，0.5为居中对齐，1为右对齐，hjust表示水平方向，vjust表示竖直方向

my_palette <- c('#FA8072', '#DCDCDC', '#7FFFD4')
library(ggplot2)
library(ggrepel)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, size = abs(log2FoldChange))) +
  geom_label_repel(data = de_top, aes(label = id)) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  labs(x = 'log2 fold change', 
       y = '-log10(pvalue)',
       title = 'volcano plot',
       size = 'log2 fold change') +
  #ylim(c(0, 15)) +
  theme_bw() +
  theme(plot.title = element_text(size = 18, hjust = 0.5))

去除图例

去除一个图例使用的函数是guides(),使映射size的的值为FALSE就可以实现对图例的取消

my_palette <- c('#FA8072', '#DCDCDC', '#7FFFD4')
library(ggplot2)
library(ggrepel)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, size = abs(log2FoldChange))) +
  geom_label_repel(data = de_top, aes(label = id)) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  labs(x = 'log2 fold change', 
       y = '-log10(pvalue)',
       title = 'volcano plot',
       size = 'log2 fold change') +
  #ylim(c(0, 15)) +
  guides(size = FALSE) +
  theme_bw() +
  theme(plot.title = element_text(size = 18, hjust = 0.5))

修改图例位置

修改图例位置也是通过theme()函数实现的，legend.position表示图例，c(0.93, 0.8)表示它在图层中所处的位置这个图层是画布图层，其可以看作是第一二象限，横坐标从0-5.0的长度是1，从0--5.0的长度是1，纵坐标从0-20的长度是1。故要在右上角大概坐标就为（0.9，0.85）.
此外就是图例的背景颜色需要修改，使用legend.background进行修改，使用的函数是element_blank

my_palette <- c('#FA8072', '#DCDCDC', '#7FFFD4')
library(ggplot2)
library(ggrepel)
ggplot(data = de_result, aes(x = log2FoldChange, y = -log10(padj))) +
  geom_point(aes(color = direction, size = abs(log2FoldChange))) +
  geom_label_repel(data = de_top, aes(label = id)) +
  geom_hline(yintercept = -log10(0.05), linetype = 'dashed') +
  geom_vline(xintercept = c(-1, 1), linetype = 'dashed') +
  scale_color_manual(values = my_palette) +
  scale_size(range = c(0.1,3)) +
  labs(x = 'log2 fold change', 
       y = '-log10(pvalue)',
       title = 'volcano plot',
       size = 'log2 fold change') +
  #ylim(c(0, 15)) +
  guides(size = FALSE) +
  theme_bw() +
  theme(plot.title = element_text(size = 18, hjust = 0.5), 
        legend.background = element_blank(), 
        legend.position = c(0.93, 0.80))