アレイデータとNGSデータをIGVで並べたい
最近、メチル化アレイのデータをNGSデータと並べて見てみたいから、IGVで可視化できないかと相談を受けました。
日頃、様々なツール間でデータをやり取りするためにファイルを変換するスクリプトをよく書いているので、つい発想が「じゃあアレイデータをbedgraphへ変換するスクリプトでも書こうかな」となってしまうのですが、よく考えたら、アレイデータをNGSデータと並べて見たい、というのは、多くの人が考えそうなこと。
つまり、誰かが先にやっているのでは……？
先人の知恵探し
ひとまず検索をかけて見ると、某バイオインフォマティシャン御用達の掲示板で、「メチル化アレイをIGVで可視化したい」という、正にいま求めていた質問があります。しかし、そこには私と同じ発想で「bedgraphに変換するのが適切じゃないかな（※スクリプトは自分で書く）」「もしくはbigWigもいいよね（※スクリプトは自分で書く）」なんて助言が……。
考えることはみんな同じなんですね！
ELMERパッケージ
そんな訓練されきったバイオインフォマティシャンのコメントをかいくぐり、メチル化アレイ→IGVで可視化可能なファイル形式へ変換するツールをいくつか見つけました。
最終的に私が使ったのは、Bioconductorで公開されているELMERパッケージです。 createBigWigDNAmetArray というわかりやすい名前の関数があります。
メチル化アレイ2種類（450K または EPIC）に対応しており、私の場合はGenomeStudioで作成できるFinalReportを読みこんでこの関数に渡しました。大変簡単でした（あえて不満を言うなら、ドキュメントは不親切でした……）。
library(ELMER) library(data.table) # FinalReportファイルの読み込み dat <- fread("FinalReport.txt", header=TRUE, stringsAsFactors=FALSE, skip=8, data.table=FALSE) # 1列目のTarget IDを行名に設定し、Target IDの列は削除する rownames(dat) <- dat[,1] dat <- dat[-1] # シグナルデータだけを取り出す。列数はサンプル数に合わせて指定 sig_dat <- dat[c(1,4,7,8)] library(ELMER) createBigWigDNAmetArray(sig_dat, genome="hg38", met.platform="EPIC", dir="IGV_tracks")
これでdir引数で指定した名前のディレクトリ（上記の例ではIGV_tracks）に、各サンプルのシグナルがbigWigに変換されて出力され、IGVで閲覧することができます。
かんたーん！
なお、data.tableパッケージの fread関数で読み込んだオブジェクトをデータフレームとして取得する（ data.table=FALSE ）のは私の習慣ですが、もしかしたらデータフレームにしなくても、createBigWigDNAmetArray関数は動くかも知れません（検証はしていません）。
