超並列シーケンシング(MPS)[論文より ref.1]のためにHi-Plex(www.hiplex.org)と呼ばれる高度に多重化されたPCRベースのターゲット濃縮システムを開発した。 Hi-Plexは、非常に正確な結果を得ることができるシンプルで低コストのプロトコルである。その主要な特徴の1つは、サイズ選択によるオフターゲット増幅の除去を容易にする均一なライブラリーサイズを定義する能力であり、ペアエンドシーケンシングと組み合わせて、各アンプリコンについてのペアエンドの完全なオーバーラップを可能にする。後者は、バリアントの検出およびシーケンスエラーによって引き起こされるアーティファクトのフィルタリング両方において高度な厳密さを可能にする。以前、Hi-Plexシーケンス [ref.2](参考ページ)によって生成されたオーバーラップリードを使うariable CallingツールであるROVERを開発し、PALB2のコード領域の遺伝的な変異のスクリーニングにHi-PlexとROVERを適用した。結果は、以前のスクリーニングによって同定された全60変異を全て検出し、偽陽性はゼロだった[ref.3、4]。
ROVERは、ターゲットとしてアンプリコン領域のゲノム座標をタブ区切り形式で記述したファイルと、リファレンスゲノムにマップされたペアエンドリードを含むBAM形式[ref.5]の1つ以上の配列ファイルを必要とする。出力はVCF形式のバリエーションのリストである[ref.6]。 ROVERは、SNVおよびsmall indelを検出することができる。バリアントは、オーバーラップするペアエンドリードの両方で同じ位置に表示される場合にのみ、ROVERによってレポートされる。
ROVERを用いてバリアントを検出するのに最も計算量が高価になる部分は、リファレンスゲノムへのマッピング(またはアライメント)である。以下で説明する1つの指標実験では、Bowtie(http://bowtie-bio.sourceforge.net/)[ref.7]とのマッピングに要した時間は、ROVERバリアントコールパイプライン全体の約78%を構成していた。リードのマッピングは全エキソームおよび全ゲノムDNAシークエンシングパイプラインの標準的な部分であるが、これまでにAmplivar [ref.8]によって実証されているように、これはmassively parallel sequencing (MPS) ベースのPCRによって避けることができる。これは、各シーケンスリードの5 '末端がPCR増幅のプライマー配列で始まり、そして既知の配列を増幅するからである。後者の情報は、プライマー設計中に決定される。プライマーペアが介在する配列のゲノム座標を同定するのでリードはリファレンスにマップする必要はない。オプションで、ペアの少なくとも1つのリードのプライマーに続く小さな配列がリファレンスと同一であることを確認することで、正しい位置にマッピングされるという確信をさらに高めることもできる。
特定のリードの開始座標を決定したら、そのシーケンスをリファレンスと比較することができる。挿入と削除検出を可能にするためには、Needleman-Wunschアルゴリズム[ref.9]のスタイルでギャップアライメントを実行する必要がある。しかしながら、このアルゴリズムの複雑さは、アライメントされたシーケンスの長さの2次関数であり、したがって、入力における各リードについて計算するコストが大きくなる。幸運なことに、入力のほとんどのリードはリファレンスと同じになるか、少数の不一致でしか違いがないため、このコストを避けることができる。このような場合、リードとリファレンスの単純な線形比較で十分である。線形の比較が失敗した場合にのみ、ギャップアライメントアルゴリズムに戻る。
Amplivarは、UNDR ROVERの基礎となる構想に似ているが、異なる仕組みを適用している。 Amplivarはプライマー配列を使用してアンプリコンとリードを関連させ、リードをグループとしてアライメントさせることで計算上のオーバーヘッドを削減する(BLAT [ref.10]を使用)。さらに、Amplivarは(SeqPrep [ref.11]を使用して)オーバーラップリードをマージするが、UNDR ROVERは両方のリードを保持して厳密なフィルタリングシステムの一部として一致するかテストする。
UNDR ROVERは、Hi-Plexターゲットシーケンスをサポートするように設計されているが、シーケンシングリードで遺伝子特異的プライマー配列を保持する他のアンプリコンベースのターゲットシーケンスシステムとも互換性がある。 AmpliSeqのデータは、例えば、遺伝子特異的プライマーがライブラリー生成の間に大部分切断されるため、適合性がない。重要な考慮点は、UNDR-ROVERは、ペアエンドリードのかなりのオーバーラップを示すシーケンシングデータを処理することを意図しているため、これを達成していないシステムでの使用はお勧めしない。
インストール
mac os 10.13のvirtualenv環境でテストした。
公式ではvirtualenvを使い、仮想環境で導入、実行することを勧めている。virtualenvを使ったインストールは以下のように行う。
virtualenv --system-site-packages undr_rover_dev source
undr_rover_dev/bin/activate
pip install git+https://github.com/bjpop/undr_rover
> undr_rover --version
$ undr_rover --version
undr_rover 2.0.0
> undr_rover
$ undr_rover
usage: undr_rover [-h] [--version] --primer_coords PRIMER_COORDS
--primer_sequences FILE
[--primer_prefix_size PRIMER_PREFIX_SIZE]
[--kmer_length KMER_LENGTH]
[--kmer_threshold KMER_THRESHOLD]
[--primer_bases PRIMER_BASES] [--proportionthresh N]
[--absthresh N] [--qualthresh N] [--overlap OVERLAP]
[--max_variants N] --reference FILE [--id_info FILE] --out
FILE [--log FILE] [--coverdir COVERDIR] [--fast]
[--snvthresh N] [--genotype] [--ploidy {1,2}]
[--error ERROR]
fastqs [fastqs ...]
undr_rover: error: too few arguments
(undr_rover_dev) kamisakBookpuro:undr_rover_dev kamisakakazuma$ undr_rover -h
usage: undr_rover [-h] [--version] --primer_coords PRIMER_COORDS
--primer_sequences FILE
[--primer_prefix_size PRIMER_PREFIX_SIZE]
[--kmer_length KMER_LENGTH]
[--kmer_threshold KMER_THRESHOLD]
[--primer_bases PRIMER_BASES] [--proportionthresh N]
[--absthresh N] [--qualthresh N] [--overlap OVERLAP]
[--max_variants N] --reference FILE [--id_info FILE] --out
FILE [--log FILE] [--coverdir COVERDIR] [--fast]
[--snvthresh N] [--genotype] [--ploidy {1,2}]
[--error ERROR]
fastqs [fastqs ...]
Find variants from fastqs via a mapping-free approach.
positional arguments:
fastqs Fastq files containing reads.
optional arguments:
-h, --help show this help message and exit
--version show program's version number and exit
--primer_coords PRIMER_COORDS
Primer coordinates in TSV format.
--primer_sequences FILE
Primer base sequences as determined by a primer
generating program.
--primer_prefix_size PRIMER_PREFIX_SIZE
Size of primer key for blocks in dictionary.
--kmer_length KMER_LENGTH
Length of k-mer to check after the primer sequence.
Defaults to 30.
--kmer_threshold KMER_THRESHOLD
Number of single nucleotide variants deemed acceptable
in kmer.
--primer_bases PRIMER_BASES
Number of bases from primer region to use in gapped
alignment.Helps with variant calling near the edges of
a block. Defaults to 5.
--proportionthresh N Keep variants which appear in this proportion of the
read pairs. For a given target region, and bin
otherwise. Defaults to 0.05.
--absthresh N Only keep variants which appear in at least this many
read pairs. Defaults to 2.
--qualthresh N Minimum base quality score (phred).
--overlap OVERLAP Minimum proportion of block which must be overlapped
by a read. Defaults to 0.9.
--max_variants N Ignore reads with greater than this many variants
observed. Defaults to 20.
--reference FILE Reference sequences in Fasta format.
--id_info FILE File containing rs ID information in VCF format.
--out FILE Name of output file containing called variants.
--log FILE Logs progress in specified file, defaults to stdout.
--coverdir COVERDIR Directory to write coverage files, defaults to current
working directory.
--fast Use gapped alignment less often, leading to faster run
time.
--snvthresh N Distance between two single nucleotide variants before
going to a gapped alignment.
--genotype Compute genotypes for SNVs. Defaults to False.
--ploidy {1,2} Ploidy for genotyping 1 = haploid, 2 = diploid.
Defaults to 2.
--error ERROR Expected base read error rate. Defaults to 0.002.
解析が終わったら"deactivate"で仮想環境を終了する。
ラン
テストラン。
git clone https://github.com/bjpop/undr_rover
cd undr_rover/example/
undr_rover --primer_coords example_coords.txt --primer_sequences \ example_primers.txt --reference hg19.fa --out example.vcf \
--genotype example_R1.fastq example_R2.fastq
出力はVCF形式に準じている。
引用
UNDR ROVER - a fast and accurate variant caller for targeted DNA sequencing.
Park DJ, Li R, Lau E, Georgeson P, Nguyen-Dumont T, Pope BJ
BMC Bioinformatics. 2016 Apr 16;17:165.