雙向富集分析圖

1.需求

上下調(diào)基因分別富集分析褪储，把結(jié)果放在一張圖上:

2.示例數(shù)據(jù)和r包準(zhǔn)備

library(ggplot2)
library(dplyr)
set.seed(2012)
data <- data.frame(
  id = paste0("ca",1:20), 
  p = runif(20,0,0.05),
  change = sample(c("up","down"),20,replace = T)
)
head(data)

##    id          p change
## 1 ca1 0.01091543     up
## 2 ca2 0.03720259   down
## 3 ca3 0.01408370   down
## 4 ca4 0.04675872     up
## 5 ca5 0.03732084   down
## 6 ca6 0.04373917     up

模擬的是kegg的富集結(jié)果，編了一個數(shù)據(jù)慧域。p是p值鲤竹。畫圖的時候用-log10(p)。

因為我們需要上下調(diào)分兩個方向昔榴，所以先把下調(diào)的整成負(fù)值log10(p)辛藻。

為了讓通路按照上下調(diào)和p值大小順序排列起來，還需要用有序因子控制一下橫坐標(biāo)的順序互订。

3.實現(xiàn)雙向畫圖

data <- data %>% 
  mutate(p2 = ifelse(change == "up", -log10(p), log10(p))) %>% 
  arrange(change,p2) 
data$id = factor(data$id,levels = unique(data$id),ordered = T)
head(data)

##     id           p change        p2
## 1 ca16 0.002428496   down -2.614663
## 2 ca19 0.005513276   down -2.258590
## 3 ca18 0.007077920   down -2.150094
## 4  ca3 0.014083697   down -1.851283
## 5 ca20 0.027069219   down -1.567524
## 6 ca15 0.032271096   down -1.491186

這里有一點繞：剛才把down對應(yīng)的那些通路-log10（p）值設(shè)置為負(fù)了吱肌，需要在圖上改回去，直接改坐標(biāo)即可仰禽。

以前都是手動修改breaks和labels參數(shù)氮墨，今天我突然想到，ggplot2的默認(rèn)的breaks是怎么生成的呢吐葵？一番搜索找到了答案规揪，茅塞頓開：https://stackoverflow.com/questions/38486102/how-does-ggplot-calculate-its-default-breaks

tmp = with(data, labeling::extended(range(p2)[1], range(p2)[2], m = 5));tmp

## [1] -3 -2 -1  0  1  2

lm = tmp[c(1,length(tmp))];lm

## [1] -3  2

這里的tmp就是默認(rèn)的breaks。lm是坐標(biāo)范圍温峭，也可作為參數(shù)設(shè)置進(jìn)去粒褒。

有兩種辦法可以畫出來雙向的直方圖

方法1 ：geom_segment

ggplot(data, aes(x=id, y=p2)) +
  geom_segment( aes(x=id, xend=id, y=0, yend=p2, color=change), size=5, alpha=0.9) +
  theme_light() +
  theme(
    panel.border = element_blank()
  ) +
  xlab("") +
  ylab("-log10(pvalue)")+
  ylim(lm)+
  scale_y_continuous(breaks = tmp,
                     labels = abs(tmp))+
  coord_flip()

方法2 ：geom_bar

ggplot(data, aes(x=id, y=p2)) + 
  geom_bar(stat='identity', aes(fill=change), width=.7)+
  coord_flip()+
  theme_light() +
  ylim(lm)+
  scale_y_continuous(breaks = tmp,
                     labels = abs(tmp))+
  theme(
    legend.position = "none",
    panel.border = element_blank()
  ) 

4.搞定kegg和go雙向富集圖的函數(shù)

氛圍到了，寫著寫著就安排上了唄诚镰。拿去用_{1.4.0以上版本的tinyarray可用奕坟，今天（2020.10.16）剛更新上去的。輸出格式是ggplot2清笨，所以ggplot2細(xì)節(jié)調(diào)整月杉，都可以添加到這張圖上去}

load("step4output.Rdata")
library(tinyarray)
pt = double_enrich(deg,n = 5)
pt

## $kp

## 
## $gp