ex3 jb - kuratanicde/20200617-nibb-omics GitHub Wiki

Build genome browser with JBrowse

共生バクテリアである Buchnera aphidicola のゲノムブラウザをJBrowseで作成してみよう。

References

開発環境の準備

UNIX環境。

テスト用の簡易ウェブサーバ。ここでは、以下の短いRuby scriptを使う。

(交換会後に補足:rubyの実行に問題がある場合は、以下のpythonのワンライナーでもウェブサーバを起動できます。)

python -m SimpleHTTPServer

今回、以下のAmazon AWSの仮想サーバでライブデモをやります。(情報交換会後はポートを閉じます。)

Set up JBrowse software

JBrowse開発者のHPから最新版のJBrowseソフトウェアをダウンロードする。

$ wget https://github.com/GMOD/jbrowse/releases/download/1.16.9-release/JBrowse-1.16.9.zip
$ unzip JBrowse-1.16.9.zip
$ mv JBrowse-1.16.9 jb1
$ cd jb1
$ ./setup.sh

web server を立ち上げ、インストールが正常に行われたかどうか確認。volvoxのテストデータが表示されればOK。 (webserver.rb はここから取得してください。コードをコピーペしてwebserver.rbの名前で保存)。

$ ruby webserver.rb

お好みのウェブブラウザで http://[IP_address]:8000 へアクセス。(手元のコンピュータの場合、http://localhost:8000)

Build reference genome

NCBI Genomes からBuchnera aphidicola のリファレンスゲノムデータを取得する。そのデータを使ってJBrowseでリファレンスゲノムを構築する。

上記 (https://www.ncbi.nlm.nih.gov/genome/170) にアクセスし、最上部の四角で囲われたリファレンス情報の中から、genomeとGFFを取得する。

  • GCF_000009605.1_ASM960v1_genomic.gff (genome annotation in GFF3 format)
  • GCF_000009605.1_ASM960v1_genomic.fna (genome seq in FASTA format)

Build reference genome sequence

jb1 ディレクトリ直下で、

$ cd jb1
$ bin/prepare-refseqs.pl --fasta ../GCF_000009605.1_ASM960v1_genomic.fna

ウェブサーバを起動し、ウェブブラウザでアクセスして確認。

次に、アノテーションgffデータをインポート。

$ bin/flatfile-to-json.pl --gff ../GCF_000009605.1_ASM960v1_genomic.gff --trackLabel NCBI_Annot --trackType CanvasFeatures

ウェブブラウザでアクセスして確認。

次に、name, id のindexingを行う。これによって、JBrowseから名前で検索が可能になる。

$ bin/generate-names.pl -v

ウェブブラウザでアクセスして確認。nameやIDで検索してみる。

Customize annotation track

trackの設定は、data/trackList.json にJSONフォーマットで管理されている。このファイルを編集することで柔軟に設定を変えたり、機能を追加することができる。

遺伝子IDからNCBIのページに飛べるようにする。

trackList.jsonを以下のように書き換える。fmtDetailValue_Name... の挿入

[trackList.json] (v0 => v1) | Full code : trackList.v1.json

<          "urlTemplate" : "tracks/NCBI_Annot/{refseq}/trackData.json"
---
>          "urlTemplate" : "tracks/NCBI_Annot/{refseq}/trackData.json",
> 	 "fmtDetailValue_Name": "function(name) { return '<a href=\"http://www.ncbi.nlm.nih.gov/protein/'+name+'\">'+name+'</a>'; }"

右クリックのサブメニューから外部サイトへ飛べるようにする。

右クリックでKEGGのデータベースにジャンプできるようにする。

[trackList.json] (v1 => v2) | Full code : trackList.v2.json

29c29,37
< 	 "fmtDetailValue_Name": "function(name) { return '<a href=\"http://www.ncbi.nlm.nih.gov/protein/'+name+'\">'+name+'</a>'; }"
---
> 	 "fmtDetailValue_Name": "function(name) { return '<a href=\"http://www.ncbi.nlm.nih.gov/protein/'+name+'\">'+name+'</a>'; }",
> 	 "menuTemplate": [
>               {
> 		  "label": "KEGG: {name}",
> 		  "iconClass": "dijitIconDatabase",
> 		  "action": "newWindow",
> 		  "url": "https://www.genome.jp/dbget-bin/www_bget?buc:{name}"
>               }
> 	 ]

遺伝子名でGooge検索の機能を実装する。

[trackList.json] (v2 > v3)| Full code : trackList.v3.json

36c36,42
<               }
---
>               },
>              {
>                  "label": "Google: {name}",
>                  "iconClass": "dijitIconSearch",
>                  "action": "newWindow",
>                  "url": "https://www.google.com/search?q={name}+Buchnera"
>              }

Add tracks

BLAST検索の結果を可視化する

E. coli の全遺伝子セットをBuchnera genomeに対してBLAST検索した結果をBEDフォーマットに変換し、それをjbrowseにインポートする。

E. coli K12株のリファレンスプロテオームデータを取得する。

$ wget ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/reference_proteomes/Bacteria/UP000000625_83333.fasta.gz
$ unzip 

BLAST search

$ makeblastdb -in GCF_000009605.1_ASM960v1_genomic.fna -dbtype nucl -parse_seqids
$ tblastn -db GCF_000009605.1_ASM960v1_genomic.fna -query UP000000625_83333.fasta -outfmt 6 -evalue 1.0e-5 -max_target_seqs 1 > UP000000625.vs.GCF_000009605.tblastn.fmt6

注)-max_target_seqs 1 は推奨されないオプションだが、今回はデモのためそのまま進める。

以下のフォーマットで出力される。 [UP000000625.vs.GCF_000009605.tblastn.fmt6]

sp|P00350|6PGD_ECOLI	NC_002528.1	68.951	467	143	2	1	466	111628	113025	0.0	684
sp|P00448|SODM_ECOLI	NC_002528.1	63.366	202	73	1	1	202	204463	205065	4.29e-92	289
sp|P00452|RIR1_ECOLI	NC_002528.1	75.362	759	187	0	1	759	195489	193213	0.0	1234
s
...

カラム左から、query acc.ver, subject acc.ver, % identity, alignment length, mismatches, gap opens, q. start, q. end, s. start, s. end, evalue, bit score

これをBED format に変換する。

blast6_to_bed6.rb

inputf = ARGV[0]

File.open(inputf).each do |l|
  a = l.chomp.split(/\t/)
  query = a[0]
  target = a[1]
  from = a[8].to_i
  to = a[9].to_i
  score = a[11].to_f

  left, right = [from, to].sort
  left = left - 1 # 0-based coordinate
  strand = ( from < to ? "+" : "-")

  [query, target, left, right, strand]

  puts [target, left, right, query, score, strand].join("\t")
end
ruby blast6_to_bed6.rb UP000000625.vs.GCF_000009605.tblastn.fmt6 >UP000000625.vs.GCF_000009605.tblastn.fmt6.bed

UP000000625.vs.GCF_000009605.tblastn.fmt6.bed ができた。これをjbrowse にインポート。

bin/flatfile-to-json.pl --bed ../UP000000625.vs.GCF_000009605.tblastn.fmt6.bed  --trackLabel Ecoli --trackType CanvasFeatures

webブラウザでアクセスして確認。

bam, vcfなど

bamやvcf のデータもインポートすることができます。 詳細は、開発者documentを見てください。

もちろんIGVの演習で使った、bam, vcfもインポート可能です。ただ、リファレンスの配列のIDが異なるので変換する必要があります。