CNV calling | DELLY

基于paired-end mapping methodssplit-read analysis,检测结构变异。

写在前面

Delly is an integrated structural variant (SV) prediction method that can discover, genotype and visualize deletions, tandem duplications, inversions and translocations at single-nucleotide resolution in short-read massively parallel sequencing data. It uses paired-ends, split-reads and read-depth to sensitively and accurately delineate genomic rearrangements throughout the genome. Structural variants can be visualized using Delly-maze and Delly-suave.

参考资料

Manual: https://github.com/dellytools/delly
FAQ: https://groups.google.com/forum/#!topic/delly-users/fOL8J3Exod4

实战

版本:Delly v0.7.8 可以一次性call 所有类型的SVs;低版本的需要通过-t指定SV类型
00 输入文件

  • Bam files need to be sorted, indexed and ideally duplicate marked. If multiple libraries are present for a single sample these need to be merged in a single bam file with unique ReadGroup tags.
  • .excl: a file include information of telomere and centromere regions and also all unplaced contigs, those regions will be removed from calling.
  • 有群体时,建议将所有个体的.bcf文件合并后,再进行genotyping。
#install
git clone --recursive https://github.com/dellytools/delly.git
cd delly/
make all
#也可以直接下载编译后的二进制文件,直接使用,不需安装 https://github.com/dellytools/delly/releases/
#Running
DB="/root/cc/DB"
i="SXY"
samtools index $i.recal.bam
~/cc/biosoft/delly/delly_v0.7.8_linux_x86_64bit call -g $DB/hg19_chr.fa -x exc1.excl -o delly_out/$i.exc1.bcf $i.recal.bam

#bcftools view $i.target.bcf >$i.target.vcf 

#Merge SV sites into a unified site list
delly merge -o sites.bcf s1.bcf s2.bcf ... sN.bcf

#Genotype this merged SV site list across all samples. This can be run in parallel for each sample.
delly call -g $DB/target-ref.fa -v sites.bcf -o s1.geno.bcf  s1.bam

#Merge all genotyped samples to get a single VCF/BCF using bcftools merge
bcftools merge -m id -O b -o merged.bcf s1.geno.bcf s2.geno.bcf ... sN.geno.bcf

#Apply the germline SV filter which requires at least 20 unrelated samples
delly filter -f germline -o germline.bcf merged.bcf
bcftools view germline.bcf >germline.vcf

01 SV统计:长度分布、在各样本中的情况等
长度统计:DEL,DUP和INV,可以用INFO:END-POS获得;INS的长度值为 INFO:INSLEN,而[POS, INFO:END] 提供的大概的插入位置信息。
这里使用作者开发的一个小工具,统计size分布、每个样本的变异数目等。

git clone --recursive https://github.com/dellytools/svprops.git
cd svprops/
./svprops germline.vcf >germline.tab #size分布,各类型
./sampleprops germline.vcf >sample.tab #每个样本在多少SV中是RR,RA,AA(基因型统计)

#绘图
Rscript R/svprops.R germline.tab
Rscript R/sampleprops.R sample.tab

head -3germline.tab
#chr    start   chr2    end id  size    vac vaf singleton   missingrate svtype  ct  precise ci  inslen  homlen  ce  refgq   altgq   gqsum   rdratio medianrc    refratio    altratio    maxaltratio PEsupport   SRsupport   supportsum
#chr22  17900897    chr22   17900914    DEL00000000 18  7   0.152174    NA  0   DEL 3to5    1   6   0   5   1.93828 75  10000   41614   0.563525    632 0   0.721311    0.787879    0   228 737
#chr22  18003439    chr22   18004197    DEL00000006 759 9   0.195652    NA  0   DEL 3to5    0   174 0   0   0   15  45  558 0.658926    107 0   0.25    0.333333    11  0   123

head -3 sample.tab
#sample missing homref  het homalt
#sample1    0   76  37  14
#sample2    0   69  47  11

germline.tab各列信息如下

  • vac: variant allele count (across all samples)
  • vaf: variant allele frequency (across all samples)
  • singleton: N/A if present in multiple samples or sample name if only present in one
  • missingrate: how many samples have a missing genotype (useful after GQ calibration)
  • ct: Delly's INFO:CT,Paired-end signature induced connection type
  • precise: 1 if INFO:PRECISE,Precise structural variation
  • ci: Delly's INFO:CIPOS,PE confidence interval around POS
  • inslen: insertion length
  • homlen: homology length
  • ce: Delly's INFO:CE,Consensus sequence entropy
  • refgq: median reference genotype quality for all SV non-carriers (GT==0/0)
  • altgq: median alternative genotype quality for het. SV carriers
  • gqsum: total GQ sum (useful to flag repetitive sites that are poorly genotyped)
  • rdratio: read-depth ratio of SV carriers to non-carriers (useful for filtering CNVs)
  • medianrc: median coverage
  • refratio: median REF support for non-carriers
  • altratio: median ALT support for carriers
  • maxaltratio: max. ALT support for carriers
  • PEsupport: total paired end support across all samples
  • SRsupport: total split-read support across all samples
  • supportsum: total depth across all samples

sample.tab各列信息如下

  • missing:该个体有多少位点的genotype为missing
  • homref:该个体有多少位点为ref的纯合位点
  • het:该个体有多少位点为杂合位点
  • homalt:该个体有多少位点为alt的纯合位点


    sample.tab

Note

  • 样本数比较少时,利用filter可能使所有SV被过滤掉,作者建议根据INFO FIELD(PASS)来过滤SV。

  • exc1.excl:排除区域(端粒等区域,DELLY github上有参考文件)

chr1  #排除整条染色体
chr22   0   10000   telomere
chr22   13000000    16000000    centromere
chr22   51294566    51304566    telomere

Error

bcf文件转为vcf时,发现里面一片空白,可能是bcftools版本过低,建议更新后再转换格式即可。

bcftools view file.bcf >file.vcf
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
平台声明:文章内容(如有图片或视频亦包括在内)由作者上传并发布,文章内容仅代表作者本人观点,简书系信息发布平台,仅提供信息存储服务。

推荐阅读更多精彩内容