【原】火山圖|給你geneList，幫我標到火山圖上

生信補給站 2021-05-24

展開全文

火山圖（Volcano Plot）常用于展示基因表達差異的分布，橫坐標常為Fold change（倍數(shù)），越偏離中心差異倍數(shù)越大；縱坐標為P value（P值），值越大差異越顯著。得名原因也許是因為結果圖像火山吧!

火山圖只標示指定基因？這需求都遇到過吧。

一載入R包，數(shù)據(jù)

library(ggplot2)
library(openxlsx)
library(dplyr)

#繪制火山圖數(shù)據(jù)
data <- read.xlsx("火山圖.xlsx", sheet = 1)
head(data) #查看數(shù)據(jù)，主要有P值，F(xiàn)old change和基因ID即可。

二繪制火山圖（標示最顯著的基因）

2.1 先根據(jù)閾值分出上調和下調基因；

data$change <- as.factor(ifelse(data$adj.P.Val < 0.01 & abs(data$logFC) > 1,
ifelse(data$logFC > 1,'UP','DOWN'),'NOT'))

2.2 標示差異顯著的基因

data$sign <- ifelse(data$adj.P.Val < 0.001 & abs(data$logFC) > 2.5,data$GENE_SYMBOL,NA)head(data)

2.3 繪制火山圖

ggplot(data = data, aes(x = logFC, y = -log10(adj.P.Val), color = change)) +
geom_point(alpha=0.8, size = 1) +
theme_bw(base_size = 15) +
theme(panel.grid.minor = element_blank(),panel.grid.major = element_blank()) +
geom_hline(yintercept=2 ,linetype=4) +
geom_vline(xintercept=c(-1,1) ,linetype=4 ) +
scale_color_manual(name = "", values = c("red", "green", "black"), limits = c("UP", "DOWN", "NOT")) +
geom_text(aes(label = sign), size = 3)

了解一下ggplot2繪圖的方式，標示的基因就是各個基因的text，然后想辦法將其賦予到一個 aes 中即可。

三標示指定基因

和上面類似，將指定基因添加到繪制數(shù)據(jù)中即可。

3.1 讀入含有geneList的文件

gene <- read.xlsx("火山圖.xlsx", sheet = 2)
gene$geneList <- gene$gene

額外生成一列相同列是為了后面合并后還有一列存在，這一列用于標示基因。（方法有點笨）

3.2 合并火山圖數(shù)據(jù)

data2 <- data %>%

left_join(gene,by = c("GENE_SYMBOL" = "gene"))

head(data2)

增加了geneList列，為了后面使用text的方式添加上基因。

3.3 標示文件中的指定基因

ggplot(data = data2, aes(x = logFC, y = -log10(adj.P.Val), color = change)) +
geom_point(alpha=0.8, size = 1) +
theme_bw(base_size = 15) +
theme(panel.grid.minor = element_blank(),panel.grid.major = element_blank()) +
geom_hline(yintercept=2 ,linetype=4) +
geom_vline(xintercept=c(-1,1) ,linetype=4 ) +
scale_color_manual(name = "", values = c("red", "green", "black"), limits = c("UP", "DOWN", "NOT")) +
geom_text(aes(label = geneList), size = 5,color = "blue")

3.4 ggrepel 解決重疊問題

如果目標標示基因太多會導致重疊，可使用ggrepal函數(shù)

library(ggrepel)

ggplot(data = data2, aes(x = logFC, y = -log10(adj.P.Val), color = change)) +
geom_point(alpha=0.8, size = 1) +
theme_bw(base_size = 15) +
theme(panel.grid.minor = element_blank(),panel.grid.major = element_blank()) +
geom_hline(yintercept=2 ,linetype=4) +
geom_vline(xintercept=c(-1,1) ,linetype=4 ) +
scale_color_manual(name = "", values = c("red", "green", "black"), limits = c("UP", "DOWN", "NOT")) +
geom_label_repel(aes(label=geneList), fontface="bold", color="grey50", box.padding=unit(0.35, "lines"), point.padding=unit(0.5, "lines"), segment.colour = "grey50")