0.关于snpEff,建议大家去官网看一下作者。逝者已去,深切缅怀。同时看下官方文档
官方文档有几处语言或者程序命令错误
-
教程地址
下载snpEff地址
解压unzip snpEff_latest_core.zip
我的路径是/home/chaim/bsa/snpEff/ - 配置玉米zm437版本的数据库
在路径/home/chaim/bsa/snpEff/snpEff/
目录下创建文件夹data
,
cd /home/chaim/bsa/snpEff/snpEff/
mkdir data
cd data
mkdir genomes
mkdir zm437
#在zm437目录存放基因组注释文件genes.gff, 蛋白库,protein.fa
#在genomes目录放置基因组参考序列 zm437.fa
注意上述的基因组注释文件是gff3格式。
修改snpEff.config的参数
添加如下内容
#maize genome,version zm437
zm437.genome:maize
回到snpEFF目录,运行命令
java -jar snpEff.jar build -gff3 -v zm437
- 对vcf格式文件进行注释:
bwa
目录存放着GATK4处理之后的文件common_filtration.vcf
在/home/chaim/bsa/bwa
目录执行下面命令
java -Xmx8g -jar /home/chaim/bsa/snpEff/
snpEff/snpEff.jar zm437 common_filtration.vcf > common.eff
会输出三个文件,
snpEff_genes.txt
snpEff_summary.html
common.eff.vcf
如果不是在安装目录运行snpEff,需要加上-c参数。
参数说明:
-c 安装snpEff的目录config文件的位置
-verbose 详细模式,方便查看出错信息
zm442 我自己生成的基因组注释信息
命令如下:
java -Xmx8g -jar /public/home/soft/snpEff/snpEff/snpEff.jar -c /public/home/soft/snpEff/snpEff/snpEff.config zm442 F8.Filter.vcf -verbose -stats F8.summary.html -csvStats F8.csv >F8.ann.vcf
- 如果想更改使用其他注释文件,
删除/home/chaim/bsa/snpEff/snpEff//data/zm437/snpEffectPredictor.bin
该文件删除即可。
重新从步骤1开始即可。
5.产出结果分析
#默认是ann
java -Xmx8g -jar /home/chaim/bsa/snpEff/
snpEff/snpEff.jar ann zm442 common_filtration.vcf > common.eff
*.ann.vcf 是一个注释结果文件,其就在vcf的info信息新添加了anno一列信息,其具体每个值含义如下:
Allele
突变之后的碱基,第一个突变位点由T碱基突变成了C碱基,对应Allel的值为C
Annotation
由sequence ontology定义的突变类型
Annotation_Impact
对变异位点有害程度的简单评估,取值有HIGH, MODERATE, LOW, MODIFIER 4种,含义如下
Impact | Meaning | Example |
---|---|---|
HIGH | The variant is assumed to have high (disruptive) impact in the protein, probably causing protein truncation, loss of function or triggering nonsense mediated decay. | stop_gained, frameshift_variant |
MODERATE | A non-disruptive variant that might change protein effectiveness. | missense_variant, inframe_deletion |
LOW | Assumed to be mostly harmless or unlikely to change protein behavior. | synonymous_variant |
MODIFIER | Usually non-coding variants or variants affecting non-coding genes, where predictions are difficult or there is no evidence of impact. | exon_variant, downstream_gene_variant |
Snp下游的分析
可以使用snpEff注释的vcf进行4DTv位点分析,然后用其构建进化树。
或者是直接使用vcf构建进化树。
两种方法构建进化树均已经实现流程自动化。Vcf2Tree github
注释之后,可以计算位点的Ka,Ks值
根据https://www.biostars.org/p/304091/的方法,提取snpeff输出的结果,从中获取Ka,Ks
sample="样本名称"
grep '^#\|missense_variant\|synonymous_variant' ${sample}.ann.vcf >${sample}.mis_syn.txt
python3 KaKs.py ${sample}.mis_syn.txt >${sample}.KaKs.txt
kaks.py
从https://github.com/MerrimanLab/selectionTools/blob/master/extrascripts/kaks.py下载,我把它修改为支持python3的模式存放在https://github.com/chaimol/bio_code/blob/master/pipline/Vcf2Tree/KaKs.py下载。
关于KaKs的值,建议是过滤掉>3的Ks的值。
更多的关于过滤Ks和Ka/Ks的标准的讨论,参考文献https://onlinelibrary.wiley.com/doi/full/10.1111/mec.15275