IMPORTANT ANNOUNCEMENT: Our testing has shown that not all combinations of snpEff/database versions produce high-quality results. Please see the Current Recommended Best Practices When Running SnpEff and Analysis of SnpEff Annotations Across Versions sections below to familiarize yourself with our recommended best practices BEFORE running snpEff.
Contents |
Introduction
Until recently we were using an in-house annotation tool for genomic annotation, but the burden of keeping the database current and our lack of ability to annotate indels has led us to employ the use of a third-party tool instead. After reviewing many external tools (including annoVar, VAT, and Oncotator), we decided that SnpEff best meets our needs as it accepts VCF files as input, can annotate a full exome callset (including indels) in seconds, and provides continually-updated transcript databases. We have implemented support in the GATK for parsing the output from the SnpEff tool and annotating VCFs with the information provided in it.SnpEff Setup and Usage
- Download the SnpEff core program. If you want to be able to run VariantAnnotator on the SnpEff output, you'll need to download a version of SnpEff that VariantAnnotator supports from this page (currently supported versions are listed below). If you just want the most recent version of SnpEff and don't plan to run VariantAnnotator on its output, you can get it from here.
- Unzip the core program
- Open the file snpEff.config in a text editor, and change the "database_repository" line to the following:
database_repository = http://sourceforge.net/projects/snpeff/files/databases/
- Download one or more databases using SnpEff's built-in download command:
java -jar snpEff.jar download GRCh37.64
- A list of available databases is here. The human genome databases have GRCh or hg in their names. You can also download the databases directly from the SnpEff website, if you prefer.
- The download command by default puts the databases into a subdirectory called data within the directory containing the SnpEff jar file. If you want the databases in a different directory, you'll need to edit the data_dir entry in the file snpEff.config to point to the correct directory.
- Run SnpEff on the file containing your variants, and redirect its output to a file. SnpEff supports many input file formats including VCF 4.1, BED, and SAM pileup. Full details and command-line options can be found on the SnpEff home page.
Supported SnpEff Versions
- If you want to take advantage of SnpEff integration in the GATK, you'll need to run SnpEff version 2.0.5 (note: the newer version 2.0.5d is currently unsupported by the GATK, as we haven't yet had a chance to test it)
Current Recommended Best Practices When Running SnpEff
These best practices are based on our analysis of various snpEff/database versions as described in detail in the Analysis of SnpEff Annotations Across Versions section below.- We recommend using only the GRCh37.64 database with SnpEff 2.0.5. The more recent GRCh37.65 database produces many false-positive Missense annotations due to a regression in the ENSEMBL Release 65 GTF file used to build the database. This regression has been acknowledged by ENSEMBL and is supposedly fixed as of 1-30-2012, however as we have not yet tested the fixed version of the database we continue to recommend using only GRCh37.64 for now.
- We recommend always running with "-onlyCoding true" with human databases (eg., the GRCh37.* databases). Setting "-onlyCoding false" causes snpEff to report all transcripts as if they were coding (even if they're not), which can lead to nonsensical results. The "-onlyCoding false" option should only be used with databases that lack protein coding information.
- Do not trust annotations from versions of snpEff prior to 2.0.4. Older versions of snpEff (such as 2.0.2) produced many incorrect annotations due to the presence of a certain number of nonsensical transcripts in the underlying ENSEMBL databases. Newer versions of snpEff filter out such transcripts.
Analysis of SnpEff Annotations Across Versions
- Analysis of the SNP annotations produced by snpEff across various snpEff/database versions: File:SnpEff snps comparison of available versions.pdf
- Both snpEff 2.0.2 + GRCh37.63 and snpEff 2.0.5 + GRCh37.65 produce an abnormally high Missense:Silent ratio, with elevated levels of Missense mutations across the entire spectrum of allele counts. They also have a relatively low (~70%) level of concordance with the 1000G Gencode annotations when it comes to Silent mutations. This suggests that these combinations of snpEff/database versions incorrectly annotate many Silent mutations as Missense.
- snpEff 2.0.4 RC3 + GRCh37.64 and snpEff 2.0.5 + GRCh37.64 produce a Missense:Silent ratio in line with expectations, and have a very high (~97%-99%) level of concordance with the 1000G Gencode annotations across all categories.
- Comparison of SNP annotations produced using the GRCh37.64 and GRCh37.65 databases with snpEff 2.0.5: File:SnpEff snps ensembl 64 vs 65.pdf
- The GRCh37.64 database gives good results provided you run snpEff with the "-onlyCoding true" option. The "-onlyCoding false" option causes snpEff to mark all transcripts as coding, and so produces many false-positive Missense annotations.
- The GRCh37.65 database gives results that are as poor as those you get with the "-onlyCoding false" option on the GRCh37.64 database. This is due to a regression in the ENSEMBL release 65 GTF file used to build snpEff's GRCh37.65 database. The regression has been acknowledged by ENSEMBL and is due to be fixed shortly.
- Analysis of the INDEL annotations produced by snpEff across snpEff/database versions: File:SnpEff indels.pdf
- snpEff's indel annotations are highly concordant with those of a high-quality set of genomic annotations from the 1000 Genomes project. This is true across all snpEff/database versions tested.
Example SnpEff Usage with a VCF Input File
Below is an example of how to run SnpEff version 2.0.5 with a VCF input file and have it write its output in VCF format as well. Notice that you need to explicitly specify the database you want to use (in this case, GRCh37.64). This database must be present in a directory of the same name within the data_dir as defined in snpEff.config.java -Xmx4G -jar snpEff.jar eff -v -onlyCoding true -i vcf -o vcf GRCh37.64 1000G.exomes.vcf > snpEff_output.vcfIn this mode, SnpEff aggregates all effects associated with each variant record together into a single INFO field annotation with the key EFF. The general format is:
EFF=Effect1(Information about Effect1),Effect2(Information about Effect2),etc.And here is the precise layout with all the subfields:
EFF=Effect1(Effect_Impact|Effect_Functional_Class|Codon_Change|Amino_Acid_Change|Gene_Name|Gene_BioType|Coding|Transcript_ID|Exon_ID),Effect2(etc...It's also possible to get SnpEff to output in a (non-VCF) text format with one Effect per line. See the SnpEff home page for full details.
Adding SnpEff Annotations using VariantAnnotator
Once you have a SnpEff output VCF file, you can use the VariantAnnotator walker to add SnpEff annotations based on that output to the input file you ran SnpEff on.There are two different options for doing this:
Option 1: Annotate with only the highest-impact effect for each variant
NOTE: This option works only with supported SnpEff versions. VariantAnnotator run as described below will refuse to parse SnpEff output files produced by other versions of the tool, or which lack a SnpEff version number in their header.The default behavior when you run VariantAnnotator on a SnpEff output file is to parse the complete set of effects resulting from the current variant, select the most biologically-significant effect, and add annotations for just that effect to the INFO field of the VCF record for the current variant. This is the mode we plan to use in our Production Data-Processing Pipeline.
When selecting the most biologically-significant effect associated with the current variant, VariantAnnotator does the following:
- Prioritizes the effects according to the categories (in order of decreasing precedence) "High-Impact", "Moderate-Impact", "Low-Impact", and "Modifier", and always selects one of the effects from the highest-priority category. For example, if there are three moderate-impact effects and two high-impact effects resulting from the current variant, the annotator will choose one of the high-impact effects and add annotations based on it. See below for a full list of the effects arranged by category.
- Within each category, ties are broken using the functional class of each effect (in order of precedence: NONSENSE, MISSENSE, SILENT, or NONE). For example, if there is both a NON_SYNONYMOUS_CODING (MODERATE-impact, MISSENSE) and a CODON_CHANGE (MODERATE-impact, NONE) effect associated with the current variant, the annotator will select the NON_SYNONYMOUS_CODING effect. This is to allow for more accurate counts of the total number of sites with NONSENSE/MISSENSE/SILENT mutations. See below for a description of the functional classes SnpEff associates with the various effects.
- Effects that are within a non-coding region are always considered lower-impact than effects that are within a coding region.
java -jar dist/GenomeAnalysisTK.jar \ -T VariantAnnotator \ -R /humgen/1kg/reference/human_g1k_v37.fasta \ -A SnpEff \ --variant 1000G.exomes.vcf \ (file to annotate) --snpEffFile snpEff_output.vcf \ (SnpEff VCF output file generated by running SnpEff on the file to annotate) -L 1000G.exomes.vcf \ -o out.vcfVariantAnnotator adds some or all of the following INFO field annotations to each variant record:
- SNPEFF_EFFECT - The highest-impact effect resulting from the current variant (or one of the highest-impact effects, if there is a tie)
- SNPEFF_IMPACT - Impact of the highest-impact effect resulting from the current variant (HIGH, MODERATE, LOW, or MODIFIER)
- SNPEFF_FUNCTIONAL_CLASS - Functional class of the highest-impact effect resulting from the current variant (NONE, SILENT, MISSENSE, or NONSENSE)
- SNPEFF_CODON_CHANGE - Old/New codon for the highest-impact effect resulting from the current variant
- SNPEFF_AMINO_ACID_CHANGE - Old/New amino acid for the highest-impact effect resulting from the current variant
- SNPEFF_GENE_NAME - Gene name for the highest-impact effect resulting from the current variant
- SNPEFF_GENE_BIOTYPE - Gene biotype for the highest-impact effect resulting from the current variant
- SNPEFF_TRANSCRIPT_ID - Transcript ID for the highest-impact effect resulting from the current variant
- SNPEFF_EXON_ID - Exon ID for the highest-impact effect resulting from the current variant
1 874779 . C T 279.94 . AC=1;AF=0.0032;AN=310;BaseQRankSum=-1.800;DP=3371;Dels=0.00;FS=0.000;HRun=0;HaplotypeScore=1.4493;InbreedingCoeff=-0.0045; MQ=54.49;MQ0=10;MQRankSum=0.982;QD=13.33;ReadPosRankSum=-0.060;SB=-120.09;SNPEFF_AMINO_ACID_CHANGE=G215;SNPEFF_CODON_CHANGE=ggC/ggT; SNPEFF_EFFECT=SYNONYMOUS_CODING;SNPEFF_EXON_ID=exon_1_874655_874840;SNPEFF_FUNCTIONAL_CLASS=SILENT;SNPEFF_GENE_BIOTYPE=protein_coding;SNPEFF_GENE_NAME=SAMD11; SNPEFF_IMPACT=LOW;SNPEFF_TRANSCRIPT_ID=ENST00000342066 1 874816 . C CT 2527.52 . AC=15;AF=0.0484;AN=310;BaseQRankSum=-11.876;DP=4718;FS=48.575;HRun=1;HaplotypeScore=91.9147;InbreedingCoeff=-0.0520; MQ=53.37;MQ0=6;MQRankSum=-1.388;QD=5.92;ReadPosRankSum=-1.932;SB=-741.06;SNPEFF_EFFECT=FRAME_SHIFT;SNPEFF_EXON_ID=exon_1_874655_874840; SNPEFF_FUNCTIONAL_CLASS=NONE;SNPEFF_GENE_BIOTYPE=protein_coding;SNPEFF_GENE_NAME=SAMD11;SNPEFF_IMPACT=HIGH;SNPEFF_TRANSCRIPT_ID=ENST00000342066
Option 2: Annotate with all effects for each variant
VariantAnnotator also has the ability to take the EFF field from the SnpEff VCF output file containing all the effects aggregated together and copy it verbatim into the VCF to annotate.Here's an example of how to do this:
java -jar dist/GenomeAnalysisTK.jar \ -T VariantAnnotator \ -R /humgen/1kg/reference/human_g1k_v37.fasta \ -E resource.EFF \ --variant 1000G.exomes.vcf \ (file to annotate) --resource snpEff_output.vcf \ (SnpEff VCF output file generated by running SnpEff on the file to annotate) -L 1000G.exomes.vcf \ -o out.vcfOf course, in this case you can also use the VCF output by SnpEff directly, but if you are using VariantAnnotator for other purposes anyway the above might be useful.
List of Genomic Effects
Below are the possible genomic effects recognized by SnpEff, grouped by biological impact. Full descriptions of each effect are available on this page.High-Impact Effects
- SPLICE_SITE_ACCEPTOR
- SPLICE_SITE_DONOR
- START_LOST
- EXON_DELETED
- FRAME_SHIFT
- STOP_GAINED
- STOP_LOST
Moderate-Impact Effects
- NON_SYNONYMOUS_CODING
- CODON_CHANGE (note: this effect is used by SnpEff only for MNPs, not SNPs)
- CODON_INSERTION
- CODON_CHANGE_PLUS_CODON_INSERTION
- CODON_DELETION
- CODON_CHANGE_PLUS_CODON_DELETION
- UTR_5_DELETED
- UTR_3_DELETED
Low-Impact Effects
- SYNONYMOUS_START
- NON_SYNONYMOUS_START
- START_GAINED
- SYNONYMOUS_CODING
- SYNONYMOUS_STOP
- NON_SYNONYMOUS_STOP
Modifiers
- NONE
- CHROMOSOME
- CUSTOM
- CDS
- GENE
- TRANSCRIPT
- EXON
- INTRON_CONSERVED
- UTR_5_PRIME
- UTR_3_PRIME
- DOWNSTREAM
- INTRAGENIC
- INTERGENIC
- INTERGENIC_CONSERVED
- UPSTREAM
- REGULATION
- INTRON
Functional Classes
SnpEff assigns a functional class to certain effects, in addition to an impact:- NONSENSE: assigned to point mutations that result in the creation of a new stop codon
- MISSENSE: assigned to point mutations that result in an amino acid change, but not a new stop codon
- SILENT: assigned to point mutations that result in a codon change, but not an amino acid change or new stop codon
- NONE: assigned to all effects that don't fall into any of the above categories (including all events larger than a point mutation)
---- note by WZY ----
this is an example script to run snpEff:
java -jar snpEff.jar download hg19
java -Xmx2G -jar snpEff.jar -v -chr chr21 -i vcf -o vcf hg19 input.vcf > output.vcf 
if you are not running the program from the directory where it is installed, you have to use:
ReplyDelete-config /home/zhengyuan/snpEff_3_1/snpEff.config