構造変化はゲノム変化の重要なソースであり、表現型の変化、遺伝性疾患、進化に関与する可能性がある。集団における構造的変異の程度は、主にNGSのおかげで、最近になって認められているに過ぎない。事実、いくつかのヒト個体のゲノムをシーケンシングすることにより、構造変化に関与するDNAを一塩基多型(SNP)よりも多く見つけることができる[論文より ref.8]。しかし、リードが短いため、これらの変異はSNPよりも同定するのがはるかに困難である。これまで提案されたほとんどの方法は、リファレンスゲノム上のリードをマッピングすることに頼っている。主なアプローチは、マップされたペアエンドのリードが、ペアエンドリードの予想されるインサートサイズと方向に関して不一致のマッピングを示す場合に、構造的なバリアントブレイクポイントをコールする[ref.7]。複雑なゲノムでは、主にリピートとマッピングエラーが原因で、高い誤陽性率が得られ、方法間のオーバーラップも小さ苦なる[ref.1]。注目すべきことに、コピー数変化にほとんどの注意と努力が当てているようであるが、inversionなどのbalanced構造変化はあまり研究されておらず、これはinversionのような構造変化がショートリードデータでは検出することがより困難であることを示唆している。
これらのアプローチはすべて、リファレンスゲノムおよび第1ステップであるマッピングに依存する。これは、利用可能なリファレンスゲノムがないか、品質が悪い、またはあまりにも遠すぎる生物の場合、強い制限となる。一方、リシーケンシングされたゲノムの全面的なアセンブリを行い、得られたアセンブリを補完することもできる[ref.6]が、デノボアセンブリは困難かつ資源集約的な作業であり、inversionに注目することで軽減できる。したがって、著者らが取り組む問題は、リファレンスゲノムやリードの完全なアセンブリを必要とせずにraw NGSリードでinversionのシグネチャを直接識別できる方法である。アセンブリのグラフで特定のパターンを検索することによって、未加工のrawリードから関心のある生物学的事象を直接コールするいくつかの方法が開発されている。それらのうちのいくつかは、SNPまたはsmall indels [ref.10,9,13]またはRNA-seqデータ[ref.11]におけるオルタナティブスプライシングイベント検出に専用されている。 Cortex var [ref.4]は、バブルパターンを生成するいかなるバリアントも検出できると主張しているが、inversionやその他のbalanced構造変化バリアントを特に対象としている訳ではない。
この論文の主な貢献は、de Bruijn graphにて、inversionによって発生したトポロジカルパターンの解析と正式なモデリングである。さらに、このようなinversionパターンを検出するアルゴリズムを提案する。このアルゴリズムはTakeABreakというツールで実装され、http://colibread.inria.fr/TakeABreak/からダウンロードできる。
(以下略)
TakeABreakに関するツイート
インストール
mac os10.12で動作テストを行った。
ビルド依存
- CMake 2.6+; see http://www.cmake.org/cmake/resources/software.html
- c++ compiler; compilation was tested with gcc and g++ version>=4.5 (Linux) and clang version>=4.1 (Mac OSX).
本体 Github
#Anaconda環境で導入、バージョン指定しないとversion2が入る。
conda install -y -c bioconda takeabreak
#またはソースからビルドする
# get a local copy of source code
git clone --recursive https://github.com/GATB/TakeABreak.git
# compile the code an run a simple test on your computer cd TakeABreak
sh INSTALL
./build/bin/TakeABreak -h
$ ./build/bin/TakeABreak -h
ERROR: Unknown parameter '-h'
[TakeABreak options]
[Input / output options]
-in (1 arg) : input read file(s) [default '']
-graph (1 arg) : input graph file (likely a hdf5 file) [default '']
-out (1 arg) : prefix for output files [default '']
[Graph options]
-kmer-size (1 arg) : size of a kmer [default '31']
-abundance-min (1 arg) : minimal abundance threshold for solid kmers [default 'auto']
-abundance-max (1 arg) : maximal abundance threshold for solid kmers [default '2147483647']
-solidity-kind (1 arg) : way to consider a solid kmer with several datasets (sum, one, or all) [default 'one']
-max-disk (1 arg) : max disk (in MBytes) [default '0']
-max-memory (1 arg) : max memory (in MBytes) [default '2000']
[Inversion options]
-max-sim (1 arg) : max similarity percentage between a and b' and between u and v' [default '80']
-repeat (1 arg) : maximal repeat size at the breakpoint (longest common suffix between a and b') [default '8']
-lct (1 arg) : local complexity threshold (LCT) [default '100']
[General options]
-help (0 arg) : help
-version (0 arg) : version
-nb-cores (1 arg) : number of cores [default '0']
-verbose (1 arg) : verbosity level [default '1']
——
k-merカウンタはDSKを使っている。
テストラン
git clone https://github.com/GATB/TakeABreak.git
cd TakeABreak/tests/data/
#実行
TakeABreak -in toy_example_reads.fasta,toy_example_with_inv_reads.fasta
- -kmer-size size of a kmer [default '31']
- -abundance-min minimal abundance threshold for solid kmers [default 'auto']
graphファイル(binary)と、breakpoint上のfasta配列が出力される(*1)。リファレンスがあるなら、この配列をblast検索すればブレイクピントは簡単に探せる。
原理上ホモの inversionは検出できないので、inversionを全て検出したいならコントロールと一緒に解析する必要があります。複数fastqを与えるには、fastqをリストにして、リストファイルを指定してランするようです。詳細はGithubのREADMEで確認してください。
引用
Mapping-Free and Assembly-Free Discovery of Inversion Breakpoints from Raw NGS Reads
Claire Lemaitre, Liviu Ciortuz, Pierre Peterlongo
AlCoB 2014 (Tarragona), Lecture Notes in Computer Science vol. 8542, pp. 119--130.
*1
k-merサイズをあげれば、出力される配列も長くなるが、長すぎると共通k-merの頻度が下がり、何も検出されなくなる。また計算時間も増加する。