IGVでの格好いい見せ方の、もとデータ

前回、IGVでPacBioを格好良く見せる方法を紹介しました。
その時例に使ったデータは何かと言うと、Sequelでヒトゲノム NA12878 を10xくらいの深度で読んだものです。
具体的には、ここのデータ

ライブラリの長さは25kb、Blue Pippinを使って15kbにサイズセレクション

• 使用したSMRT Cell 1M の数：10
• トータルラン時間：60時間
• 出力塩基数：32.8 Gb
• リード数：340万本
• リード長のN50 ：11.823 bp

このとき使用した試薬は、旧バージョン、v.1.2 のもの
なので今ならセルあたりの出力はもっと多いはず。

とにかくこれでヒトゲノムの10倍のデータが出た。

このデータをヒトゲノムにマップするのですが、ここで使ったツールは、NGM-LR　＋　PBHoney
PBHoneyは構造変異を検出するツールです。

NGM-LR って何？　という方、これはロングリード用のマッピングツールです。
Next-Gen Mapping tool for Long Read、だったかな？何かそんな名前。
Githubにもあるので、興味のある方はここからどうぞ。

PacBioリードは1本が長いので、例えば 1kb 程度の挿入・欠損をまたいで読むことが可能。
しかし通常のマッパーでは 1kb の変異を考慮してほかの配列を綺麗にゲノムにマップすることができなかった。NGM-LRは、二箇所に分かれてマップするような、ロングリード独特な性質をフルに発揮できるマッピングツール。

BWAとNGM-LRのマッピング結果　Aaronのスライドより

さて、こうしてマッピングした結果は、もちろん参照できます。

先ず、DNAnexusのデモアカウントと作りましょう！
いえいえ、決して私はDNAnexusの手先ではありません。
仕方無いんです。ここにアクセスした方が、データ参照が楽だから。

私はアカウント持っているのですが、ロングインするとこんな感じです。

左下の、"PacBio Sequel Data" というところをクリックします。
これが例の10カバレッジのSequelデータ

"Sequel Data" を開きます

NA12878.reads.ngm.bamというファイルが、マッピングファイルですが、20Gbもあって大きいです。
そこで、indexsession ファイルをダウンロードします。
これは、IGVに取り込むと、DNAnexusのサーバにアクセスしてデータを表示してくれるインデックスファイルです。
IGVはこのように、必ずしもローカルに大きなサイズのマッピングファイルを持っておく必要がありません。

さ、IGV を開きましょう。前回のあれ、ですよ。わからないひとはちょうどこの前の記事をチェック！

ゲノム配列が "Human hg19" であることをチェックして（違ったら Human hg19 を選ぶ）

File > Open Sessions 

今ダウンロードした　indexsession ファイルを選択。

何か聞かれるけどOKをする。

ゲノムのポジションを入力する場所に、試しに、

chrX:116453100-116453795

と入れてGo！

いぇーーーーい !!

ほかにもチラッと見てみたいポジションは、こんなところかな？

chrX:116454160-116454859
chr10:92213800-92216245

InDelの変異箇所は、先の　DNAnexusのアカウントから、

NA12978_Output を選び、ここから ... del.bed  ins.bed ファイルをダウンロードしてきて、取り込むと面白いかも！

いかがだったでしょうか？

IGVでの格好いい見せ方

Integrative Genomics Viewer (IGV) といえば、NGSやっているひとなら一度は聞いたことのある、ゲノムビューワーですね。
フリーで使えてさくっと見るのにはとても便利。
先日、本社とのウェブミーティングで、このIGVを使って、PacBioのデータを格好良く見せる方法を教わりました。
皆さんも知っていると便利なツールだったのでシェアします！

まず、普通にIGV（ここではv.2.3.88）を使ってPacBioのアライメントファイルを開くとこう見えます。

ぜんぜん綺麗じゃない！　
SubreadをアラインしているのでInDelエラーが目立ってるんです。

これをある方法を使うと、このように見ることができます。

とてもすっきりしてますねえ。

これ、Development Snapshot版を使っています。

先ず、バイナリを落としてきます（上図カーソルの場所から）。
Zipを解凍して、中に作られる、igv.batファイルをダブルクリックして起動させます（テスト環境はWindows 7&10）

最初は、これまでIGVで表示していた、デフォルトで、PacBioアライメントが表示されるかもしれません。ここで、View > Preferences　を選択

Alignmentsタブを開き、とりあえず以下のように数字を入れてみて下さい。

特に、

• Label indels > 1：　2塩基以上のindelに表示が出る
• Hide indels < 20：　19塩基以下のindelは表示されない

としてチェックを入れて、OKします。

ちなみに、普通のIGV(v.2.3.88)でも、これに似た画面はありますが、上記オプションはありませんね。

これはv.2.3.88のIGV

さてこれで、先ほどのような綺麗な表示ができるようになりました！

こんな感じに。

しかし例えば、下図のような場所があったとします。

大きなDeletionがある領域です。

そのほかにも、小さなInDelが点在するようですね。

では、Allele Frequencyが小さいものはエラーの可能性が高いので、除くことにします。

例えば、35％以上のアレルだけを残したい、そんな場合は、再び View > Preferences > Alignments から、Coverage allele-fraction threshold: を、0.35としてみる

そうするとこんなにすっきりします！

大きなDeletionがあるところの上流側に、ハプロタイプ（ここでは黄緑のA、リファレンスはT）があるのに気がつきます。

この黄緑色Aを右クリックして、Group alignments by > base at [塩基の場所]　を選択

するとこうなります

大きなDeletionと上流SNV（T -> A）がリンクしているのがわかると思います。

IGVはあくまでビューワーですので、これ自体が解析するソフトではありませんが、データを見せる方法としてはとても使えるツールだと思います。

尚、この機能は、PacBioのプログラマー、Aaronさんによって作られました。

さて、次回は、この例に使用したヒトゲノムデータについてお話ししたいと思います。

多分誰でも使えるデモデータ、のはず。。

10X Genomics のデータだけに頼ってはいけない理由

 10X Genomics （10XG）という会社をご存知ですか？
Synthetic Long Read、またの名をLinked Read と呼ばれるデータを出力し、Scaffoldを作ります。
Synthetic Long Readというのは訳すと、合成ロングリード？擬似的ロングリード？でしょうか。
でも、これは本当の意味でロングリードじゃない！　という意見があったのか、最近は Linked Read、リンクしたリード？、という呼び方がされるようです。

http://www.10xgenomics.com/products/より

バンクーバーのASHGでも、ランチョンやその他のワークショップでも存在感を出していた10XGですが、この会社の対象マーケットは、大きく、1)デノボアセンブリ、2)シングルセル発現解析、の2つがメインだと言えます。
シングルセル発現解析のほうは既存の技術に比べて何十倍ものハイスループット解析ができる、といった優位性があると思います。
ですが、もうひとつのデノボアセンブリは、どうかなあ・・・
10XGに頼って解析を進めていると、大きな情報を見落としていることに気づかないことがあるでしょう。今日はそういう話。

ここから先は、PacBioのJasonさん（Falconなどの開発者）が解析してくれた結果。
彼はバイアスをかけない解析をしますので、私は、とても客観的なデータだと思います。


まず、NA19240サンプルを、PacBioでアセンブルした結果と、10XG v.1.3 でアセンブルした結果があります。
Contiguityを比較するのは大人気ない、というかそもそもContigとScaffoldを比較するのはあまり意味が無い。
そこで、アセンブルした後の結果の、構造変異を検出した場所（Segmental Duplication以外の場所で）を見てみます。

• 10xGアセンブリ、PacBioアセンブリ共に、数千～1万強の、挿入、欠損部位を検出した
• 10xGアセンブリの方が多くの欠損変異を検出し、5,573個は、PacBioアセンブリでは検出していなかった
• 反対に、2,693個の欠損変異は、PacBioアセンブリでは検出され、10xGアセンブリでは検出されなかった

【欠損：10XGで検出されたがPacBioで検出されなかった5,573のうちの１例】

これは、10XGのアセンブリでは176bpの欠損を示していた箇所。
リファレンスには確かに176bpあり、PacBioのアセンブリでも、アセンブリ前の補正後リード（p-reads）でも、確かに176bpがあった例。
この箇所は、GCの連続配列で、10xGが使用しているショートリードでは読めない。
なので間違ってコールされてしまった例。

ほかランダムに抽出した９例

最初の１例を合わせた10例のうち、8例は擬陽性だった。

【欠損：PacBioアセンブリでは検出されたが、10xGでは検出されなかった2,693例のうちの1例】

この配列も、GCが多い場所で、10xGでは読めていない。しかし、この場所には、リファレンスには無い58塩基の欠損領域が、PacBioアセンブリから検出できている。
しかもそれは、ヘテロ欠損。
研究者なら、このようなヘテロな場所の変異に、より興味が湧くのでは深いのでは？

ほかにランダムに抽出した9例

先の1例を除き、全てで10xGのアセンブリでは、読めていない箇所に、実際は欠損変異が存在した。

ということで結論、
読めない配列が原因でコールされた変異は擬陽性の可能性が高い
読めない場所から変異をコールするのは不可能
10xGのデータだけに頼って解析していては、重要な変異を取りこぼすことになる

今回は欠損変異だけに注目しました。（スライドはJasonさん作）

どんな技術もそうですが、それだけに頼って解析していては、見落とすものは必ずあります。　なのでバリデーションは大事。
特に、比較的新しい技術に関しては、一見素晴らしいように見えても、欠点もちゃんと意識して使わないと、レビューワーから突っ込みを受けることになりますね。

10xGユーザの方は、上記のような擬陽性の可能性もあることを考えて、一度、PacBioで読んでバリデートしましょう！！

アメリカ人類遺伝学会 ASHG PacBioネタ

ちょっと報告が遅れましたが、ASHG＠バンクーバーのPacBioネタです。

初日のマクロジェン社のランチョンでは、10Kアジア人ゲノムプロジェクトの話がありました。当然PacBioのにとも何人か参加していましたよ。もちろん私も。
韓国人ゲノムAK1は、最近論文になりましたよね。
これからはアジア人を10000人、読みまくる！という威勢のいい話でした。
もちろん全部をPacBioで読むわけではないでしょう。

二日目はPacBioのランチョンセミナー
これはもう、YouTubeで紹介されていますので、興味のあるひとはここを是非チェック！

「ヒトゲノム＋構造解析」は、これからのPacBioの大きなテーマです。
今までは、ゲノムアセンブリに注力してきたPacBioですが、アセンブリの分野ではもうグローバルスタンダード、になりました。
これ以上長く読める、技術的にも安定したシークエンサーは、今のところPacBio以外存在しない！　というのは過言ですが、客観的に言ってもPacBioのSMRTテクノロジーは素晴らしい技術です。
巷ではナノポア技術が競合とみなされているようですけれど、まだまだ実績が違う。
PacBioで読んだ結果は、ちゃんとこれまで1700本の論文になりましたから。
PacBio は、サイエンスが議論できるレベルのシークエンサーです！


さて、学会期間中に新商品やサービスをアナウンスするのは企業によくあるパターンですが、PacBioもやってくれました。

その名も、新試薬、バージョン 1.2.1 for Sequel
て、地味だ。地味すぎる・・・

いまのSequel 試薬が1.2で、それの新しい版が1.2.1
やっぱり地味だ。

でもこの試薬のすごいところは、Sequelの最初の必要DNA量は変わらないけれど、そこからランできるSMRT Cellの数が増えたこと。
ゲノムDNA20ugくらいからスタートして、今までSequel は数セルしかランできなかったけれど、10セル程度かそれ以上ランできるように、ライブラリ作製効率が上がった（注：この辺の数字は大きくサンプル依存）。
具体的には、プレート上にロードするDNA濃度は、RSIIと同じ、6.5～40 pMになった。
プライマーも新しくなった。
その他、試薬ではないけれどロボット系の動きも少し変わった。

ベータテストの結果、20kb以上の長いライブラリ（w/ サイズセレクション）を読んで、出力は1セルあたり5～7Gb
平均リード長は10kb、11kbくらい（平均ではなくてN50 ならもっと長いですけどね）
バクテリア、植物、ヒト細胞から同じような出力結果が出ているようで。

さらに来年春のアップデートではもっと大きく改善する予定。
と、その時に名前は 1.2.2 とかではなくて、もっとインパクトがあるものにして欲しい！

PacBio v.s. Oxford Nanopore 特許侵害訴訟に！ うーん、仲良くやりたいんだけどなあ

11月3日に飛び込んできたニュース
この話は書こうかどうか迷ったけれど、これは一応個人のブログだし、PacBio社も2日経っても公式に何も言っていないので、好き勝手書きます。

PacBioがOxford Nanopore Technologiesを、特許侵害で訴えた
ニュース記事はこちら

これはですね、PacBioが2013年4月に出願していた特許「U.S. Patent No. 9,404,146」が、今年の8月に公開され、その中に記載されている請求項のいくつかが、ONT社が使っているテクノロジーに抵触している、ということ。
で、この特許･･･146は、すごく広範囲を抑えているんですね。私も全部読んだわけではないですが、キーは、「二本鎖の一分子DNAを、Sense・Anti-sense両方読む」というところが、どうやらONTの2Dシークエンス技術が侵害している、ということらしいです。

PacBioは、ダンベル型のライブラリを作って読むので、二本鎖DNAの両方のStrandを、片側ずつ何度も読むことができるCCSという方法があります。
ONTも、片側ダンベルアダプターをつけて、二本鎖DNAの両方のStrandを、片側１回ずつ読む、2Dシークエンスという方法がある。

また、この特許･･･146は、一分子DNAをそのまま読む、という方法論についてもカバーしている。
つまり一分子シークエンサーは、ことごとく抵触してしまう！

この特許は、国際出願されていないらしいので、とりあえずアメリカマーケットが問題になるのでしょうね。でもアメリカは世界最大のマーケット。日本への影響は多大でしょう。

この宣戦布告の発表があったその日、PacBioの第３四半期の発表もありました。
絶好調とまでは言わないまでも、これまでにSequelを75台も売り、アップデートの予定など今後の見通しも明るい、というものでした。
でーも、特許侵害のニュースのインパクトの方が大きかったのか、翌日の株価は20%の下落！
まあ、「トランプ大統領候補が勝つかも知れない不安」のせいで、Nasdaqも全体的に少し下げたんですが、20%は無い。


私は、個人的には、仲良くやって欲しいというのが本音です。
知っているひとは知っていると思いますが、前からの友人もいますし。
立場上、ライバルなので言うことは言いますけどね。
でも一緒に業界を盛り上げていこう、という気持ちです。

これからどうなるのか。