本期学习笔记的内容是在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
© 素材来源于网络,侵权请联系后台删除
往期推荐: