2017年07月01日(土)

細胞と抗体の名前付きでENCODEのeCLIPのデータをダウンロードする方法

投稿者:

札幌はガトーキングダムで開かれた内藤カンファレンス。自宅から車で40分ほどにもかかわらず3泊4日もすることになり家族が悲しむかと思いきやむしろ喜んでいるようだったのでちょっぴり悲しくなりつつも、今回のテーマはずばり、"Noncoding RNA: Biology, Chemistry, & Diseases"。GordonやKeystoneに引けを取らない豪華な海外ゲスト。超ウルトラスーパーデラックスエキサイティングでグレートでプロジーヂャスなミーティングは非常に熱気に満ちていて、興味のある話はいくら聞いていても飽きないということを本当に実感した4日間でした。

今回のミーティングの圧巻はHoward ChangさんとGene Yeoさんのトーク。NGSでできることをやりつくす、というか焼き尽くす、と言った趣で、彼らと同じことをやっていてももうだめだ、ということをつくづく思い知らされたような気もします。まあ、逆に言えばやらなくても良いことは何かということが気づけたということで、自分の立ち位置を見直す良い機会になった、ということもできるでしょうか。いずれにせよ、RNAの二次構造全部決めちゃいました、70種類のRNA結合蛋白質の結合サイト全部決めちゃいました、というのはリソースとしては非常に重要で、特に70種類ものRNA結合蛋白質の結合サイトに関する情報は、今後自分の実験をデザインするために非常に有用なツールになると思われます。データを独占せずにコミュニティー全体が使えるようにしてくれている。それによってコミュニティーが活性化されてまた有用なデータが蓄積される。その好循環が回っているのがこの業界のすごいところだと思います。GeneさんのCLIPはENCODEのプロジェクトなのでそのうちにUCSCブラウザやその他のプラットフォームでも見ることができるようになるでしょうが、しばらくかかるようなので、当面は自分でダウンロードして使うことになりそうです。ENCODEのページに行くとfastqやらbamやら全てのCLIPデータをダウンロードできるページはあるのですが、そのままダウンロードしてもファイル名に細胞や抗体の名前が含まれていないので、はあ?みたいな感じで非常にわかりにくいのと、全部ダウンロードするととんでもないサイズになりそうなので、IGVで表示させるのに便利なbigWigファイルとbedファイルを抗体と細胞の名前つきでダウンロードするシェルスクリプト作ったので置いときます。hg19とGRC38にマップしたファイルがあるようなので、それぞれ別のシェルスクリプトでダウンロードできるようにしてあります。

GRC38のbedファイルをダウンロードするシェルスクリプトはこちら
GRC38のbigWigファイルをダウンロードするシェルスクリプトはこちら
hg19のbedファイルをダウンロードするシェルスクリプトはこちら
hg19のbigWigファイルをダウンロードするシェルスクリプトはこちら

これらのshファイルをダウンロードして(注:ダウンロードすると.txtの拡張子がつくみたいなのでそれを削除してください)、Macのターミナルでファイル名を指定すれば、カレントディレクトリに次々とダウンロードされるはずです。例えばデスクトップにENCODEフォルダを作り、そこにスクリプトのファイルをおいてGRC38のbedをダウンロードするのなら

$ cd Desktop/ENCODE
$ ls
Encode_download.R		download_bed_hg19.sh		download_bigWig_hg19.sh
download_bed_GRCh38.sh		download_bed_hg19.txt		metadata.tsv
download_bed_GRCh38.txt		download_bigWig_GRCh38.sh
$ sh download_bed_GRCh38.sh

みたいな感じ。パッと見た感じではeCLIP実験の再現性は非常に高いようなので、ブラウザで見るだけならduplicateは一つだけダウンロードするので良いかもしれません。この辺り、ダウンロードURLを含めたENCODE CLIPの全てのmetadataが入ったファイルがありますので、Rでいろいろいじれば好きなセットをダウンロードするファイルができると思います。ちなみに、ダウンロード用のシェルスクリプトを作る時に使ったRのコマンドもつけときます。中身は下記の通りです。恥ずかしいぐらい単純なことしかしてませんし、エクセルのコピペでも同じようなことすぐできますが、あれ、これどうしたっけな、とかとんでもないミスも出るので、やはりRは便利ですね。

#read metadata
meta<-read.delim("metadata.tsv")
head(meta)

#remove mock data
mock<-grep("mock",meta$Experiment.target)
meta<-meta[-mock,]
head(meta)

#select required column only
meta1<-meta[,c(2,3,7,17,30,42,43)]
head(meta1)

#change names
meta1$File.format<-sub("bed narrowPeak","bed",meta1$File.format)
meta1$Output.type<-sub("minus strand signal of unique reads","minus",meta1$Output.type)
meta1$Output.type<-sub("plus strand signal of unique reads","plus",meta1$Output.type)
meta1$Experiment.target<-sub("-human","",meta1$Experiment.target)
head(meta1)

#select bigWig
bigWig<-subset(meta1,meta1$File.format=="bigWig")
head(bigWig)

#select bed
bed<-subset(meta1,meta1$File.format=="bed")
head(bed)

#create bed download command
bed$download<-paste("wget ",bed$File.download.URL," -O ",bed$Biosample.term.name,"_",bed$Experiment.target,"_",bed$Biological.replicate.s.,".bed.gz",sep="")
head(bed)

#create bigWig download command
bigWig$download<-paste("wget ",bigWig$File.download.URL," -O ",bigWig$Biosample.term.name,"_",bigWig$Experiment.target,"_",bigWig$Output.type,"_",bigWig$Biological.replicate.s.,".bigWig",sep="")
head(bigWig)

#select genome assembly<
bigWig_hg19<-subset(bigWig,bigWig$Assembly=="hg19")
bigWig_GRCh38<-subset(bigWig,bigWig$Assembly=="GRCh38")
bed_hg19<-subset(bed,bed$Assembly=="hg19")
bed_GRCh38<-subset(bed,bed$Assembly=="GRCh38")
head(bigWig_hg19)
head(bigWig_GRCh38)
head(bed_hg19)
head(bed_GRCh38)

#create download command file
write(bigWig_hg19$download,"download_bigWig_hg19.txt")
write(bigWig_GRCh38$download,"download_bigWig_GRCh38.txt")
write(bed_hg19$download,"download_bed_hg19.txt")
write(bed_GRCh38$download,"download_bed_GRCh38.txt")

これでできた.txtファイルの最初の行に#!/bin/shとつけて保存し、拡張子を.shに変えればターミナルのshコマンドで走らせられます。IGVで眺めているといろいろ実験のアイデアが湧いてきそうです。NGSの素晴らしさを手軽に体感できますので皆さんも是非!

(追記)

中川 真一

北海道大学 薬学研究院 教授
▶ プロフィールはこちら

ブログアーカイブ

ログイン

サイト内検索