如何用ggplot2绘制基因功能富集气泡图?

简介: 如何用ggplot2绘制基因功能富集气泡图?

本期学习笔记的内容是在R语言中绘制功能富集泡泡图,最终能实现下面这种结果,包含差异表达基因数据的转换整理和绘图参数调整等内容,所有数据和代码已整理放在文末链接中。

学习绘图之前,首先学一下怎么看这张图:

  • 横轴表示差异基因占比,数值越大表示富集程度越大
  • 纵轴表示富集的通路,越靠上的通路越显著
  • 气泡的颜色表示pvalue的大小(越红表示pvalue的负对数越大,pvalue越小)
  • 气泡大小表示差异基因个数,越大则数量越多

功能富集分析泡泡图用来展示一组基因参与哪些功能调节通路,有利于理解基因的功能和潜在意义,这次跟着生信宝典的R语言学习教程,记录一下学习笔记。

准备阶段

富集矩阵数据

  • Description:通路描述
  • GeneRatio:该通路的差异基因占总差异基因比例
  • qvalue:富集的显著程度
  • Count:该通路差异基因的个数
  • Type:区分样本

原始数据

整理数据

对富集矩阵进行处理,读取为表格形式,设置分隔符和行名,并输出整理后的数据。

e_data <- read.table(text=data,header = T,row.names = NULL,
                     sep = "\t",quote="")
head(e_data)

整理后的数据

加载R包

library(plyr)
library(stringr)
library(grid)
library(ggplot2)

需要用到如上所述R包,其中ggplot2用于绘图。

单样本绘制

原始数据包含两个样本,这里先取出数据中的一个样本。

e_data_1=droplevels(e_data[e_data$Type=="Baodian_up",])

构造一个函数,使分数转化为小数。

mixedToFloat <- function(x){
  x <- sapply(x, as.character)
  is.integer  <- grepl("^-?\\d+$", x)
  is.fraction <- grepl("^-?\\d+\\/\\d+$", x)
  is.float <- grepl("^-?\\d+\\.\\d+$", x)
  is.mixed    <- grepl("^-?\\d+ \\d+\\/\\d+$", x)
  stopifnot(all(is.integer | is.fraction | is.float | is.mixed))
  numbers <- strsplit(x, "[ /]")
  ifelse(is.integer,  as.numeric(sapply(numbers, `[`, 1)),
         ifelse(is.float,    as.numeric(sapply(numbers, `[`, 1)),
                ifelse(is.fraction, as.numeric(sapply(numbers, `[`, 1)) /
                         as.numeric(sapply(numbers, `[`, 2)),
                       as.numeric(sapply(numbers, `[`, 1)) +
                         as.numeric(sapply(numbers, `[`, 2)) /
                         as.numeric(sapply(numbers, `[`, 3)))))
}

进行测试一下,看看刚才构造的函数有什么用处,可以发现该函数能将输入的分数转化为小数输出。

> mixedToFloat(c('2 3/4', '2/3', '11 1/4'))
[1]  2.75000000.666666711.2500000

转化数据

利用刚才构造的函数对数据进行转化,将基因的比例和计数转化为小数格式,并设置矩阵的列名;qvalue转换为负对数,并添加为新的一列。

> e_data_1$GeneRatio=mixedToFloat(e_data_1$GeneRatio)
> e_data_1$Count=mixedToFloat(e_data_1$Count)
> log_name <- paste0("negLog10_","qvalue")
> col_name_e_1 <- colnames(e_data_1)
> col_name_e_1 <- c(col_name_e_1,log_name)
> e_data_1$log_name <- log10(e_data_1$qvalue) * (-1)
> colnames(e_data_1) <- col_name_e_1

转换后的数据

排序统计

获取数据中通路的出现次数并进行排序,获得Description和ID信息。

> e_data_1_freq <- as.data.frame(table(e_data_1$Description))
> colnames(e_data_1_freq) <- c("Description","ID")
> head(e_data_1_freq)
                                       Description ID
1                            ERK1 and ERK2 cascade  1
2   negative regulation of Notch signaling pathway  1
3                         neuron apoptotic process  1
4             positive regulation of cell adhesion  1
5 positive regulation of cytoskeleton organization  1
6                        regulation of cell growth  1

将两个矩阵合并到一起,以Description列为基准,然后对合并后的矩阵排序(按照ID、GeneRatio、负对数值为序)

e_data_2 <- merge(e_data_1,e_data_1_freq,by="Description")
e_data_3 <- e_data_2[order(e_data_2$ID,
                           e_data_2$GeneRatio,
                           e_data_2$negLog10_qvalue),]

unique删除重复的Description行,然后转换成因子,并进行排序。

t_order <- unique(e_data_3$Description)
e_data_1$Description <- factor(e_data_1$Description,
                               levels = t_order,ordered = T)

绘图阶段

设置颜色

color_1 <- c("green","red")

初始化坐标系

设置x轴和y轴,增加标签x="GeneRatio",y="GO description"

p <- ggplot(e_data_1,aes(x=GeneRatio,y=Description)) +
  labs(x="GeneRatio",y="GO description") + labs(title="")

设置点状图的大小和颜色信息,scale_color_gradient用于生成一组渐变色。

p <- p + geom_point(aes(size=Count,color=negLog10_qvalue)) +
  scale_color_gradient(low = color_1[1],high=color_1[2],name="negLog_qvalue")

对纵轴的标签进行处理,str_wrap用于设置字符串单行长度不超过60.

p <- p + scale_y_discrete(labels=function(x) str_wrap(x,width = 60))

当前结果图像

设置绘图区域参数

初始化参数如下:

top="top"
bottom="bottom"
left="left"
right="right"
none="none"
legend_pos_par=right

使用判断语句来自动估算图形的长宽数据,设置平面大小。

uwid=0
vhig=12
if (uwid == 0 || vhig == 0) {
  x_len = length(unique(e_data_1$Description))
  if(x_len<10){
    vhig = 10
  } elseif(x_len<20) {
    vhig = 10 + (x_len-10)/3
  } elseif(x_len<100) {
    vhig = 13 + (x_len-20)/5
  } else {
    vhig = 40
  }
  uwid = vhig 
  if(legend_pos_par %in% c("left", "right")){
    uwid = 1.5 * uwid
  }
}

调整气泡图

使用如下代码调整图像,步骤:设置图例的位置——设置边界和背景——设置坐标轴参数——输出plot。

p <- p + theme(legend.position = legend_pos_par)
p <- p + theme(panel.grid=element_blank(),
               panel.border = element_blank(),
               legend.background = element_blank(),
               axis.line.x=element_line(size=0.4,colour = "black",linetype = "solid"),
               axis.line.y=element_line(size=0.4,colour = "black",linetype = "solid"),
               axis.ticks=element_line(size=0.4))
p

最终结果图片

参考资料:

http://www.ehbio.com/Bioinfo_R_course

END

© 素材来源于网络,侵权请联系后台删除

往期推荐:

文献丨群体转录组分析锁定关键转录因子

文献丨转录组RNA seq——青年阶段!

笔记丨ggplot2热图入门学习笔记

笔记丨PCA分析基本知识和数学原理

相关文章
|
13天前
|
数据可视化
绘制GGPLOT2双色XY区间面积图组合交叉折线图数据可视化
绘制GGPLOT2双色XY区间面积图组合交叉折线图数据可视化
|
15天前
R语言中绘制箱形图的替代品:蜂群图和小提琴图
R语言中绘制箱形图的替代品:蜂群图和小提琴图
21 0
|
5月前
|
存储 数据可视化
使用 plotly 绘制旭日图
使用 plotly 绘制旭日图
74 0
|
8月前
|
C++ Python
Python绘制比例直方柱状比重图
Python绘制比例直方柱状比重图
145 1
Python绘制比例直方柱状比重图
|
9月前
|
数据可视化 数据挖掘 Linux
数据可视化丨优雅的绘制带显著性标记的箱线散点图,主要使用ggsignif和ggplot2
数据可视化丨优雅的绘制带显著性标记的箱线散点图,主要使用ggsignif和ggplot2
|
11月前
|
数据挖掘
ggplot2| 绘制KEGG气泡图
ggplot2| 绘制KEGG气泡图
201 0
|
11月前
|
数据可视化
ComplexHeatmap|根据excel表绘制突变景观图(oncoplot)
ComplexHeatmap|根据excel表绘制突变景观图(oncoplot)
157 0
|
11月前
R绘制多彩的森林图:基于ggplot2
R绘制多彩的森林图:基于ggplot2
142 1
|
11月前
|
人工智能 数据可视化
跟SCI学umap图| ggplot2 绘制umap图,坐标位置 ,颜色 ,大小还不是你说了算
跟SCI学umap图| ggplot2 绘制umap图,坐标位置 ,颜色 ,大小还不是你说了算
725 1
|
12月前
|
机器学习/深度学习 算法 数据可视化
R绘图案例|基于分面的面积图绘制
R绘图案例|基于分面的面积图绘制
10830 0