CLC Genomics WorkbenchでPacBioデータ

11月はあっという間に過ぎていきました。
サンフランシスコで人類遺伝学会があり、本当はその報告もしたいのですが、日常の仕事に追われ・・・。
PacBioのセミナーは、既に論文やプレゼンで知っている内容だったので、私にとっては驚きは少なかったのも影響してるのかな？　

アカデミックな報告ではありませんが、Linked-Inでつながっていてまだあったことの無かったひとと、会場で会い、ランチをご馳走になりました。
前ブログ「ショートリードの憂鬱」での、PacBioの最初の紹介が、Expression Analysisというボストンの受託会社のウェブセミナーを見たあとだったのですが、当時からその会社の担当者とはLinkedInでやり取りしていて、今回やっと会えました。　CEOにも会ってきました。
そこでの会話の内容はお伝えできませんが、同じPacBioユーザ同士、共感できるものがありました。

また展示会場のCLC-Bio社のブースで、前職にいたとき知り合った、CLCデンマーク本社のひとと1年ぶりくらいに会いました。　懐かしかった。
というわけで今日はCLCの話。

Genomics Workbenchは、使い勝手の良いソフトですが、まだPacBio用の入口がありません。
しかし、諦めるのはまだ早い！　
まず、Pacのリード（filtered_subreads.fastq）は、Import＞Illumina ( !! )からインポートします。

Quality value はNCBI/Sanger を選びます。　Paired-endとかのチェックは全て外します。

PacBioのデータはHiSeqとかに比べると桁違いに小さいので、Windows版のDesktop Editionでも十分使えるでしょう。



インポートが無事終わると、このようになります。　シーケンスの塩基数に注目！　10000bp 越えもある！

 

今のはfilterd_subreads、この意味は、SMRT PortalというPac付属のソフトで、

1. Read長が5000 base 以上
2. Sub-read長が5000 base以上
3. Read Qualityが75以上

などという基準でフィルタリングした後のsubreadということです。

まだエラーコレクションも何もしていません。

 

Pacのエラーコレクションで、pacBioToCAというものがあります。　これを流すと、CCSなどのHighクオリティショートリードでCLRの精度が高まります。

こいつのOutputが、fasta + qual ファイルです。

fasta + qual形式ファイルの「エラー修正subread」を取り込むには、Import > Roche 454 ( !! ) 

fasta と qual ファイルを選んだらNextボタンでインポート

というふうにするとちゃんと入ります。　

まあ、fasta と qual があれば、fastqファイルが作れますけどね。
既に誰かが作ったfasta+qual --> fastq converterはいろいろあると思いますが、私はSEQanswersのこのスレッドに出ていた、このスクリプトはいいのでは？と思います。Perlで動きます。
メモリに全部展開するのが嫌？なひともいるでしょうけど。

さて、リードを入れて何が面白いかというと、CLCのGenomics Workbenchはレポート機能が充実しているからそれをやってみると、Pacのデータの特徴が簡単に見えるのです。GCの割合や、リード長の分布、なども1つのレポートにまとまっているから便利。
下図の左はエラーコレクション前、右はエラーコレクション後の、クオリティ（Phredスコア）の分布図です。　例えばこんなグラフを見て、確認しながら、LinuxのPacサーバに戻ってアセンブる。

 

Pacサーバでアセンブるわけは、Genomics Workbenchのアセンブラーが、de Bruijn Graph 式アセンブラーだからです。
ハイスループットのショートリードには向いていますが、長さがバラバラの超ロングリードには向いていないかも（あくまで「かも」ですが）しれません。
現バージョン5.5.1では、k-mer を64以上には増やせないようですので、これ以上のk-mer では試して無いのですが、どうでしょう？
k-mer=64 ではアセンブルの成績は良くありませんでした。


Genomics Workbenchでリードが取り込めれば、リファレンスがあればマッピングも、できます。

上側がエラーコレクション後、下がコレクション前です。E.coliのデータです。
マッピングパラメータはデフォルトでやったのですが、この例では、エラーコレクション前はSubread全体の37%しかマッピングできませんでした。それがエラーコレクションをした後は、99%以上がマッピングされました！　


SMRT Portal/Analysis というPacBio付属のアセンブラーやマッピングツールも勿論よくできています。　Pacデータをガッツリ解析するにはこっちの方が向いているでしょう。
GUI以外にコマンドラインも充実しています。
また、最初のsubreadフィルタリングなどには必須です。

でもCLCのGenomics Workbenchは、レポート機能が充実している、というか、手軽にデータを確認できるという点で優れていると思います。