macでインフォマティクス

macでインフォマティクス

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

Varscan2 の解析の流れ

修正 不確かな情報を削除

2019 2/17 誤字修正

 

Using VarScan 2 for Germline Variant Calling and Somatic Mutation Detection(Daniel C. Koboldt et al., 2013)より

  

 シングルヌクレオチド変異(SNV)および小さな挿入/欠失(indels)のようなバリアントの同定は、リシークエンシングデータの解析における重要なステップである。次世代シークエンシング(NGS)機器の驚異的なスループットは、従来のキャピラリーベースのシークエンシング(Mardis 2008)のコストのほんの一部かつ比較的短時間で、数千人のエクソームまたはゲノムのシーケンシングを可能にした。しかし、カバレッジの変動、シークエンシングエラー、アライメントアーチファクト、およびその他の問題が原因で、正確なバリアントコールが依然として課題となっている(Koboldt、Ding他2010)。プールされたサンプルまたは腫瘍 - ノーマルペアのシーケンシングのようなNGSの多くの応用は、二倍体サンプルのために設計された確率的方法を使用して検出することがより困難な低頻度で変異体の検出を必要とする。

 著者らは、以前に個人またはプールサンプルのSNV / indelコールを行うVarScan(http://varscan.sourceforge.net)を開発した(Koboldt、Chen et al、2009)。その後のリリースでは、tomor- normal比較(Koboldt、Zhang et al。2012)におけるsomatic変異(図1)およびコピー数変化(CNA)の検出能力、およびfamily trioにおけるde novo 突然変異検出まで拡張している。

 ここでは、VarScan 2を使用した生殖細胞系列、体細胞系列およびtrioの変異検出戦略に関する推奨プロトコルを提示する。基本プロトコル1は、VarScanを用いて単一の個体または一緒に(プール)シーケンシングされた個体のSNVおよびindelsをコールする方法を説明する。代替プロトコル1はまた、個体群のコホート、またはクロスサンプルバリアントコールしても知られている生殖系列のバリアントコールも記述する。基本プロトコル2は、体細胞突然変異およびLOHイベント(wiki)の検出を含むtumor -normalペアの分析のためのステップを概説する。代替プロトコル1は、tumor -normalペアにおいて体細胞性CNAを同定するためのexomeベースのコピー数変化コールを記述する。基本プロトコール3は、family trioにおける生殖系列変異およびデノボ変異の同定について記載する。サポートプロトコール1では、生殖系列、体細胞および新生児の突然変異のバリアントコールセットに対して推奨される偽陽性のバリアントフィルタリング戦略について概要を説明する。

 

発表から時間は経ってますが、複数のケースに対応できる高感度なバリアントコーラーで、tumor検出などで今でも使われているツールです。他のパイプラインに組み込まれることも多く、引用回数も多くなっています。

 

f:id:kazumaxneo:20180820214531p:plain

3つの解析例。論文より転載。

  

公式ページ(NGSの基本ツールへのリンクあり)

http://varscan.sourceforge.net/using-varscan.html#v2.3_pileup2snp

wiki

https://wiki.bits.vib.be/index.php/Varscan2

Varscan2に関するツイート

 インストール

mac os10.13に導入した。

Anaconda環境ならcondaで導入できる。

conda install -c bioconda varscan

https://bioconda.github.io/recipes.html#recipes

> varscan

$ varscan

VarScan v2.4.3

 

***NON-COMMERCIAL VERSION***

 

USAGE: java -jar VarScan.jar [COMMAND] [OPTIONS] 

 

COMMANDS:

pileup2snp Identify SNPs from a pileup file

pileup2indel Identify indels a pileup file

pileup2cns Call consensus and variants from a pileup file

mpileup2snp Identify SNPs from an mpileup file

mpileup2indel Identify indels an mpileup file

mpileup2cns Call consensus and variants from an mpileup file

 

somatic Call germline/somatic variants from tumor-normal pileups

copynumber Determine relative tumor copy number from tumor-normal pileups

readcounts Obtain read counts for a list of variants from a pileup file

 

filter Filter SNPs by coverage, frequency, p-value, etc.

somaticFilter Filter somatic variants for clusters/indels

fpfilter Apply the false-positive filter

 

processSomatic Isolate Germline/LOH/Somatic calls from output

copyCaller GC-adjust and process copy number changes from VarScan copynumber output

compare Compare two lists of positions/variants

limit Restrict pileup/snps/indels to ROI positions

 

condaでバージョン指定せずに導入すると2.4.3が入った。

dockerイメージも複数上がっている。

https://hub.docker.com/search/?isAutomated=0&isOfficial=0&page=1&pullCount=0&q=varscan&starCount=0

オフィシャルイメージのpullは"docker pull jeltje/varscan2"

 

ラン

1、Germline(Basic Protocol 1)

VarScanはシーケンスカバレッジ、サポートするリード数、バリアントアレル頻度(VAF)、Fisherの正確検定のp値、に関してユーザー定義の最小基準を満たすバリエーションを探す。

germlineのバリアント検出には3つのコマンドが使用できる。

  • mpileup2snp - calls single nucleotide polymorphisms (SNPs)
  • mpileup2indel - calls insertions and deletions (indels)
  • mpileup2cns - calls a consensus genotype (reference, SNP, or indel) 

パラメータ

samtools mpileup –B –q 1 –f ref.fa sample.bam > sample.mpileup

#sambambaを使うなら
sambamba mpileup sample.bam -t 12 --samtools –B –q 1 -f ref.fa > sample.mpileup

#SNV
varScan mpileup2snp sample.mpileup –-min-coverage 10 –-min-var-freq \
0.20 –-p-value 0.05 > sample.varScan.snp

#indel
varScan mpileup2indel sample.mpileup -–min-coverage 10 –-min-var-freq \
0.10 -–p-value 0.10 > sample.varScan.indel

#Filter SNP calls to remove those near indel positions
varScan filter sample.varScan.snp –-indel-file sample.varScan.indel \
–-output-file sample.varScan.snp.filter

#Filter indels to obtain a higher-confidence set
filter sample.varScan.indel -–min-reads2 4 –-min-var-freq 0.15 \
–-p-value 0.05 –-output-file sample.varScan.indel.filter

 

2、Alternate Protocol 1

 複数サンプルのパイルアップ(mpileup)出力を生成するSAMtoolsの機能により、多数のサンプルセットにわたってバリアントを同時にコールすることができる。クロスサンプルバリアントコールと呼ばれるこのアプローチは、サンプルグループ内のすべてのバリアント位置を特定し、それぞれのサンプルすべてに対して genotype callsを提供する。各サンプルがすべてのバリアントの位置でジェノタイピングされなければならないので、出力サイズは含まれるサンプルの数とともに指数関数的に増加する。しかし、コホート内のすべてのバリアントコールを含む単一のマスターVCFファイルは、ダウンストリーム分析や共有に便利なのでよく使われる。コマンドとパラメータの設定は、基本プロトコル1の場合と同じである。

#pileup(*1)
samtools mpileup –B –q 1 –f reference.fasta sample1.bam sample2.bam sample3.bam … > cohort.mpileup


#SNV
varscan mpileup2snp cohort.mpileup –-vcf-sample-list samples.txt \
–-min-coverage 10 –-min-var-freq 0.20 –-p-value 0.05 –-output-vcf 1 \
>cohort.varScan.snp.vcf

#indel
varscan mpileup2indel cohort.mpileup –-vcf-sample-list samples.txt \
-–min-coverage 10 –-min-var-freq 0.10 -–p-value 0.10 –-output-vcf 1 \
>cohort.varScan.indel.vcf

 

3、Somatic(Basic Protocol 2)

 腫瘍サンプルと対照となる(一致した)正常サンプルの次世代シーケンシングは、ガンの開始、発達、の遺伝的基礎を研究するための強力な戦略である。獲得された体細胞突然変異(バリアントの0.01%以下)を先天的な生殖系列変異( 99.99%以上のバリアント)と区別することができるので、患者からの一致する正常サンプルの重要性は過小評価することはできない。分析に影響を与える可能性のあるバッチ効果を最小限に抑えるため、腫瘍および正常サンプルは同一実験条件下でシーケンシングすることが推奨される。Tumor-normal ペアのNGSデータについても、VarScanは、十分なカバレッジであれば、正常および腫瘍データを直接比較することにより、体細胞突然変異、生殖系列変異およびLOHイベントを同定することができる。 VarScanの体細胞コマンドは、正常および腫瘍サンプルからのmpileup入力を(この順序で)受け入れる。最小カバレッジ要件を満たすあらゆる位置で、可能なバリアントを識別するために両方のサンプルを独立にコールする。いずれかのサンプルでバリアントがコールされた場合、両方のサンプルの遺伝子型および支持データを比較して、体細胞状態(生殖細胞系、LOH、または体細胞性)を推測する。

 具体的には、正常および腫瘍におけるリファレンスおよび変異対立遺伝子を支持するリード数のフィッシャーexact testは、バリアントの証拠が正常と腫瘍とで異なる(および統計的に有意である)かどうかを決定する。このアプローチは、サブクローナル集団または低い腫瘍細胞性に起因する低頻度突然変異のような、正常と腫瘍との間の微妙な違いの検出を可能にする。

パラメータ

#tumorとnormalをこの順序で指定してpileup
samtools mpileup –B –q 1 –f ref.fa normal.bam tumor.bam > normal-tumor.mpileup

#Run VarScan in somatic mode
varscan somatic normal-tumor.mpileup output \
–min-coverage 10 –min-var-freq 0.08 –somatic-p-value 0.05

#Run the processSomatic subcommand to divide the output into separate files based on somatic status and confidence
varscan processSomatic output.basename.snp
varscan processSomatic output.basename.indel

#上記のコマンドを実行すると6つのファイルに分けて出力される。そのうち.hcの拡張子のファイルは empiricalな基準でhigh-confidenceと判断されたバリアントが出力されている。例えば、高い信頼性を有する体細胞突然変異は、tumor VAF> 15%、normal VAF <5%、およびsomatic p value <0.03を有する。ユーザーが調節可能。

#high-confidenceのsomaticバリアントに対して追加のフィルタリングを実行。
#indelsの近くのアラインメントの問題のために誤検出される可能性のある体細胞突然変異を除去。上記のsnpとindelの出力を指定する
varscan somaticFilter output.snp.Somatic.hc –indel-file output.indel –output-file output.snp.Somatic.hc.filter

 

4、Support Protocol 1

著者らは、ストランド、mismatch qualityの合計、リードポジションバイアス、マッピングクオリティ、および他の基準(Koboldt、Zhang et al、2012)に基づく、false positiveの可能性の高いバリアントをフィルタリングする戦略を開発: bam-readcountを開発した(紹介)。bam-readcountはBAMファイルを使用する。体細胞変異の場合は、tumorのBAMファイルを使用する。生殖細胞およびLOHイベントについては、normalのBAMファイルを使用する。

varscanのtab出力ファイルを入力とする。

bam-readcount –q 1 –b 20 –f ref.fa –l varScan.variants input.bam > variants.readcounts

perl fpfilter.pl varScan.variants variants.readcounts –output-basename variants.filter

 

5、Family Trios(Basic Protocol 3)

 ファミリーpedigreesの次世代シーケンシングは遺伝性疾患を研究する強力なアプローチを提供する。家族のトリオ(母親、父親、および病に冒された子供)を研究して、伝わった対立遺伝子または新規の変異を同定し、疾患の感受性に関与する可能性のある遺伝子を絞り込む。著者らは、バリアントコール精度を向上させ、明らかなMendelian Inheritance Errors(MIE)を特定し、信頼性の高いde novo突然変異を検出するために家族関係を活用するVarScanのTrio解析モジュールを開発した。家族の全ゲノムシーケンシングからのデータは、ヒトにおける新規の突然変異率は、約1.1 × 10^-8 per haploid genome、であることを示している(Roach、Glusman et al、1000 Genomes Project Consortium 2010)。この推定により、個体の二倍体ゲノムは平均して32億塩基対のうちの約64のde novo変異を有する。現在のヒトコンセンサスコード領域配列(CCDS)のサイズ(約3,400万bp)を考えると、二倍体個体あたりのデノボコード変異は1つ未満であると予想される。

 この極端な希少性のために、デノボ変異は保守的に呼ばれるべきである。これに対処するために、VarScanは、緩和されたパラメータを用いて各親の明らかなde novo突然変異を再評価し、一方または両方の親の生殖系列変異としてのいくつかの証拠を有するものを再分類する。同様の方法で、VarScanは見かけ上のMendelian Inheritance Errorsを調整しようとする。Trioサブコマンドの出力は単一のVCFで、すべてのバリアントは生殖細胞系(送信または非送信)、de novo、またはMIEのいずれかに分類されている。

#Run SAMtools mpileup on father, mother, and child BAM files in that order
samtools mpileup –B –q 1 –f ref.fa father.bam mother.bam child.bam > trio.mpileup

#Run VarScan mpileup2snp to call SNVs. 000747096 for definitions and recommended values of command-line parameters
varscan trio trio.mpileup trio.mpileup.output --min-coverage 10 –-min-var-freq 0.20 –-p-value 0.05 –adj-var-freq 0.05 –adj-p-value 0.15

SNPsの trio.mpileup.output.snp.vcfと、indelの trio.mpileup.output.indel.vcfが出力される。

 

 

論文では、これ以外にsomatic CNVの解析例(Alternate Protocol 2)などが説明されています。また、論文のStrategic Planning ~ Basic Protocolチャプターでは、Varscan2の基本的な使い方からより実践的プロトコルまで丁寧に説明されています。出力のvarscan nativeとVCF準拠フォーマットについての表もあるので確認してみてください(リンク)。

 

 

類似ツール

  引用

Using VarScan 2 for Germline Variant Calling and Somatic Mutation Detection
Daniel C. Koboldt, David E. Larson, Richard K. Wilson

Curr Protoc Bioinformatics. Author manuscript; available in PMC 2014 Dec 29.

 

VarScan 2: somatic mutation and copy number alteration discovery in cancer by exome sequencing
Koboldt DC1, Zhang Q, Larson DE, Shen D, McLellan MD, Lin L, Miller CA, Mardis ER, Ding L, Wilson RK

Genome Res. 2012 Mar;22(3):568-76.


VarScan: variant detection in massively parallel sequencing of individual and pooled samples
Daniel C. Koboldt,* Ken Chen, Todd Wylie, David E. Larson, Michael D. McLellan, Elaine R. Mardis, George M. Weinstock, Richard K. Wilson, and Li Ding

Bioinformatics. 2009 Sep 1; 25(17): 2283–2285.

 

A review of somatic single nucleotide variant calling algorithms for next-generation sequencing data.

https://www.ncbi.nlm.nih.gov/pubmed/29552334

 

Detailed comparison of two popular variant calling packages for exome and targeted exon studies.

https://www.ncbi.nlm.nih.gov/pubmed/25289185

 

参考

*1

論文にはusers may encounter I/O issues when processing many (>100) BAM files simultaneously、とある

 

ほか資料

 

 ガン変異コーラーの被引用数

https://cell-innovation.nig.ac.jp/surfers/cancer_caller.html

  

 


上記とは直接関係ありませんが、pileup2snpやpileup2indelを走らすと、vcf出力にエラーがあるらしく、--output-vcf 1をつけてもvcfのヘッダーが出力されませんでした。ヘッダーがないとbzip圧縮してvcf-toolsなどのフィルタリングができなくなります。ヘッダーが正常に出力されるmpileup2vcfコマンドを走らせ、得られたvcfをSNVとindelにパースすることで、エラーは回避できます。

VarScan / Discussion / Help:VCF header is missing in VARSCAN2 output