構造変化（SV）はドナーゲノムの構造変化をもたらすゲノム変異である。Indels、コピー数変動（CNV）およびゲノム再編成はすべてSVのサブクラスである。多くの研究は、SVが正常なヒト集団[論文より ref.1,2]ならびに癌ゲノム[ref.3-5]において広く広がっていることを明らかにした。ハイスループットシーケンシング（HTS）はSV検出のための主要なプラットフォームとなっており、HTSデータを用いたSV検出のための多くのアルゴリズムが開発されている[ref.6]。 HTSデータに基づくゲノム研究において、重要な問題は、異なるシナリオにおける様々なアルゴリズムの性能をベンチマークすることである。アルゴリズムの性能は、その設計、その実装、ならびにシーケンシングデータの品質および機能に依存する。ベンチマークのための理想的な方法は、シーケンシングとそれに続く実験的検証によるものだが、この方法は高価で、労働集約的で時間のかかるものである。したがって、HTSデータのシミュレーションは、ゲノム変動検出アルゴリズムをベンチマークするための強力かつ費用対効果の高い代替方法となる。
通常、シミュレーションによるHTSベースのSV検出アルゴリズムのベンチマーキングには、（1）シミュレートされたSVを含むゲノムの生成、および（2）シミュレートされたSVを有するゲノムに基づくHTSショートリードのシミュレーションが含まれる。実際のHTSデータを近似するために、生成されたゲノムは実際のゲノムと類似しているはずであり、シミュレートされたHTSデータは様々なシークエンシングエラーおよびバイアスを含むべきである。腫瘍サンプルはしばしば正常な汚染および異種サブクローナル腫瘍細胞を含むので、腫瘍ゲノムのシミュレーションデータは正常な汚染および/または複数のサブクローンを含むべきである。
RSVsim [ref.7]、SCNVsim [ref.8]、VarSim [ref.9]、IntSIM [ref.10]およびSInC [ref.11]を含むいくつかのSVシミュレーションパッケージが利用可能である。これらのパッケージは、コミュニティに大きなリソースを提供する。しかし、著者らが知る限り、SVベンチマークのためのデータシミュレーションにおけるすべての問題を体系的に解決したパッケージはない。例えば、RSVsimは広範囲のSVをシミュレートすることができるが、CNVを生成することはできず、normal contaminationとsubclonesで腫瘍データを生成することはできない。 SCNVsimは、aneuploidy、normal contamination、multiple-subclonesなどの腫瘍データをシミュレートする上での問題を考慮し体細胞SV / CNVでは腫瘍遺伝子を生成できたが、SNVのような事象は発生しない。 VarSimは、包括的なクラスのゲノム変異をシミュレートすることができるが、aneuploidy、normal contamination、multiple-subclonesで腫瘍データを生成することもできない。 IntSimおよびSInCは、生殖系列および体細胞変異の両方をシミュレートすることができるが、SNVとCNVのみシミュレートする。さらに、GCバイアスはHTSデータに広がっていることはよく知られている。シミュレータpIRS [ref.12]はシミュレーションデータにGCバイアスを導入することができるが、そのGCバイアスプロファイルは1つのデータセットで訓練され、ユーザーは基本的に1つのタイプのGCバイアスのみを生成できる。 Cancer Genome Atlasと1000ゲノムプロジェクトのデータから得られた何百ものシーケンシングデータを分析したところ、GCバイアスはさまざまな形を取っていた[ref.13]。また、pIRSはこれらのGCバイアスをシミュレートするほど柔軟ではない。
この論文では、SVシミュレーション用のPysim-svを紹介する。他のHTSデータシミュレーションパッケージと比較して、Pysimには主に3つの利点がある。
１、SNVsだけでなくSVの全スペクトルをシミュレートできる。
２、Aneuploidy、normal contamination、multiple-subclonesのサブクローンを有する腫瘍データのシミュレーションを可能にする。
３、任意の形式のGCバイアスを用いてHTSデータを生成できる。
ワークフロー。論文より転載。
インストール
依存
- ART（紹介）
本体 Github （マニュアルPDFも含まれる）
GitHub - xyc0813/pysim: A python package to simulate structural variation
git clone https://github.com/xyc0813/pysim.git
cd pysim/
> python specify_ploidy.py -h
$ python specify_ploidy.py -h
2018-04-26 13:18:19
simulation begins at:0.041636
Usage: specify_ploidy.py -i <file> -c <config file> -o <out_fasta>
specify_ploidy
Author: Yuchao Xia
Description: specify the number of ploidy for different chromesome
Options:
-h, --help show this help message and exit
-i FILE, --inFile=FILE
A reference fasta file.
-c FILE, --config=FILE
the number of ploidy of different chromesome
-o file, --output=file
output fasta file
——
> python pysim_main.py -h
$ python pysim_main.py -h
simulation begins at:0.053883
Usage: pysim_main.py -c ini.config
Author: Yuchao Xia
Description: simulate germline and somatic SVs.
Options:
-h, --help show this help message and exit
-c FILE, --config=FILE
config file
——
> python run_art.py -h
$ python run_art.py -h
2018-04-26 13:55:38
simulation begins at:0.048871
Usage: run_art.py -i <file> -l <read length> -f <coverage> -m <insert size> -s <sd of insert size> -o <out_fastq>
sim_SV v1.1
Author: Yuchao Xia
Description: run art
Options:
-h, --help show this help message and exit
-i FILE, --inFile=FILE
A reference fasta file.
-l 100 paired-end read length
-f 30 read coverage
-m 350 library insert size
-s 50 sd of insert size
-o file, --output=file
output fastq file
> python mixture_v0.5.py -h
$ python mixture_v0.5.py -h
2018-04-26 14:12:26
simulation begins at:0.074645
Usage: mixture_v0.5.py -i <file>
simulate_copy_number v1.1
Author: Yuchao Xia
Description:
Options:
-h, --help show this help message and exit
-i FILE, --inFile=FILE
config file.
-o file, --output=file
output fastq file
——
Usage: GC_bias.py -i <file> -g <GC_file> -b <name_sort_sam file> -o <out_fastq> -m <mode> -k <y/n>
GC v1.0.1
Author: Yuchao Xia
Description: accord GC content to filter paired-end reads from name sorted simulate bam file.
Options:
-h, --help show this help message and exit
-i FILE, --inFile=FILE
A reference fasta file.
-g FILE, --GC=FILE GC_function file,splited by ,the filst col is GC
content
-b FILE, --bam=FILE the name sorted sam/bam format file generated by ART
-o file, --output=file
output fastq file
-k y/n, --keepsam=y/n
generate sam file
-m int, --mode=int the function of GC content
ラン
流れは導入したい変異によって変わってくる。ここでは５つのコマンドを記載する。
１、ploidyを考慮したFASTAの調整。
倍数性のゲノムをシミュレートする時に行う。ランにはconifgファイルを使う。configはタブ区切りで以下のようなフォーマットにする。１行目は説明で"#"を使いコメントアウトし、２行目に情報を記載する。１列目はFASTAのヘッダー名、２列目がploidyで、diploidなら"2"。
複数クロモソームを指定するなら、同じフォーマットで３行以降に追記していく。
コマンドを走らせる時は、FASTAファイルと作成したconfigファイルを指定する。
python specify_ploidy.py -i input.fasta -c config1 -o output_ploid.fa
出力はFASTAファイルである。ploidyの数に応じてchromosome配列が複製される。上記のconfig条件では2なので、出力のFASTAのファイルサイズは２倍になる。
２、SNPs、SVを導入したgermlineとsomaticのFASTAの発生。
germlineのFASTAにはSVや新規indelは導入されない。一方、dbSNPのSNPsはgermlineとsomatic両方に導入される。またsomaticのcloneから生じたsubcloneには、追加の変異が導入される。
ランにはパラメータの詳細を記したconifgファイルを使う。
PDFより転載。
２行目のSVのconfigは以下のようなフォーマットにして別ファイル保存し、そのパスを上記のconfigに記載する。SVのconfigにはSVの種類とサイズ、導入回数をタブ区切りで記載する。SVの種類はdel、inv、translocation、tandem duplicationのいずれかを指定する。
３行目のリファレンスは一の出力。３行目のdbSNPのVCFは、sbSNPのSNPs情報。このSNPsはgermlineもsomatyicも導入される。それ以外のパラメータは上記のコメント部分を参照。
configファイルを指定してランする。
pysim_main.py –c SV.config
- -i A reference fasta file.
- -l 100 paired-end read length
- -f 30 read coverage
- -m 350 library insert size
- -s 50 sd of insert size
- -o utput fastq file
FASTAの他に、SVやSNPsの場所を記載したファイルができる。出力のprefixは "out_prex=test_chr22_new_"。subはサブクローン。クローンからさらに追加の変異が発生した細胞のゲノムのfastaも指定数発生する（sub_clone=）。サブクローンに新規出現した変異も別ファイルで保存される。
３、fastqのシミュレーション。例えば2で作ったFASTAを指定する。
python run_art.py -i test_chr22_new_somatic.fa -l 100 -f 30 -m 350 -s 50
- -i A reference fasta file.
- -l 100 paired-end read length
- -f 30 read coverage
- -m 350 library insert size
- -s 50 sd of insert size
- -o utput fastq file
出力はARTによるシミュレーションと同じ。
Tumorのsubclone とsomatic (normal)、またはsomaticとgermlineの比較などを行うならば、それぞれのFASTAを使い、順番にリードを発生させる必要がある。
４、mixしたfastqの作成。
python mixture_v0.5.py –i config
割合を示したconfigファイルを指定する。それぞれのsublconeの割合を指定する。２つsubclone FASTAをstep2で作り、それを体細胞のnormalと混ぜている。混ぜるなら、3の作業を行う必要はない。
５、GCバイアスの導入。分布はベルヌーイ分布に従う。step3またはstep4で出力されるsamを使うので、先にstep3かstep4を行なう。
python GC_bias.py -i input.fa -b sim_SV_genome.sam -o output.fa -m 2 -k n
- -i a reference fasta file
- -g GC_function file,splited by ,the filst col is GC content
- -b the name sorted sam/bam format file generated by ART
- -o output fastq file
- -m the function of GC content
- -k y/n, --keepsam=y/n
samを元に、biasがかかったfastqが出力される。
dbSNPsのSNPs。
上からオリジナルゲノム、germline、somatic、subclone。dbSNPsのSNPsは全ての細胞に発生している。
somatic deletion。
上からオリジナルゲノム、germline、somatic、subclone。somaticのSVはsomatic、subcloneだけに発生している。
somatic subclone deletion。
上からオリジナルゲノム、germline、somatic、subclone。subcloneのSVはsubcloneだけに発生している（ポジションは、somaticのfastaを元に設定されているので、元のゲノムからジャンプするとずれがある）。
マニュアルPDFに、germline、somatic、tumorのsubclone変異を発生させる例が載っています。
引用
Pysim-sv: a package for simulating structural variation data with GC-biases
Xia Y, Liu Y, Deng M, Xi R
BMC Bioinformatics. 2017 Mar 14;18(Suppl 3):53
