macでインフォマティクス

macでインフォマティクス

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

ロングリードのマッピングからSVを検出する Sniffles

structural variations (SV) Pacbio Nanopore long read

2019 3/10 リンク追記

2019 7/11 インストール修正. -o pacbioは正しくは-x pacbioです。修正しました。

2020 7/13  ツイート追記、構成修正

2023/07/13 ツイート追記

 

 

　SnifflesはロングリードのSV caller。Githubの説明によれば、主にPacBioのリード用に設計されているが、Oxford Nanoporeのリードにも使用できるとされる。ターゲット SVは、ゲノム上の構造変化（例えば、欠失、重複、挿入、逆位および転座）である。 Snifflesは、これらのタイプのすべてを検出することができ、nested SV（e.g. inversion flanked by deletions or an inverted duplication）を検出することができるとされる。 ロングリードのアライナーNGMLRと一緒にPreprintが公開されており、現在論文準備中と思われる。

 

 

wiki

https://github.com/fritzsedlazeck/Sniffles/wiki

 

 

 

 

インストール

cent OSに導入した。

本体 Github

https://github.com/fritzsedlazeck/Sniffles

wget https://github.com/fritzsedlazeck/Sniffles/archive/master.tar.gz -O Sniffles.tar.gz 
tar xzvf Sniffles.tar.gz 
cd Sniffles-master/ 
mkdir -p build/ 
cd build/ 
cmake .. 
make 

cd ../bin/sniffles* 
./sniffles

#bioconda(link)
mamba install -c bioconda sniffles
mamba install -c bioconda ngmlr

$ ./sniffles -h

 

Usage: sniffles [options] -m <sorted.bam> -v <output.vcf> 

 

Input/Output:

    -m <string>,  --mapped_reads <string>

        (required)  Sorted bam File

    -v <string>,  --vcf <string>

        VCF output file name

    -b <string>,  --bedpe <string>

         bedpe output file name

    --Ivcf <string>

        Input VCF file name. Enable force calling

    --tmp_file <string>

        path to temporary file otherwise Sniffles will use the current directory.

 

General:

    -s <int>,  --min_support <int>

        Minimum number of reads that support a SV. [10]

    --max_num_splits <int>

        Maximum number of splits per read to be still taken into account. [7]

    -d <int>,  --max_distance <int>

        Maximum distance to group SV together. [1000]

    -t <int>,  --threads <int>

        Number of threads to use. [3]

    -l <int>,  --min_length <int>

        Minimum length of SV to be reported. [30]

    -q <int>,  --minmapping_qual <int>

        Minimum Mapping Quality. [20]

    -n <int>,  --num_reads_report <int>

        Report up to N reads that support the SV in the vcf file. -1: report all. [0]

    -r <int>,  --min_seq_size <int>

        Discard read if non of its segment is larger then this. [2000]

    -z <int>,  --min_zmw <int>

        Discard SV that are not supported by at least x zmws. This applies only for PacBio recognizable reads. [0]

    --cs_string

        Enables the scan of CS string instead of Cigar and MD.  [false]

 

Clustering/phasing and genotyping:

    --genotype

        Enables Sniffles to compute the genotypes. [false]

    --cluster

        Enables Sniffles to phase SVs that occur on the same reads [false]

    --cluster_support <int>

        Minimum number of reads supporting clustering of SV. [1]

    -f <float>,  --allelefreq <float>

        Threshold on allele frequency (0-1).  [0]

    --min_homo_af <float>

        Threshold on allele frequency (0-1).  [0.8]

    --min_het_af <float>

        Threshold on allele frequency (0-1).  [0.3]

 

Advanced:

    --report_BND

        Report BND instead of Tra in vcf output.  [false]

    --report_seq

        Report sequences for indels in vcf output. (Beta version!)  [false]

    --ignore_sd

        Ignores the sd based filtering.  [false]

 

Parameter estimation:

    --skip_parameter_estimation

        Enables the scan if only very few reads are present.  [false]

    --del_ratio <float>

        Estimated ration of deletions per read (0-1).  [0.0458369]

    --ins_ratio <float>

        Estimated ratio of insertions per read (0-1).  [0.049379]

    --max_diff_per_window <int>

        Maximum differences per 100bp. [50]

    --max_dist_aln_events <int>

        Maximum distance between alignment (indel) events. [4]

 

 パスの通ったディレクトリに移動しておく。 macでビルドするには、cmakeのときにコンパイラのg++とgccのパスを与える必要がある（Github説明より）。

 

 

ラン

1、NGMLRでロングリードをマッピングする。 ONTのリードをマッピングするときは-x ont、pacbioのリードをマッピングするときは-x pacbioをつけ、パラメータを一括設定する。

#出力をsamtoolsに渡してsort.bamにする。
ngmlr -t 20 -r ref.fa -q input.fastq -x ont | samtools sort -@ 20 -O BAM -o mapped.sort.bam -
samtools index -@ 20 mapped.sort.bam

２、SnifflesでSVsをコール。

sniffles -s 20 -m mapped.sort.bam -v output.vcf
  • -m <string>  (required)  Sorted bam File 
  • -v <string>      VCF output file name []
  • -s <int>　  Minimum number of reads that support a SV. Default: 10

出力はVCFに準拠している。

#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  /home/uesaka/mapped.sort.bam

CP023201.1      0       0       N       <DUP>   .       PASS    PRECISE;SVMETHOD=Snifflesv1.0.9;CHR2=CP023201.1;END=4748374;STD_quant_start=1.200961;STD_quant_stop=1.640825;Kurtosis_quant_start=8.970050;Kurtosis_quant_stop=11.824786;SVTYPE=DUP;SUPTYPE=SR;SVLEN=4748374;STRANDS=-+;RE=312  GT:DR:DV        ./.:.:312

CP023201.1      15406   1       N       <INV>   .       PASS    IMPRECISE;SVMETHOD=Snifflesv1.0.9;CHR2=CP023201.1;END=15495;STD_quant_start=12.186058;STD_quant_stop=2.097618;Kurtosis_quant_start=-0.240186;Kurtosis_quant_stop=1.566116;SVTYPE=INV;SUPTYPE=SR;SVLEN=89;STRANDS=++;RE=10

       GT:DR:DV        ./.:.:10

CP023201.1      15397   2       N       <DEL>   .       PASS    PRECISE;SVMETHOD=Snifflesv1.0.9;CHR2=CP023201.1;END=24137;STD_quant_start=7.355270;STD_quant_stop=8.185353;Kurtosis_quant_start=4.369543;Kurtosis_quant_stop=3.125435;SVTYPE=DEL;SUPTYPE=SR;SVLEN=8740;STRANDS=+-;RE=16 GT:DR:DV        ./.:.:16

CP023201.1      15397   3       N       <DEL>   .       PASS   

 

 

Githubでは、SURVIVOR（検出されたSVを1000genomeのデータセットなどと比較して、SVにアノテーションをつける方法論）も使いわけ、複数サンプルの結果をマージして統合したVCFを出力する方法についても記載されています。

https://github.com/fritzsedlazeck/Sniffles/wiki/SV-calling-for-a-population

引用

Accurate detection of complex structural variations using single molecule sequencing

View ORCID ProfileFritz J Sedlazeck, Philipp Rescheneder, Moritz Smolka, Han Fang, Maria Nattestad, Arndt von Haeseler, Michael Schatz

doi: https://doi.org/10.1101/169557 bioRxiv preprint

https://www.biorxiv.org/content/early/2017/07/28/169557

 

関連