snpeff-variant-annotation

SnpEff + SnpSift — Variant Annotation and Filtering

Overview

SnpEff annotates variants in VCF files by predicting their functional consequences: impact level (HIGH, MODERATE, LOW, MODIFIER), affected gene and transcript, amino acid change, and HGVS notation. SnpSift is the companion tool for filtering, sorting, and enriching annotated VCFs with external databases such as ClinVar and dbSNP. Together they form a fast, self-contained pipeline for going from raw variant calls to biologically interpretable, filtered variant sets. Both tools are Java-based and are invoked from the command line or Python subprocess; pre-built genome databases (hg38, GRCh37, mm10, and 100+ others) are downloaded with a single command.

When to Use

• Annotating VCF files from GATK, DeepVariant, bcftools, or other callers with predicted gene-level functional consequences before manual review or downstream filtering
• Prioritizing clinically relevant variants by filtering to HIGH-impact stop-gain, frameshift, and splice-site variants for rare disease or cancer gene panel analysis
• Adding ClinVar pathogenicity classifications and dbSNP rsIDs to a variant set for cross-study comparison or clinical reporting
• Extracting structured, tab-delimited fields (gene, protein change, AF, ClinSig) from annotated VCFs into pandas DataFrames for statistical analysis
• Identifying candidate de novo variants in trio analysis by combining allele frequency thresholds, impact filters, and parent VCF exclusion
• Use ANNOVAR instead when comprehensive annotation from multiple databases (gnomAD, CADD, SpliceAI) in a single run is required
• Use Ensembl VEP instead when REST API access or VEP-specific plugins (CADD, LOFTEE, SpliceRegion) are needed

Prerequisites

• Java: Java 11+ (required for SnpEff 5.x)
• SnpEff JAR: downloaded from the SnpEff releases page or via conda/bioconda
• Python packages (optional): cyvcf2, pandas, matplotlib, seaborn for Python-side parsing and visualization
• Reference genome database: downloaded once per assembly (e.g., hg38, GRCh37, mm10)

Check before installing: The tool may already be available in the current environment (e.g., inside a pixi / conda env). Run command -v snpEff first and skip the install commands below if it returns a path. When running inside a pixi project, invoke the tool via pixi run snpEff rather than bare snpEff.

# Download SnpEff JAR
wget https://snpeff.blob.core.windows.net/versions/snpEff_latest_core.zip
unzip snpEff_latest_core.zip
# JAR is at snpEff/snpEff.jar and snpEff/SnpSift.jar

# Or via conda (recommended for reproducibility)
conda install -c bioconda snpeff

# Verify
java -jar snpEff/snpEff.jar -version
# SnpEff 5.2a (build 2024-02-06)

# Install Python packages for downstream parsing
pip install cyvcf2 pandas matplotlib seaborn

Quick Start

# 1. Download the hg38 genome database (one-time setup, ~1 GB)
java -jar snpEff/snpEff.jar download hg38

# 2. Annotate variants
java -jar snpEff/snpEff.jar \
    -v hg38 \
    input.vcf.gz \
    > annotated.vcf

# 3. Filter to HIGH-impact variants
java -jar snpEff/SnpSift.jar filter \
    "ANN[*].IMPACT = 'HIGH'" \
    annotated.vcf \
    > high_impact.vcf

echo "Done. Check snpEff_summary.html and snpEff_genes.txt for QC stats."

Workflow

Step 1: Install SnpEff and Download Genome Database

Download and verify the pre-built genome database for your target assembly. Databases include transcript models from Ensembl or UCSC and are built into SnpEff's local cache directory (~/.snpEff/data/ or snpEff/data/).

# List all available databases (grep for your assembly)
java -jar snpEff/snpEff.jar databases | grep -i "GRCh38\|hg38"

# Download hg38 (Homo sapiens, GRCh38)
java -jar snpEff/snpEff.jar download hg38

# Download GRCh37 (older assemblies still widely used in clinical pipelines)
java -jar snpEff/snpEff.jar download GRCh37.75

# Download mouse mm10
java -jar snpEff/snpEff.jar download mm10

# List installed databases
ls ~/.snpEff/data/

Step 2: Annotate VCF with Functional Effects

Run SnpEff annotation to add ANN INFO fields to each variant. The ANN field encodes pipe-separated annotations per transcript: Allele|Effect|Impact|Gene|GeneID|Feature|FeatureID|BioType|Rank|HGVS.c|HGVS.p|cDNA_pos|CDS_pos|Protein_pos|Distance|Errors.

# Annotate a gzipped VCF (outputs to stdout, pipe or redirect)
java -Xmx8g -jar snpEff/snpEff.jar \
    -v \
    -stats snpeff_summary.html \
    hg38 \
    input.vcf.gz \
    > annotated.vcf

# Compress and index the output for downstream tools
bgzip annotated.vcf
tabix -p vcf annotated.vcf.gz

echo "Annotated variants in annotated.vcf.gz"
echo "QC report: snpeff_summary.html"
echo "Gene table: snpEff_genes.txt"

Step 3: Filter HIGH-Impact Variants with SnpSift

SnpSift filter evaluates boolean expressions over INFO and FORMAT fields. The ANN[*] syntax iterates over all transcript annotations for each variant — any matching transcript qualifies the variant.

# Filter for HIGH-impact variants only (stop-gain, frameshift, splice-site)
java -jar snpEff/SnpSift.jar filter \
    "ANN[*].IMPACT = 'HIGH'" \
    annotated.vcf.gz \
    > high_impact.vcf

# Filter HIGH or MODERATE impact and allele frequency < 1%
# (requires AF in INFO field, e.g., from GATK or gnomAD annotation)
java -jar snpEff/SnpSift.jar filter \
    "(ANN[*].IMPACT = 'HIGH' | ANN[*].IMPACT = 'MODERATE') & (AF < 0.01)" \
    annotated.vcf.gz \
    > rare_functional.vcf

# Count filtered variants
grep -v "^#" high_impact.vcf | wc -l

Step 4: Add ClinVar and dbSNP Annotations

SnpSift annotate transfers INFO fields from a reference VCF (ClinVar, dbSNP) to the target VCF by matching on chromosome and position. Download ClinVar and dbSNP VCFs from NCBI FTP before running.

# Download reference databases (one-time setup)
# ClinVar (GRCh38)
wget https://ftp.ncbi.nlm.nih.gov/pub/clinvar/vcf_GRCh38/clinvar.vcf.gz
wget https://ftp.ncbi.nlm.nih.gov/pub/clinvar/vcf_GRCh38/clinvar.vcf.gz.tbi

# dbSNP (GRCh38, b156 or later)
wget https://ftp.ncbi.nlm.nih.gov/snp/latest_release/VCF/GCF_000001405.40.gz
wget https://ftp.ncbi.nlm.nih.gov/snp/latest_release/VCF/GCF_000001405.40.gz.tbi

# Annotate with ClinVar CLNSIG and CLNDN fields
java -jar snpEff/SnpSift.jar annotate \
    -info CLNSIG,CLNDN,CLNREVSTAT \
    clinvar.vcf.gz \
    annotated.vcf.gz \
    > annotated_clinvar.vcf

# Annotate with dbSNP rsIDs (adds RS field to INFO)
java -jar snpEff/SnpSift.jar annotate \
    GCF_000001405.40.gz \
    annotated_clinvar.vcf \
    > annotated_full.vcf

echo "ClinVar + dbSNP annotations added to annotated_full.vcf"

Step 5: Extract Fields to Tab-Delimited Output

SnpSift extractFields converts VCF annotations into a flat table. Use ANN[0] to take the first (most severe) transcript annotation, or ANN[*] to expand all transcripts.

# Extract key fields: CHROM, POS, REF, ALT, impact, gene, protein change, AF, ClinVar sig
java -jar snpEff/SnpSift.jar extractFields \
    -s "," \
    -e "." \
    annotated_full.vcf \
    CHROM POS REF ALT \
    "ANN[0].GENE" \
    "ANN[0].EFFECT" \
    "ANN[0].IMPACT" \
    "ANN[0].HGVS_P" \
    "ANN[0].HGVS_C" \
    "ANN[0].FEATUREID" \
    AF \
    CLNSIG \
    CLNDN \
    > variants_table.tsv

# Preview
head -3 variants_table.tsv
# CHROM  POS   REF  ALT  ANN[0].GENE  ANN[0].EFFECT            ANN[0].IMPACT  ...
# chr1   69511 A    G    OR4F5        synonymous_variant        LOW            ...
# chr1   925952 G   A    SAMD11       missense_variant          MODERATE       ...

Step 6: Parse Annotated VCF in Python with cyvcf2

Use cyvcf2 (or PyVCF2) for programmatic VCF parsing in Python. cyvcf2 wraps htslib and is significantly faster than pure-Python parsers.

import re
import pandas as pd
from cyvcf2 import VCF

records = []

for variant in VCF("annotated_full.vcf.gz"):
    ann