こんにちは。ほしのはやしです。
生物を対象にした研究をする人は見たことがあるVolcano plotについて、初心者向けの解説とRを用いた書き方について説明します!

Volcano plotとは?
Volcano Plot(火山プロット)は、生物学や医学研究、特に遺伝子発現解析やプロテオミクス解析で広く使用される可視化手法の一つです。
このプロットは、【統計的有意性(p値)と変動の大きさ(fold change)】を同時に示すことで、実験群と対照群の間で顕著に変化している遺伝子やタンパク質を素早く特定するのに役立ちます。
Volcano plotの基本
Volcano plotを書くためには、以下の情報が必要です。
Fold Change(変動の大きさ)
⇒実験群と対照群の発現量の比を表す指標
p値(統計的有意差)
⇒実験群と対照群の比較において、発現量の差が偶然によるものかどうかを示す指標
遺伝子やタンパク質の名前
⇒どのデータに着目するかに役立つ
例えばRNA-seq解析では、リードカウントのデータからDESeq2などを用いて以下のようなデータシートが作成されます。

基本的には、大量にデータを解析することになるので、多重性の問題がでてきます。
そのためp値は【adj.P.Val】など補正されたp値を使うことになります。
多重性についての詳しい説明は以下のページをご覧ください!

Volcano plotの構造
Volcano Plot は散布図の一種であり、以下のような軸で描画されます。
x軸(横軸):Log2 Fold Change (log2FC)
実験群と対照群の比較において、遺伝子やタンパク質の発現量の変化の大きさを示します。
y軸(縦軸):-log10(p値)
p値が小さいほどその変化が有意であることを意味します。
-log10(p値) を取ることで、p値が小さいほど(統計的に有意であるほど)上の方にプロットされます。
しきい値の設定
有意な変動を示すデータポイント(遺伝子など)は、特定のしきい値を超えるものとしてハイライトしてわかりやすくします。
【例】
統計的有意性のしきい値:p値 < 0.05(-log10(0.05) ≈ 1.3)
fold change のしきい値:|log2FC| > 1(2倍以上の変動)
Volcano Plotの用途
Volcano Plotは、主に以下のようなデータ解析に用いられます。
RNA-seq やマイクロアレイデータの解析
RNA-seq解析では、遺伝子の発現量を比較し、どの遺伝子が統計的に有意な変化を示しているかを調べるのに使われます。例えば、がん細胞と正常細胞の遺伝子発現を比較し、発現量が大きく変化しているがん関連遺伝子を特定することができます。
プロテオミクスデータの解析
質量分析 (Mass spectrometry) を用いたプロテオミクス解析において、タンパク質の発現量を比較し、疾患特異的なバイオマーカーの特定に役立ちます。
メタボローム解析
代謝産物の比較分析において、異なる条件間で有意に変化する代謝物を特定するために使用されます。
CRISPRスクリーンのデータ解析
CRISPR-Cas9を用いたスクリーニングデータの解析にも使用され、特定の遺伝子ノックダウンの影響を評価するのに役立ちます。
Volcano plotの解釈
基本的には、以下の4つの視点で概略を掴むことが大事です。
右上の点(高いlog2FCかつ有意なp値)
→ 実験群で顕著に発現上昇した遺伝子(upregulated genes)
左上の点(低いlog2FCかつ有意なp値)
→ 実験群で顕著に発現低下した遺伝子(downregulated genes)
中央の点(log2FCが0付近)
→ 変動が少なく、有意ではない遺伝子
x軸付近の点(p値が大きい = -log10(p値) が小さい)
→ 統計的に有意でない変動

ぱっとみてどの遺伝子/タンパクが動いてるか直観を得るのに有用なのだ!
R studioを使った実際の書き方
まずは必要な項目だけに絞った練習用のデータセットを作成します。
# 必要なパッケージの読み込み
library(tidyverse)
# サンプルデータの作成
set.seed(42)
num_genes <- 500 # 遺伝子数を増やす
data <- data.frame(
gene_name = paste0("Gene", 1:num_genes),
logFC = rnorm(num_genes, mean = 0, sd = 1.5), # より現実的な変化
adj.P.Value = runif(num_genes, min = 0.0001, max = 0.2) # 0.05以下の有意な値を含める
)


続いて、p値をlog10で対数変換して、視覚的に見やすくするためのFold change, p値のしきい値を決めます。
またその遺伝子のみ表示するようにします。
# p値の-log10変換
data$neg_log10_pval <- -log10(data$adj.P.Value)
# ラベルを付ける遺伝子のフィルタリング(|logFC| > 1.5 かつ adj.P.Value < 0.05)
data$label <- ifelse(data$adj.P.Value < 0.05 & abs(data$logFC) > 1.5, data$gene_name, "")


今回は、LogFCが1.5より大きいかつP値<0.05をラベルするようにしました!
ご自身のデータに応じてしきい値は修正してください!
# 必要なパッケージの読み込み(ラベルの重なりを防ぐパッケージ)
library(ggrepel)
# Volcano plotの作成
ggplot(data, aes(x = logFC, y = neg_log10_pval)) +
geom_point(aes(color = adj.P.Value < 0.05), size = 2, alpha = 0.6) + # 有意な遺伝子を赤色
scale_color_manual(values = c("gray", "red")) +
geom_hline(yintercept = -log10(0.05), linetype = "dashed", color = "blue") + # 有意水準線
geom_vline(xintercept = 0, linetype = "dashed", color = "black") + # 中心線
geom_text_repel(aes(label = label), size = 4, max.overlaps = 15, box.padding = 0.5) + # ラベル追加
labs(title = "Volcano Plot", x = "log2 Fold Change", y = "-log10 Adjusted P-Value") +
theme_minimal()


これでVolcano plotを書くことができました!
別のパッケージを使ってよりビジュアル系に書く方法についてもご説明します!
「EnhancedVolcano」パッケージを以下のコードをコピペしてインストールします。
if (!requireNamespace('BiocManager', quietly = TRUE))
install.packages('BiocManager')
BiocManager::install('EnhancedVolcano')
途中にupdateするかと聞かれれば以下のように「a」と入力してアップデートしましょう。


先ほどのデータシートを使って以下のコードでグラフを書きます。
library(EnhancedVolcano)
EnhancedVolcano(data,
lab = data$gene_name,
x = 'logFC',
y = 'adj.P.Value',
title = 'Test data',
pCutoff = 10e-2,
FCcutoff = 0.5,
pointSize = 3.0,
labSize = 6.0)


開発者のページはこちらになりますので、細かい修正はこちらをご覧ください!
まとめ
Volcano plotの簡単な説明とRを使った作り方について解説しました!
少しでも皆様のお役に立てたなら幸いです!
コメント