macでインフォマティクス

macでインフォマティクス

NGS関連のインフォマティクス情報についてまとめています。

(ヒトゲノム)高速かつ精度の高いロングリードのSVコーラー cuteSV

 

 構造変化(SV)とは、欠失、挿入、逆位、重複、転座などのゲノムリアレンジメントで、その大きさが50 bpを超えるものを指す。ヒトゲノム上で最大のdivergencesとして、SV はヒトの疾患(遺伝性疾患やガンなど)、進化(遺伝子欠損やトランスポゾン活性など)、遺伝子制御、その他の表現型(交配や本質的な生殖隔離など)と密接に関連している。

 ショートリードベースのSVコールアプローチを開発するための努力がなされてきた。それらの多くは、リードデプス、不一致リードペア、スプリットリードアラインメント、ローカルアセンブリ、またはそれらの組み合わせなどの方法を用いており、1000 Genomes Projectのような大規模なゲノム研究において重要な役割を果たしてきた。しかし、これらのツールは比較的短いリード長であるため、高感度なSV検出には限界があり、偽陽性も存在する。

 Pacific Bioscience (PacBio)や Oxford Nanopore Technology (ONT)などのロングリードシークエンシング技術の急速な発展に伴い、ロングレンジスパニング情報は、より高分解能でより包括的に SV を検出する機会を提供している。しかし、高いシーケンシングエラー率(通常5~20%)とリード長(平均10kbp以上)に対応するためには、新たな計算アプローチが必要である。これまでの研究では、主に2種類のアプローチが採用されてきた。

 デノボアセンブリーベースのアプローチ は、より長いゲノム配列(コンティグやスキャホールド)にリードをアセンブリし、アセンブリした配列とリファレンス配列とのアラインメントからSVを発見することを目的としている。このようなアプローチは、リードアライメントに基づくアプローチに比べて、リファレンスの影響を受けにくく、特にリードアライメントのアーチファクトがない。しかし、これらのアプローチは通常計算量が多く、大規模なゲノムのハプロタイプ配列の再構成にはまだいくつかの困難がある。

 リードアライメントに基づくアプローチは、リードをリファレンスに対して直接アライメントし、アライメント結果を解析することでSVを検出する。このようなアプローチは、感度を欠くことなく計算リソースに対してより費用対効果が高く、ロングリードベースのSVコーラーで広く使用されている。いくつかのリードアライメントベースのSVコーラーが提案されており、PB-Honey [ref.31]、SMRT-SV [ref.32]、Sniffles [ref.33]、PBSV (https://github.com/PacificBiosciences/pbsv)、およびSVIM [ref.34]などがある。これらの手法では、発散性の高いアラインメントを持つ局所的なゲノム領域の同定、クリッピングされたリード部分の局所的なアセンブリと再アラインメント、SV-スパンニングシグネチャクラスタリングなど、リードアラインメントによって暗示されるSVの証拠を見つけるために様々な方法を使用している。さらに、通常、BLASR [ref.36]、NGMLR [ref.33]、Minimap2 [ref.37]、PBMM2 (https://github.com/PacificBiosciences/pbmm2) などの最先端のロングリードアライナーがリードアライメントに使用されていた。

 しかし、リードアライメントに基づくSVコールは、依然として非自明である。シーケンシングエラーが多く、複雑なSVが存在する状況下では、SVブレークポイント周辺のリードのアライメントはキメラ的で不均質であり、通常は感度や精度が低い。そのため、リードアラインメントから得られるSVシグネチャは非常に複雑であり、様々な種類のSVに対して感度の高い検出を行うためには、それらを収集・解析することが困難である。主に、最先端のツールでは、以下のような技術的な課題がある。(1) 全体的に感度がまだ満足のいくものではない。(2) いくつかのアプローチ(rMETL [ref.38]、rCANID [ref.39]、npInv [ref.40]など)は、特定の設計のために、サブセットまたは特定のクラスのSVしか検出できない。(3) いくつかのアプローチ(PBSVやSMRT-SVなど)は、まだ時間がかかり、スケーリング性能が良くないため、多くの大規模データセットには適していない可能性がある;そして(4) いくつかのアプローチ(SMRT-SVやPB-Honeyなど)は、1種類のシーケンシングデータしかサポートしていない(例えばPacBioのリードのみ)ものもある。このような欠点が、ロングリードシーケンシングデータを広く利用する上でのボトルネックとなっている。 

 ここでは、いくつかの利点を持つ多用途なSV検出アプローチであるcuteSVを紹介する。(1) cuteSVは、最先端のSVコーラーよりも高いSV検出率を持っている。特に、カバレッジの低いデータに対しても、精度を落とさずに高感度に検出することができる。(2) cuteSVは、主流のロングリードシーケンシングプラットフォームで作成された様々なエラー率のデータセットをサポートしており、様々なタイプのSV(欠失、挿入、重複、逆位、転座を含む)を検出することができる。(3) cuteSVは、RAMの使用量が少なく、最先端のアプローチと同等の速度を持っている。さらに重要なことは、CPUのスレッド数に応じてほぼ直線的に高速化できるスケーラビリティを持っている。これらの特徴により、cuteSVは大規模なデータ解析に適しており、最先端のゲノミクス研究への可能性を秘めている。

cuteSVの概要
ソートされたBAMファイルを入力として、cuteSVは、アライメントされたリード中の大きな挿入/欠失や分割アライメントをSVシグネチャクラスターとして抽出し、解析してSVをコールする。このアプローチには大きく分けて以下の3つのステップがある。

ステップ1:cuteSVは複数のシグネチャ抽出手法を用いて、様々なタイプのSVのシグネチャを網羅的に収集する。さらに、挿入と欠失をヒューリスティックに組み合わせて、脆弱なアライメントから本物のSVの証拠を回収する。

ステップ2: cuteSVは、特別に設計されたクラスタリングと精緻化のアプローチを用いて、キメラアラインメントされたリードを局所領域にクラスタリングし、さらにクラスタを精緻化して、SVのシグネチャヘテロ接合型SVと正確に区別する。

ステップ3:cuteSVは、いくつかのオーダーメイドのルールを使用して、洗練されたSVシグネチャクラスターに基づいてSVコールとジェノタイピングを実行する。

 

 

インストール

依存

  • python3
  • pysam
  • Biopython
  • cigar
  • numpy

Github

#bioconda (link)
conda install -c bioconda cutesv -y

#pip
pip install cuteSV

cuteSV -h

$ cuteSV -h

usage: cuteSV [-h] [--version] [-t THREADS] [-b BATCHES] [-S SAMPLE]

              [--retain_work_dir] [-p MAX_SPLIT_PARTS] [-q MIN_MAPQ]

              [-r MIN_READ_LEN] [-md MERGE_DEL_THRESHOLD]

              [-mi MERGE_INS_THRESHOLD] [-s MIN_SUPPORT] [-l MIN_SIZE]

              [-L MAX_SIZE] [-sl MIN_SIGLENGTH] [--genotype]

              [--gt_round GT_ROUND]

              [--max_cluster_bias_INS MAX_CLUSTER_BIAS_INS]

              [--diff_ratio_merging_INS DIFF_RATIO_MERGING_INS]

              [--diff_ratio_filtering_INS DIFF_RATIO_FILTERING_INS]

              [--max_cluster_bias_DEL MAX_CLUSTER_BIAS_DEL]

              [--diff_ratio_merging_DEL DIFF_RATIO_MERGING_DEL]

              [--diff_ratio_filtering_DEL DIFF_RATIO_FILTERING_DEL]

              [--max_cluster_bias_INV MAX_CLUSTER_BIAS_INV]

              [--max_cluster_bias_DUP MAX_CLUSTER_BIAS_DUP]

              [--max_cluster_bias_TRA MAX_CLUSTER_BIAS_TRA]

              [--diff_ratio_filtering_TRA DIFF_RATIO_FILTERING_TRA]

              [BAM] output work_dir

 

Long read based fast and accurate SV detection with cuteSV.

 

Current version: v1.0.6

Author: Tao Jiang

Contact: tjiang@hit.edu.cn

 

Suggestions:

 

For PacBio CLR/ONT data:

--max_cluster_bias_INS 100

--diff_ratio_merging_INS 0.2

--diff_ratio_filtering_INS 0.6

--diff_ratio_filtering_DEL 0.7

For PacBio CCS(HIFI) data:

--max_cluster_bias_INS 200

--diff_ratio_merging_INS 0.65

--diff_ratio_filtering_INS 0.65

--diff_ratio_filtering_DEL 0.35

 

 

 

positional arguments:

  [BAM]                 Sorted .bam file form NGMLR or Minimap2.

  output                Output VCF format file.

  work_dir              Work-directory for distributed jobs

 

optional arguments:

  -h, --help            show this help message and exit

  --version, -v         show program's version number and exit

  -t THREADS, --threads THREADS

                        Number of threads to use.[16]

  -b BATCHES, --batches BATCHES

                        Batch of genome segmentation interval.[10000000]

  -S SAMPLE, --sample SAMPLE

                        Sample name/id

  --retain_work_dir     Enable to retain temporary folder and files.

 

Collection of SV signatures:

  -p MAX_SPLIT_PARTS, --max_split_parts MAX_SPLIT_PARTS

                        Maximum number of split segments a read may be aligned

                        before it is ignored.[7]

  -q MIN_MAPQ, --min_mapq MIN_MAPQ

                        Minimum mapping quality value of alignment to be taken

                        into account.[20]

  -r MIN_READ_LEN, --min_read_len MIN_READ_LEN

                        Ignores reads that only report alignments with not

                        longer than bp.[500]

  -md MERGE_DEL_THRESHOLD, --merge_del_threshold MERGE_DEL_THRESHOLD

                        Maximum distance of deletion signals to be merged.[0]

  -mi MERGE_INS_THRESHOLD, --merge_ins_threshold MERGE_INS_THRESHOLD

                        Maximum distance of insertion signals to be

                        merged.[100]

 

Generation of SV clusters:

  -s MIN_SUPPORT, --min_support MIN_SUPPORT

                        Minimum number of reads that support a SV to be

                        reported.[10]

  -l MIN_SIZE, --min_size MIN_SIZE

                        Minimum size of SV to be reported.[30]

  -L MAX_SIZE, --max_size MAX_SIZE

                        Maximum size of SV to be reported.[100000]

  -sl MIN_SIGLENGTH, --min_siglength MIN_SIGLENGTH

                        Minimum length of SV signal to be extracted.[10]

 

Computing genotypes:

  --genotype            Enable to generate genotypes.

  --gt_round GT_ROUND   Maximum round of iteration for alignments searching if

                        perform genotyping.[500]

 

Advanced:

  --max_cluster_bias_INS MAX_CLUSTER_BIAS_INS

                        Maximum distance to cluster read together for

                        insertion.[100]

  --diff_ratio_merging_INS DIFF_RATIO_MERGING_INS

                        Do not merge breakpoints with basepair identity more

                        than [0.2] for insertion.

  --diff_ratio_filtering_INS DIFF_RATIO_FILTERING_INS

                        Filter breakpoints with basepair identity less than

                        [0.6] for insertion.

  --max_cluster_bias_DEL MAX_CLUSTER_BIAS_DEL

                        Maximum distance to cluster read together for

                        deletion.[200]

  --diff_ratio_merging_DEL DIFF_RATIO_MERGING_DEL

                        Do not merge breakpoints with basepair identity more

                        than [0.3] for deletion.

  --diff_ratio_filtering_DEL DIFF_RATIO_FILTERING_DEL

                        Filter breakpoints with basepair identity less than

                        [0.7] for deletion.

  --max_cluster_bias_INV MAX_CLUSTER_BIAS_INV

                        Maximum distance to cluster read together for

                        inversion.[500]

  --max_cluster_bias_DUP MAX_CLUSTER_BIAS_DUP

                        Maximum distance to cluster read together for

                        duplication.[500]

  --max_cluster_bias_TRA MAX_CLUSTER_BIAS_TRA

                        Maximum distance to cluster read together for

                        translocation.[50]

  --diff_ratio_filtering_TRA DIFF_RATIO_FILTERING_TRA

                        Filter breakpoints with basepair identity less than

                        [0.6] for translocation.

 

実行方法

bamファイルを指定する。

mkdir work_dir
cuteSV -t 20 input.bam output.vcf work_dir

vcf4.2 formatのSVコールが出力される。

引用

Long-read-based human genomic structural variation detection with cuteSV

Tao Jiang, Yongzhuang Liu, Yue Jiang, Junyi Li, Yan Gao, Zhe Cui, Yadong Liu, Bo Liu & Yadong Wang
Genome Biology volume 21, Article number: 189 (2020)