ゲノム解読論文は publication に値するか？
参考文献：
The sequence is dead: long live the genome.
Nature Biotechnology 29, 463 (2011) doi:10.1038/nbt.1901
Published online 07 June 2011
http://www.ncbi.nlm.nih.gov/pubmed/21654648
つい先日，ナツメヤシのゲノム解読論文が Nat. Biotech. に掲載された。
http://www.ncbi.nlm.nih.gov/pubmed/21623354
ナツメヤシのゲノムサイズは，~380Mb。
シークエンシング＆アセンブリ手法の概要は下記の通り：
we used the Genome Analyzer IIx to generate sequences 36–84 bp in length from genomic fragments of ~170 bp or ~370 bp. We assembled the genome by using the SOAPdenovo genome assembler...（中略）We used the SOAP Correction Tool to correct sequence reads before
assembly and closed gaps where possible with the SOAP GapCloser.
得られた contig / scaffold のサマリーは下記の通り：
Contig N50 size = 6,441bp
Scaffold N50 size = 9,339bp (NGS data only)
Scaffold N50 size = 30,480bp (NGS data with restriction map)
今回紹介する記事は，このナツメヤシの論文が，何故 Nat. Biotech. という高 IF ジャーナルへの publication に値するのか，に対しての Nat. Biotech. 誌の見解を述べている。
ボクがラフに読んだ限りでの要点を列挙する：
- ゲノム解読はもはやルーチンワークになりつつある。単一の高等生物ゲノムを解読しただけの論文は，今後，高 IF ジャーナルから姿を消してゆくことだろう。
- 一方で，ゲノム配列情報がその後の研究を促進する側面は否定できない。そのため，Nat. Biotech. 誌は，今後もゲノム解読論文を掲載する。
Nat. Biotech. 誌のゲノム解読論文に対する姿勢について，もう少し詳しく紹介する。
For all these reasons, Nature Biotechnology fully expects to continue publishing whole genome sequence papers. These papers will be of interest to a broad biotech audience because they reveal biologically or commercially relevant insights―insights that could have arisen only from a genome-wide analysis―or they report technical breakthroughs (e.g., the sequencing of complex genomes that are highly polyploid or outbred/heterogeneous in nature). Expect to see more papers that describe not one but many genomes where the value arises from comparisons of different sequences.
生物学的または商業的にオモシロい生物のゲノム解読論文は掲載しますよ。複数種よんで比較ゲノム解析しているゲノム解読論文はウェルカムですよ。ということらしい。
Nat. Biotech. 誌同様，Nat. Genet. 誌もまた，最近ゲノム解読論文を掲載している印象がある。
NGS業界の市場規模
NGS業界の市場規模って，どのくらいなのだろう？この問いに，興味はあっても，具体的な数字は知らなかった。そんなボクがふと，次の文章を読んだところで，驚いてしまった。イルミナの四半期売上が，$266.7 million。これは凄い。単純に４を掛けると，年単位の売上が１０億ドルに届くではないか。しかも，５４パーセントの伸び。これも凄い。NGS業界が，市場規模，成長率ともにこんなに魅力的だったとは知らなかった。ライフサイエンティストなら，NGS業界に投資するのもアリかもしれない。と真剣に考えてしまう。どの企業が伸びそうかなんて，予測がつきそうなもんだ。（株主優待で試薬が届いたりして。いや，個人ゲノム解読サービスの割引とかかな）第三世代，第四世代のシーケンサーを開発しているベンチャーの株を持ってみるのも楽しそうだし。（そもそも上場してないか）なーんてゆう，思考が膨らむ。いや，貧乏ポスドクに株を買う余裕があるのか？という現実もあるのですが。
The company booked $282.5 million in total revenues for the quarter, a 47 percent jump over the first quarter of 2010. Product sales generated $266.7 million, a 54 percent increase over the prior-year quarter that was driven by the sequencing business.
dbGaP への登録
SRA には controlled access なデータがあり，それらのデータは，dbGaP（Genotype and Phenotype; http://www.ncbi.nlm.nih.gov/sites/entrez?db=gap）にアクセス許可を申請することでダウンロードができるようだ。
dbGaP への申請は，このページからできる：
https://dbgap.ncbi.nlm.nih.gov/aa/wga.cgi?adddataset=phs000178&page=login
NIH scientists and staff でない人は，こちらのフォームから eRA commons アカウントの登録をする：
https://commons.era.nih.gov/commons/registration/registrationInstructions.jsp;
フォームの最初のページは，ある１点を除いては，埋めるのに困らない。ボクが困ってしまった１点とは，"DUNS Number" というところ。どうやら，これは，企業コードの世界標準らしく，次のページから検索ができるようだ。
このページに所属大学を入力すると，その後（社員が調べて）メールで送ってくれるらしい。って，今GWですよねー。この仕事は，しばらくお預けってことですかな。
fastq-dump.1.1.1 -> fastq-dump.1.2.0
久々に SRA からデータをダウンロードしようとしたら，ちょっとしたエラーが出てきた。正確には，ダウンロードのところではなく，sra-lite ファイルから fastq ファイルへ dump するところで。
エラーメッセージは，こんな感じ。
fastq-dump -A SRR101400 -M 0 /Path/To/SRA/SRRXXXXXX.lite.sra
2011-05-02 11:19:11 fastq-dump WARNING: path not found while selecting metadata within database module - physical column 'ALTREAD' failed to resolve
2011-05-02 11:19:11 fastq-dump WARNING: path not found while selecting metadata within database module - physical column 'QUALITY' failed to resolve
2011-05-02 11:19:11 fastq-dump WARNING: path not found while selecting metadata within database module - physical column 'READ' failed to resolve
2011-05-02 11:19:11 fastq-dump WARNING: path not found while selecting metadata within database module - physical column 'READ_LEN' failed to resolve
2011-05-02 11:19:11 fastq-dump WARNING: path not found while selecting metadata within database module - physical column 'READ_START' failed to resolve
2011-05-02 11:19:11 fastq-dump ERROR: type not found while resolving type within virtual database module - requested column name '.READ_TYPE' was found in metadata but failed to resolve type 'INSDC:SRA:xread_type'.
2011-05-02 11:19:11 fastq-dump ERROR: type not found while resolving type within virtual database module - failed to resolve column 'read_type', type { 49, 1 }
2011-05-02 11:19:11 fastq-dump ERROR: type not found while resolving type within virtual database module - failed to resolve column 'cs_key', type { 37, 1 }
2011-05-02 11:19:11 fastq-dump ERROR: type not found while resolving type within virtual database module - failed to resolve column 'CS_KEY', type { 37, 1 }
2011-05-02 11:19:11 fastq-dump ERROR: type not found while resolving type within virtual database module - failed to open '/Path/To/SRA/SRRXXXXXX.lite.sra'
このエラーの原因は，fastq-dump のバージョンが古くなっていたことだった。このエラーが出たときは，（たぶん）1.1.1 を使っていたのだけれども，1.2.0 にバージョンアップしたら，エラーが消えてくれた。
RNA-Seq Blog | RNA-Seq Analysis Software – Commercial Packages
もはや自分のためのメモでしかないのだけれども，RNA-Seq Blog にリストアップされていたコマーシャルツールを転載する。
- Avidis NGS – http://www.avadis-ngs.com/applications/rna_seq
- CLC Genomics Workbench – http://www.clcbio.com/index.php?id=1240
- DNASTAR QSeq – http://www.dnastar.com/t-products-qseq.aspx
- GeneSifter – http://www.geospiza.com/Products/default.shtml
- Geomatix Genome Analyzer – http://www.genomatix.de/en/produkte/genomatix-genome-analyzer.html
- Ingenuity IPA 9.0 – http://www.ingenuity.com/products/pathways_analysis.html
- Partek Genomics Suite (Partek GS) – http://www.partek.com/ngs#rnaseq
- RNA Baser Sequence Assembler – http://www.rnabaser.com/
- SeqSolve – http://www.integromics.com/ngs
- Softgenetics NextGENe – http://www.softgenetics.com/NextGENe_11.html
はは，上２つしか知らんかったなぁ。
参考文献：http://rna-seqblog.com/data-analysis/rna-seq-analysis-software-%E2%80%93-commercial-packages/
Annovar を使ってみる（マウス編）
え！？今まで使ったことないの？と言われてしまいそうであるが，そう，使ったことがなかったのだ。だって，多型解析やるの初めてなのですもの。実は。
ということで，Annovar を初めて使ってみた。ウェブページを読みつつ，作業工程を書いてみることにする。
0.1: Annovar のダウンロード
$ wget http://www.openbioinformatics.org/annovar/download/annovar.latest.tar.gz
...
$ tar zxvf annovar.latest.tar.gz
$ ls
annotate_variation.pl* example/ summarize_annovar.pl*
auto_annovar.pl* humandb/
convert2annovar.pl* retrieve_seq_from_fasta.pl*
ダウンロードには，2-3時間掛かった。
ミラーサイトからダウンロードした方が良かったのかしら。
0.2: "Quick start guide" をやってみる
http://www.openbioinformatics.org/annovar/annovar_startup.html
このページに書いてあることをやってみた。待ち時間は長かったけど，問題なくできた。
1: ここからが本番。まずは，マウスデータをダウンロードするところから。
$ annotate_variation.pl -downdb -buildver mm9 gene annovardb
...
NOTICE: the FASTA file http://www.openbioinformatics.org/annovar/download/mm9_refGeneMrna.fa.gz is not available to download but can be generated by the ANNOVAR software. PLEASE RUN THE FOLLOWING TWO COMMANDS CONSECUTIVELY TO GENERATE THE FASTA FILES:annotate_variation.pl --buildver mm9 --downdb seq annovardb/mm9_seq
retrieve_seq_from_fasta.pl annovardb/mm9_refGene.txt -seqdir annovardb/mm9_seq -format refGene -outfile annovardb/mm9_refGeneMrna.fa
...
どうやら，mm9_refGeneMrna.fa.gz をダウンロードできなかったらしい。そこで，新たな２つのステップが求められている。
とりあえず，素直にしたがってみる。
$ annotate_variation.pl --buildver mm9 -downdb seq annovardb/mm9_seq
...
$ retrieve_seq_from_fasta.pl annovardb/mm9_refGene.txt -seqdir annovardb/mm9_seq -format refGene -outfile annovardb/mm9_refGeneMrna.fa
...
WARNING: 56 regions occur more than once with discordant sequence length (for example, NR_037254, NM_001168334, NM_001123367, NM_011022, NR_037261)
ゲノム配列をダウンロードして，mm9_refGene.txt に書いてある遺伝子位置情報から，mm9_refGeneMrna.fa を作成している模様。
ここで，Ensembl のアノテーションを使いたいと思っていたことを思い出す。-downdb の引数を変更して，再度ダウンロード。
$ annotate_variation.pl -downdb -buildver mm9 ensGene annovardb/
...
（先ほどと同じエラーメッセージ）
$ retrieve_seq_from_fasta.pl annovardb/mm9_ensGene.txt -seqdir annovardb/mm9_seq -format ensGene -outfile annovardb/mm9_ensGeneMrna.fa
うむ，できた。region 関連の情報は，今回は使わなそうなので，ダウンロードしなくてよさげ。
2: 入力ファイルを用意する
入力フォーマットは，こちらに分かり易く書いてある。
ANNOVAR takes text-based input files, where each line corresponds to one variant. On each line, the first five space- or tab- delimited columns represent chromosome, start position, end position, the reference nucleotides and the observed nucleotides. Additional columns can be supplied and will be printed out in identical form.
Pileup format，SOAPsnp format など種々のフォーマットからの変換スクリプトが用意されているみたいだが，今は独自フォーマットで変異データを持っているので，自前でスクリプトを書いて，入力ファイルを用意することにした。
3: いざ，アノテーション
$ annotate_variation.pl -buildver mm9 annovar/hoge.annovar annovardb
...
$ cat annovar/hoge.annovar.variant_function
exonic Prrc2b 2 32069928 32069928 A C
...
$ cat annovar/hoge.annovar.exonic_variant_function
line1 synonymous SNV Prrc2b:NM_001159634:exon16:c.A3897C:p.A1299A, 2 32069928 32069928 A C
おぉ，できた！思ったよりもずっと簡単にできて，すてきでした。
PacBio社が第三世代シーケンサー用二次解析ツールSMRTをオープンソースで配布
As Pacific Biosciences prepares to ship its PacBio RS single-molecule real-time sequencing system later this quarter, the company has released SMRT Analysis, an open-source, secondary-analysis software suite designed to handle the system's long read data.
Pacific Biosciences 社が，SMRT Analysis というオープンソースツールを配布し始めた模様。今四半期から販売を開始する PacBio RS system のデータ解析をサポートするための準備であり，主に de novo assembly などの二次解析とその図示化に主眼をおいて開発された。このツールをオープンソースにした理由は，第三世代シーケンサー用のソフトウェア開発を促進させる狙いがある。それもそのハズで，今の第二世代シーケンサーのシェアを奪ってゆくためには，第三世代シーケンサー用のサードパーティによるソフトウェア開発は重要な要因である，との判断なのだろう。
SMRT Analysis は，APIが用意されており，新たに開発したソフトウェアを，DevNetと呼ばれる開発者ネットワークと統合・シェアすることができる。この点は非常に画期的だと，ボクは思う。なぜならば，シーケンサーメーカが積極的に，サードパーティ開発チーム間のコラボレーションをサポートしようとしているからだ。現に，DevNetからは，APIやソースコードだけでなく，シーケンスデータやシーケンサーに関してのドキュメントが揃っているようで，サードパーティによる開発を促すことに力を入れていることが分かる。すでに，DevNetには800〜1000人のユーザが登録をしているようだ。（すごい人数だ。きっと，今既にツール開発の競争が繰り広げられており，そう遠くない未来に，論文として表に出てくるのであろう。）
These features include long read lengths, expected to range between 850 and 1,500 bases; high granularity, which makes it possible to run multiple samples at a time; two new sequencing modes in addition to standard sequencing; circular consensus and strobe sequencing; as well as kinetic information, which the firm says can provide data about modifications in DNA and RNA gathered during the sequencing process.
PacBio RS system の特徴は，第一に 850-1500 bp とリードが長いこと。circular consensus によるエラー検出，ストロボシーケンスによるスキャッフォールディングは，それぞれ新しいソフトウェアを必要とするだろう。また，kinetic 情報を用いることで，メチル化などの DNA/RNA 修飾を検出するための手法も開発が望まれる。
SMRT Analysis は，SMRT Pipe と呼ばれる Python-based framework と，SMRT View というゲノムブラウザから構成されている。SMRT View は単にブラウザとしての機能だけでなく，ジョブ投入のためのGUIとしても機能する。SMRT Pipe からは，BLASR（mapping tool），Allora（de novo assembly tool），AHA（hybrid de novo assembly tool），EviCons（マッピング結果からコンセンサス塩基を推定する）などが利用できる。hybrid de novo assembly は，PacBio RS からのデータと，第二世代シーケンサーからのデータを組み合わせることにより，より長いスキャッフォールドを作ることを目的としている。
BLASRについての情報：
Sorenson said the algorithm takes currently used approaches, such as suffix arrays and dynamic programming, and "puts them together in a way that hasn't been done before."
The initial set of candidate alignments is found by querying a pre-computed index of the reference genome, and then refined until only high scoring alignments are retained. The base assignment in alignments is optimized and scored using all available quality information, such as insertion and deletion quality values. Because alignment approximates an exhaustive search, alignment significance may be computed by comparing optimal alignment score to the distribution of all other significant alignment scores.
Sorenson said BLASR performed better than other long-read alignment algorithms, such as BWA-SW, which, under specified parameters, can align PacBio data. He also said that BLASR outperforms some well-known sequence alignment algorithms including Blast, MUMmer, Exonerate, and Blat, and noted that the firm plans to publish a paper providing specific benchmark details for BLASR as compared to other methods.
Alloraについての情報：
Another tool, Allora, short for "a long read assembler," is PacBio's de novo assembly algorithm. Based on the open source assembly software package AMOS as well as other components tailored to PacBio’s long reads and error profile, Allora uses an overlap-layout-consensus approach to iteratively assemble raw reads into contigs and then outputs them as Fasta sequence and cmp.h5 files.
de Bruijn graph でなく，overlap-layout-consensus アプローチをとっている。（これは，今後 PacBio のスループットが向上したときにも，動くのだろうか）
EviConsについての情報：
A final component of the suite, EviCons, produces consensus sequences from multiple sequence alignments generated from resequencing reads or contigs. The tool uses probabilities and a likelihood ratio test to separate alignments into regions of certainty and uncertainty and then uses base quality values and the Steiner framework to produce the best estimate of the local consensus sequence for uncertain regions.
