統計学と疫学と時々、助教生活

疫学を専門とする助教の研究に関する備忘録的ページ。

R package 'qqman'の基本的な使用法

今回はRのqqmanというパッケージの基本的な使い方をまとめる。
【ページ下部に追記あり(20160922)】

概要

現在、ゲノム全体10万から50万ほどのSNPsと疾患との関連を探索するGWASが可能になっている。そこでその結果を可視化する際に用いるのがマンハッタンプロットである。マンハッタンプロットは横軸に染色体の番号を、縦軸に-log10P-valueを取る図である。このマンハッタンプロットを簡単に描出できるのがqqmanというRのパッケージである。

コード

library(qqman)
#qqmanのデモデータ(gwasResults)を表示
head(gwasResults)
#   SNP CHR BP         P
# 1 rs1   1  1 0.9148060
# 2 rs2   1  2 0.9370754
# 3 rs3   1  3 0.2861395
# 4 rs4   1  4 0.8304476
# 5 rs5   1  5 0.6417455
# 6 rs6   1  6 0.5190959
#早速、gwasRsultsの結果をmanhattan plot(デフォルト)として描画する。
manhattan(gwasResults)

f:id:ryosukefujii0320:20160118140219p:plain
赤色のラインはp値が5.0×10-8のレベルであり、GWASでの有意レベルにあることを示している。青色のラインはp値が1.0×10-5のレベルであり、suggestive level(すなわち示唆を与えるレベル=有意ではないが、なんらかの関連がありそうなグレーゾーン)であることを示している。
結果を見ると、Chromosome 3に幾つかGWAS Significant levelに到達しているSNPsを認めることがわかる。

以下、マンハッタンプロットの条件を変更する補足的なコードを参考に記す。

#プロットの大きさをcexで指定する。
manhattan(gwasResults, cex=0.8)

#染色体番号ごとに色を変える。
manhattan(gwasResults, col=c("blue4", "orange3"))

#suggestive levelの線を消す。
manhattan(gwasResults, suggestiveline=FALSE)

#先ほど全体のマンハッタンプロットで有意になった3番染色体のみを描画する。
chr3 <- subset(gwasResults, CHR==3)
manhattan(chr3) 

#興味のあるSNPsについて分かりやすく描画する。
head(snpsOfInterest)
manhattan(gwasResults, highlight=snpsOfInterest)

マンハッタンプロットに加えて、関連testのp値をQ-Qplotで描画することもできる。

qq(gwasResults$P, main="Q-Q plot of P-values")

f:id:ryosukefujii0320:20160118142300p:plain

まとめ

非常に簡単に使用できるパッケージであり、データの作成・成形をgwasResultsのようにできれば誰でも描出できる。

参考ページ

  1. qqman: an R package for visualizing GWAS results using Q-Q and manhattan plots. (http://biorxiv.org/content/biorxiv/early/2014/05/14/005165.full.pdf)
  2. Package'qqman'-CRAN- (https://cran.r-project.org/web/packages/qqman/qqman.pdf)

追記および注意点(2016/9/22)

  1. plinkのGWASファイルに関しては、そのままmanhattan(データ名)のみで描くことができる。つまり、上記のgwasResultsのように成型する必要性は一切ないことを確認。
  2. エラーメッセージとして、「plot.window(...) でエラー: 有限な 'ylim' の値が必要です」と出てきた場合には、Pの列が0から1の値以外のものがないか確認し、それを除外することで解決する(特にNAなど)。


20160118
RF