变异注释工具SnpEff详解

介绍如何使用SnpEff程序做变异注释，并衔接maftools包 。

SnpEff的优点：

1. 基于Java环境，便于安装和运行(More tiny than VEP)；2. 注释结果丰富(More HGVS than ANNOVAR)；3. 可自定义数据库。

HGVS(链接)是国际通用的变异命名标准(HGVS符号)，这里指代符合用其标准书写的变异位点。

安装和使用：

# SnpEff程序安装：
mkdir -p  ~/wgs/public
cd  ~/wgs/public
# 下载、解压后即可使用（提前具备java运行环境）
wget -c https://sourceforge.net/projects/snpeff/files/snpEff_latest_core.zip/
unzip snpEff_latest_core.zip
java -jar  ~/wgs/public/snpEff/snpEff.jar  download GRCh38.86
# 保持网络畅通，耐心等待 ...

# 使用方法：
java -jar snpEff.jar hg38 yourSampleName.normalize.vcf \
  -hgvs1LetterAa -canon > yourSampleName.snpeff.vcf

参数

-i

Input format [ vcf, bed ]. Default: VCF.

-csvStats

创建CSV摘要文件。

-hgvs

氨基酸子字段使用HGVS标注。 默认值: true

-lof

添加功能丧失(loss of function, LOF)和无义变异介导的衰变(Nonsense mediated decay, NMD) 标签。 

-canon

Only use canonical (典型/典范/标准) transcripts. 即使用典型的转录本。

-noshifthgvs

不根据HGVS表示法，来移动变异(most 3prime end, 即3'原则)。Do not shift variants according to HGVS notation (most 3prime end).

-interval <file>

在TXT/BED/BigBed/VCF/GFF文件中使用自定义间隔(可以多次使用这个选项)  

-hgvs1LetterAa

在HGVS符号中使用一个字母的氨基酸代码。默认值: FASLE （即默认不使用）。

密码子表、氨基酸单字母及单字母对照表。三字母更易读懂，单字母更简洁

其它参数：

常用选项

-chr <string>: 在染色体名称前添加'string'(例如: 'chr1'而不是'1')。 只在TXT输出。 

-classic: 使用旧风格的注释代替Sequence Ontology and Hgvs。 

-csvStats <file>: 创建CSV摘要文件。 

-download: 如果没有参考基因组，则下载。 默认值: true  

-fileList: 输入实际上要处理的文件的列表。 

-s, -stats, -htmlStats: 创建HTML摘要文件。 默认是“snpEff_summary.html”  

-noStats: 不创建stats(摘要)文件

-t: 使用多个线程(意味着'-noStats')。 默认的“off”  

结果过滤

-no-downstream: 不显示下游变化

-no-intergenic: 不显示基因间变化  

-no-intron: 不显示内含子变化  

-no-upstream: 不显示上游的变化

-no-utr: 不显示5_PRIME_UTR或3_PRIME_UTR的更改  

-no <effectType>: 不显示“effectType”。 这个选项可以多次使用。  

-fi， -filterInterval <file>: 只分析与此文件中指定的间隔(intervals)相交的变化(可以多次使用此选项)  Only analyze changes that intersect with the intervals specified in this file (you may use this option many times)

注释选项

-cancer: 进行“癌症”比较 (Somatic vs Germline)。 默认值: false

-cancerSamples <file>: 两列TXT文件，定义“original \t derived”样本。 

-formatEff: 使用与旧版本兼容的'EFF'字段(而不是'ANN')。 

-geneId: 使用基因ID代替基因名(VCF输出)。 默认值: false

-hgvsOld: 使用旧的HGVS符号。 默认值: false

-hgvsTrId: 在HGVS表示法中使用transcript ID。 默认值: false

-noHgvs: 不添加HGVS注释。 

-noLof: 不添加LOF和NMD注释。 

-oicr: 在VCF文件中添加OICR tag。 默认值: false

-sequenceOntology: 使用序列本体术语。 默认值: true

通用选项

-c, -config: 指定配置文件  

-configOption name=value: 覆盖配置文件选项  

-d， -debug: Debug mode (very verbose). 

-dataDir <path>: 覆盖配置文件中的data_dir参数。  

-nodownload: 如果本地没有SnpEff数据库，不要下载。 

-h， -help: 显示帮助并退出  

-noLog: 不向服务器报告使用统计信息  

-q, -quiet: 安静模式(不显示任何消息或错误)  

-v, -verbose: 冗长提示模式  

-version: 显示版本号并退出  

数据库选项

-canonList <file>: 只使用典型的(canonical)转录本；指定一个文件<file>，并使用其中的'gene_id transcript_id'条目替换一些转录本。 

-interaction: 使用交互进行注释(需要交互数据库)。 默认值: true

-maxTSL <TSL_number>: 只使用转录本支持性级别低于<TSL_number>。Only use transcripts having Transcript Support Level lower than <TSL_number>.

-motif: 使用Motif进行注释(需要Motif数据库)。 默认值: true

-nextProt: 使用NextProt注释(需要NextProt数据库)。 

-noGenome: 不加载任何基因组数据库(例如，使用自定义文件进行注释)。 

-noExpandIUB: 在输入变量中禁用IUB码扩展  

-noInteraction: 禁用交互注释  

-noMotif: 禁用主题注释。 

-noNextProt: 禁用NextProt注释。 

-onlyReg: 只使用调控tracks。Only use regulation tracks.

-onlyProtein: 只使用蛋白质编码转录本。 默认值: false

-onlyTr <file.txt>: 只使用此文件中的转录本。 格式:每行一个笔录ID。 

-reg <name>: 要使用的调控track(Regulation track)(该选项可以使用多次添加)。 

-ss, -spliceSiteSize <int>: 设置碱基中剪接位点(donor和acceptor)的大小。 默认值: 2  

-spliceRegionExonSize <int>: 设置外显子内，剪接位点区域大小。 默认值: 3 bases

-spliceRegionIntronMin <int>: 设置内含子内，剪接位点区域的最小碱基数。 默认值: 3 bases  

-spliceRegionIntronMax <int>: 设置内含子内，剪接位点区域的最大碱基数。 默认值: 8 bases 

-strict: 只使用“验证过的”转录本(即序列已被检查过)。 默认值: false

-ud, -upDownStreamLen <int>: 设置上游下游间隔(interval)长度(以bases为单位) 

SnpEff变异注释术语

查看SnpEff的注释结果汇总，可打开这个文件：

    snpEff_summary.html

结果解读参考资料：

    http://snpeff.sourceforge.net/VCFannotationformat_v1.0.pdf

其它技术文档：

    https://pcingola.github.io/SnpEff/se_inputoutput/

    http://snpeff.sourceforge.net/VCFannotationformat_v1.0.pdf

SnpEff注释结果中，可能对蛋白有害的术语包括：

missense, exon_loss, rare_amino_acid_variant, 

frameshift, inframe_deletion, inframe_insertion, disruptive_inframe_insertion, disruptive_inframe_deletion, 

splice_acceptor_variant, splice_donor_variant, splice_region_variant, 

start_lost, start_retained, initiator_codon_variant, 

stop_gained, stop_lost, stop_retained, 

protein_protein_contact, structural_interaction_variant,   

duplication, inversion, 

feature_ablation, gene_fusion, bidirectional_gene_fusion, 

conserved_intergenic_variant, conserved_intron_variant

另外，相比CADD或VEP，SnpEff的注释结果包含结构互作 (Structural interaction) 变异。

---

Additionally, we check whether SnpEff annotated the variant as protein_protein_contact or structural_interaction_variant.

Notes: These annotations are not present in the manually built SnpEff databases.

The score from this part (functional_part_score) is:

• 1.0 for protein changing variants localized within the functional region without benign variation;
• 1.0 for protein changing variants predicted by SnpEff to change protein residue essential for structural conformation or protein-protein interactions (this annotation is PDB based);
• 0.0 for other variants (different mutation effect, not in the functional region, or within the functional region with known benign variants).

来源：

https://intelliseq.com/variants-classification-according-to-the-acmg-amp-criteria/

根据ACMG和AMP标准的变异分类

---

SnpEff变异注释结果与maftools的衔接工具

工具的链接：

https://github.com/tsy19900929/snpeffToMaf

或工具的百度云链接：https://pan.baidu.com/s/15bypEB2QwdO6zaJwT20Erg

提取码：58ff

代码：

############################ in shell
  perl  $public/bin/snpeffToMaf-master/snpeffToMaf.pl \
  ${result}/Genotype.cohort.dbSNP.anno.g.vcf 20 0.1
  # Germline: 20 for minimum depth, 0.1 for minimum allele frequency
  # Tumour: 100 for minimum depth, 0.05 for minimum allele frequency
  # warn: you may need edit code once if vcf created by other variant callers. please issue me

# maf格式(突变注释格式)结果：
ls *.maf

cat *.maf | awk '!/Hugo_Symbol/ || NR==1' > all.maf
# concatenate all of yourSampleName.maf, but not necessary

############################ in R
library(maftools)
syn <- c("synonymous_variant","start_retained","stop_retained_variant")
df <- data.table::fread("all.maf")
vc <- names(table(df$Variant_Classification))
nonSyn <- setdiff(vc,syn)
colors <- rainbow(length(nonSyn))
names(colors) <- nonSyn
maf <- read.maf("all.maf", vc_nonSyn = nonSyn)
plotmafSummary(maf, rmOutlier = TRUE, addStat = 'median', dashboard = TRUE, titvRaw = FALSE, color = colors)

 其它变异注释方法

暂不过多介绍，例如使用annovar注释

convert2annovar.pl -format xx.vcf > xx.annovar
annotate_variation.pl -buildver hg38 --geneanno --outfile xx.anno xx.annovar   ANNOVAR/annovar/humandb/
annotate_variation.pl -buildver hg38 --dbtype knownGene --geneanno --outfile xx.anno xx.annovar  ANNOVAR/annovar/humandb/

概念与文献

数据库

分析技术