【R】ggplot2メモ（volcano plot）

・お題：ggplot2でvolcanoplotを良い感じに描くコマンドをメモしておきたい。

・EnhancedVolcanoでいいじゃないかといわれるとそうなのだけれど、ggplot2で描きたい。

・絶妙にうまいこと自分でカスタムしたいなぁと思っていたところ、以下の記事を見つけた。とても分かりやすいので、気になる方はぜひそちらをご確認いただきたい。

erikaduan.github.io

・せっかくなので、やってみる。

if (!require("pacman")) install.packages("pacman")
pacman::p_load(here,
tidyverse,
janitor,
scales,
ggrepel)

・最初この時点で何をやっているのか分からなかった。pacmanというのは、パッケージをうまいことロードできるようにするパッケージらしい。p_loadでロードすると（githubはどうか知らないが）既にインストール済のライブラリは普通にロードして、インストールしていないライブラリは勝手にインストールしてくれるらしい。しかもまとめてロードできる。このコマンドさえ打てば、初めてのPCでも勝手に環境が構築されるということらしい。ちなみに、hereはディレクトリ、janitorはデータフレームのネーミング修正、scalesはlog軸取るのに便利、ggrepelはggplotの描画に注釈をつけるのに使うっぽい。

・データを読み込む。今回は参考サイトの著者が提供しているデータをそのままカレントディレクトリに入れた。そのまま列名を修正する。

diseased_vs_healthy <-
read.delim("limma-voom_luminalpregnant-luminallactate.txt",
header = TRUE,
sep = "\t") %>%
janitor::clean_names()

・次に、少しデータを弄って解析に使いやすいようにする。

diseased_vs_healthy <- diseased_vs_healthy %>%
mutate(fold_change = 2^log_fc) %>%
select(entrezid,
symbol,
fold_change,
adj_p_val)

・DEGを図示する際に基準になる判定列を作成。

diseased_vs_healthy <- diseased_vs_healthy %>%
mutate(gene_type = case_when(fold_change >= 2 & adj_p_val <= 0.05 ~ "up",
fold_change <= 0.5 & adj_p_val <= 0.05 ~ "down",
TRUE ~ "ns"))

・ここまででデータ本体のしら準備は終了。

・カラーセットを作成する。好みの問題で、少し色を変えた。

cols <- c("up" = "orange", "down" = "blue", "ns" = "grey")

・ラベル付けで強調したい遺伝子を指定する。
sig_il_genes <- diseased_vs_healthy %>%
filter(symbol %in% c("Il15", "Il34", "Il24"))

up_il_genes <- diseased_vs_healthy %>%
filter(symbol == "Il24")

down_il_genes <- diseased_vs_healthy %>%
filter(symbol %in% c("Il15", "Il34"))

・あとはひたすら図示。やっていることは散布図描いて補助線引いたり色付けたり注釈付けたりしているのだけれど、とてもggplotらしくて面白いと感じた。

final_plot <-
#データ読み込み
ggplot(data = diseased_vs_healthy,
aes(x = log2(fold_change),
y = -log10(adj_p_val))) +
#基本のプロット
geom_point(aes(colour = gene_type),
alpha = 0.4,
shape = 16,
size = 1) +
#ラベル付けしたいUp遺伝子の点
geom_point(data = up_il_genes,
shape = 21,
size = 2,
fill = "orange",
colour = "black") +
#ラベル付けしたいDown遺伝子の点
geom_point(data = down_il_genes,
shape = 21,
size = 2,
fill = "blue",
colour = "black") +
#横の破線
geom_hline(yintercept = -log10(0.05),
linetype = "dashed") +
#縦の破線
geom_vline(xintercept = c(log2(0.5), log2(2)),
linetype = "dashed") +
#ラベル付け
geom_label_repel(data = sig_il_genes,
aes(label = symbol),
force = 2,
nudge_y = 1) +
#色を指定する
scale_colour_manual(values = cols) +
#X軸を指定する。
scale_x_continuous(breaks = c(seq(-10, 10, 2)),
limits = c(-10, 10)) +
#ラベルを指定する。
labs(title = "Gene expression changes in diseased versus healthy samples",
x = "log2(fold change)",
y = "-log10(adjusted P-value)",
colour = "Expression \nchange") +
#テーマを設定する。
theme_bw() + # Select theme with a white background
theme(panel.border = element_rect(colour = "black", fill = NA, linewidth= 0.5),
panel.grid.minor = element_blank(),
panel.grid.major = element_blank())