macでインフォマティクス

macでインフォマティクス

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

DuplicationとdeletionのSVコールから偽陽性の可能性が高いコールをフィルタリングする duphold

2019 5/2 論文追記

 

 構造変異(SV)は、重複、欠失、逆位、挿入、および転座を含む広範な種類の変異である。 SVは、一塩基変異および挿入欠失変異よりも高精度で検出することがより困難であることが知られている。そのため、偽陽性率が高くなる可能性があることからノイズと関心のある真の変異を分離することは困難である。最も一般的に使用されている構造変異コーラーref.1–5は、2種類のシーケンスアラインメントを使用して構造のバリエーションを発見する。異常な方向またはインサートサイズを持つペアエンドリード(いわゆる "不一致ペア")のアライメントと、ゲノムの異なる部分にアライメントされるsplitアライメントである。これらのメソッドはうまく機能するが、アライメントされたデプスの深さを使用しない。これは重要な制限である。例えば、二倍体でシーケンスカバー率50%は真のヘテロ接合の欠失と予想できる。構造変異の正確さを評価することを可能にするSVPlaudit(ref.6)による筆者らの何千もの候補SVの評価経験から、信頼性の高い欠失と重複のコールと偽陽性を区別する2つの一貫したパターンがあることに気づいた。第一に、明らかなデプスの減少または増加のないイベントは、人間の目には本当のイベントとして現れる可能性ははるかに低い。第二にブレイクポイントで(または近くで)デプスが急激に変化するイベントはもっともらしい。明らかに偽陽性のコールは、これらのシグナルのどちらかまたは両方を欠いている。残念ながら、ほとんどのSV発見ツールはこれらの信号を発見するのに必要とされるカバレッジ情報を調べない。そのため、本著者らはdupholdという手動の作業で得られた観察結果を強化し、SVコール結果にすばやくアノテーションを付け高品質のバリアントコールを優先させる方法を開発した。

 dupholdは、前述のhts-nim(ref.7)ライブラリを使用して、BAMまたはCRAMファイルからカバレッジ情報をすばやく抽出して配列(arrayのこと)に収納し、そこで直接問い合わせることができる。これらのデプスプロファイルは、アライメントのBAMまたはCRAMファイルから計算されたカバレッジを持つ構造変異のVCFファイルにすばやくアノテーションを付けるために使用される。

(以下略)

 

twitter


インストール

依存

  • htslib(HTSのCライブラリ)のlibhts.soが必要。htslibのビルドに必要なライブラリやコンパイラ等は以前紹介したsamtoolsのインストールを参照(リンク)。
  • lumpyのコール
#htslib1.9のダウンロードとビルド
wget https://github.com/samtools/htslib/releases/download/1.9/htslib-1.9.tar.bz2
tar xjf htslib-1.9.tar.bz2
cd htslib-1.9/
./configure --prefix=/where/to/install
make -j 8
make install
#LD_LIBRARY_PATHをhtslibに通す。
export LD_LIBRARY_PATH=/where/to/installed/htslib-1.9/

#lumpyはcondaを使えばワンライナーで入る。
conda install -c bioconda lumpy-sv

本体 Github

binaryをダウンロードして実行権をつける。

#v1.02をダウンロード
wget https://github.com/brentp/duphold/releases/download/v0.1.2/duphold
chmod a+x duphold

> ./duphold -h

# ./duphold -h

version: 0.1.2

 

  Usage: duphold [options]

 

Options:

  -v --vcf <path>           path to sorted SV VCF/BCF

  -b --bam <path>           path to indexed BAM/CRAM

  -f --fasta <path>         indexed fasta reference.

  -s --snp <path>           optional path to snp/indel VCF/BCF with which to annotate SVs. BCF is highly recommended as it's much faster to parse.

  -t --threads <int>        number of decompression threads. [default: 4]

  -o --output <string>      output VCF/BCF (default is VCF to stdout) [default: -]

  -d --drop                 drop all samples from a multi-sample --vcf *except* the sample in --bam. useful for parallelization by sample followed by merge.

  -h --help                 show help

 

実行方法

 1、dupholdでSVのVCFにカバレッジ情報を追加する。

bam/cram、リファレンスfastaSNPs/indelコールのvcf、lumpyのSVコールのVCFを指定してランする。VCFの代わりにBCFを使うとずっと速く動作する。SNPs/indelのVCF/BCFは複数サンプルをジョイントしたものでもよい。

duphold -s gatk.vcf -t 4 -v SV.vcf -b input.bam -f ref.fasta -o output.bcf
  • -s    optional path to snp/indel VCF/BCF with which to annotate SVs. BCF is highly recommended as it's much faster to parse
  • -v    path to sorted SV VCF/BCF
  • -b    path to indexed BAM/CRAM
  • -f     indexed fasta reference
  • -t     number of decompression threads. [default: 4]

bcfの方が高速動作するため、BCFに変換してから使用する事が強く推奨されている。

SV VCFのFORMAT fieldに以下のいずれかの情報が追加される。

  • DHFC: fold-change for the variant depth relative to the rest of the chromosome the variant was found on
  • DHBFC: fold-change for the variant depth relative to bins in the genome with similar GC-content
  • DHFFC: fold-change for the variant depth relative to Flanking regions

 Githubより

 

さらにINFOフィールドにはコールのGC fractionが追加される。 また、SNP/indelのVCFも指定されると以下の情報も追加される。

  • DHET: counts of SNP heterozygotes in the SV supporting: [0] a normal heterozygote, [1] a triploid heterozygote. for a DUP, we expect most hets to have an allele balance closer to 0.33 or 0.67 than to 0.5. A good heterozygous DUP will have larger values of [1], than [0], though it's possible there are no HETs in small events

  • DHHU: counts of [0] Hom-ref, [1] Hom-alt, [2] Unknown variants in the event. A heterozygous deletion may have more hom-alt SNP calls. A homozygous deletion may have only unknown SNP calls

 Githubより

 

2、bcftools viewにより高品質のSVのみ検出する。

bcftools view -i '(SVTYPE = "DEL" & FMT/DHFFC[0] < 0.7) | (SVTYPE = "DUP" & FMT/DHBFC[0] > 1.3)" $svvcf

 

引用

duphold: scalalable, depth-based annotation and curation of high-confidence structural variant calls
Brent S. Pedersen, Aaron R. Quinlan,

bioRxiv preprint first posted online Nov. 8, 2018

 

追記

Duphold: scalable, depth-based annotation and curation of high-confidence structural variant calls
Brent S Pedersen Aaron R Quinlan
GigaScience, Volume 8, Issue 4, Published: 24 April 2019

 

関連(lumpyのコール)