知名的拷贝数变异分析工具几乎都是为人类变异检测开发,对于动植物重测序分析有些尴尬。不过好在植物群体研究不必那么精细,用同样的工具也可做分析。
地址:https://github.com/abyzovlab/CNVnator
1.安装
建议直接用conda。
conda create -n cnv cnvnator
conda activate cnv
查看帮助:
$ cnvnator
Not enough parameters.
CNVnator v0.4.1
Usage:
cnvnator -root out.root [-genome name] [-chrom 1 2 ...] -tree file1.bam ... [-lite]
cnvnator -root out.root [-genome name] [-chrom 1 2 ...] -merge file1.root ...
cnvnator -root file.root [-genome name] [-chrom 1 2 ...] -vcf [file.vcf.gz | file.vcf] [-rmchr] [-addchr]
cnvnator -root file.root [-genome name] [-chrom 1 2 ...] -idvar [file.vcf.gz | file.vcf] [-rmchr] [-addchr]
cnvnator -root file.root [-genome name] [-chrom 1 2 ...] -mask strict.mask.file.fa.gz [-rmchr] [-addchr]
cnvnator -root file.root [-genome name] [-chrom 1 2 ...] [-d dir | -fasta file.fa.gz] -his bin_size
cnvnator -root file.root [-genome name] [-chrom 1 2 ...] -baf bin_size [-hap] [-useid] [-nomask]
cnvnator -root file.root [-chrom 1 2 ...] -stat bin_size
cnvnator -root file.root -eval bin_size
cnvnator -root file.root [-chrom 1 2 ...] -partition bin_size [-ngc]
cnvnator -root file.root [-chrom 1 2 ...] -call bin_size [-ngc]
cnvnator -root file.root -genotype bin_size [-ngc]
cnvnator -root file.root -view bin_size [-ngc]
cnvnator -pe file1.bam ... -qual val(20) -over val(0.8) [-f file]
cnvnator-root file.root [-chrom 1 2 ...] -cptrees newfile.root
cnvnator-root file.root -ls
Valid genomes (-genome option) are: NCBI36, hg18, GRCh37, hg19, mm9, hg38, GRCh38
2.测试
首先准备好数据,再利用一个小数据集,比如这里用一条染色体来测试一下流程。
准备基因组数据。需要将基因组按染色体/scaffold拆分成单条序列,放在一个目录下。
mkdir genome;cd genome
faSplit byname genome.fa genome
# faSplit可用conda安装,或者自己写脚本拆分
测试脚本,先用一条染色体试试:
#从bam文件中提取比对上的reads信息
cnvnator -root file.root -tree sample-1.rmdup.bam -chrom 1
#生成read depth分布图
cnvnator -root file.root -his 1000 -d genome/ -chrom 1
#计算统计结果
cnvnator -root file.root -stat 1000 -chrom 1
#RD信号分割
cnvnator -root file.root -partition 1000 -chrom 1
#拷贝数变异检测
cnvnator -root file.root -call 1000 -chrom 1 > cnv.call.txt
#转化为vcf,如果是conda安装,没有这个脚本,需要从GitHub上下载
/biosoft/CNVnator/cnvnator2VCF.pl cnv.call.txt genome >test.vcf
拷贝数cnv.all.txt结果:、
表头CNV_type coordinates CNV_size normalized_RD e-val1 e-val2 e-val3 e-val4 q0
- CNV_type有deletion和duplication两种类型;
- CNV_size 位于染色体区域;
- normalized_RD 矫正后的read depth;
- e-val1 t检验后的evalue值,通常该值越小,代表分析的结果越准确;
- q0 比对的质量值为0的reads占比,通常该值越大,代表分析的结果越不准确。
vcf结果:
没有自动加上样品名,需要自己修改。加上--prefix参数也只是改变CNV ID。
更多结果解读,请查看官方文档或网上其他教程。
3.动植物群体检测CNV
正式分析。
ls -l /project/gvcf/*.rmdup.bam |awk -F' ' '{print $8}' >sample.info
cat sample.info |while read id;do
sample=`basename $id |sed 's/.rmdup.bam//'`
echo $sample
cnvnator -root file.root -tree $id
cnvnator -root file.root -his 1000 -d genome/
cnvnator -root file.root -stat 1000
cnvnator -root file.root -partition 1000
cnvnator -root file.root -call 1000 > cnv.call.txt
/biosoft/CNVnator/cnvnator2VCF.pl cnv.call.txt genome >${sample}.cnv.vcf
sed -i "22s/cnv/${sample}/" ${sample}.cnv.vcf
bgzip ${sample}.cnv.vcf
tabix -p vcf ${sample}.cnv.vcf.gz
done
得到各个样本的拷贝数vcf文件,将它们合并成一个。
vcf-merge sample-1.cnv.vcf.gz sample-2.cnv.vcf.gz ...>merge.vcf
或对不同类型群体进行合并。
https://blog.csdn.net/yangl7/article/details/114656482
https://www.jianshu.com/p/98542359df20
https://blog.csdn.net/weixin_43569478/article/details/108079613