FunciSNP: an R/bioconductor tool integrating functional non-coding data sets with genetic association studies to identify candidate regulatory SNPs

Abstract

Single nucleotide polymorphisms (SNPs) are increasingly used to tag genetic loci associated with phenotypes such as risk of complex diseases. Technically, this is done genome-wide without prior restriction or knowledge of biological feasibility in scans referred to as genome-wide association studies (GWAS). Depending on the linkage disequilibrium (LD) structure at a particular locus, such tagSNPs may be surrogates for many thousands of other SNPs, and it is difficult to distinguish those that may play a functional role in the phenotype from those simply genetically linked. Because a large proportion of tagSNPs have been identified within non-coding regions of the genome, distinguishing functional from non-functional SNPs has been an even greater challenge. A strategy was recently proposed that prioritizes surrogate SNPs based on non-coding chromatin and epigenomic mapping techniques that have become feasible with the advent of massively parallel sequencing. Here, we introduce an R/Bioconductor software package that enables the identification of candidate functional SNPs by integrating information from tagSNP locations, lists of linked SNPs from the 1000 genomes project and locations of chromatin features which may have functional significance.

Availability: FunciSNP is available from Bioconductor (bioconductor.org).

Figure 1.

Summary of all known GWAS SNPs across a number of different cancer types. GWAS SNPs were annotated using known genomic features supplied by HOMER (version 3.9) (10), using build hg19 as reference. GWAS SNPs were extracted from the UCSC Genome table browser, track name ‘gwasCatalog’ with a P-value cut-off of 9e⁻⁰⁶ (1).

Figure 2.

Schematic flowchart to describe FunciSNP. Purple boxes represent process before integration with biofeature. Red boxes represent information after integration with biofeature.

Figure 3.

(A) Distribution of R² values of all YAFSNPs. Each marked bin contains the total number of YAFSNPs. The sum of all the counts would total the number of correlated SNPs. (B) Distribution of R² values of all YAFSNPs divided by the tagSNP and by its genomic location. (C) Histogram distribution of R² value for all 1kgSNP extracted and overlaps PolII. R² values are determined by its association to the tagSNP. (D) Scatter plot of the R² and distance to tagSNP for all 1kgSNP extracted and overlap PollI. (E) Stacked bar chart summarizing all correlated SNPs for each of the identified genomic features: exon, intron, 5UTR, 3UTR, promoter, lincRNA or in gene desert. R² cut-off at 0.5. This plot is most informative if used with a rsq value. (F) Heatmap of the number of 1kgSNPs by relationship between tagSNP and biofeature. Total number of YAFSNPs is listed within each quadrant to represent the number of potential candidate functional SNPs overlapping a biofeature (y-axis) which are in LD to the original tagSNP (x-axis).

Figure 4.

FunciSNP results viewed in UCSC genome browser. Tracks are ordered in the following manner: known GWAS hits, dbSNP135, FunciSNP result, biofeatures, refseq genes and known lincRNA. TagSNP is highlighted in the FunciSNP result track and each YAFSNP is color coded to reflect the number of biofeatures which it overlaps. The color ranges from blue (low number of biofeature overlap) to red (high number of overlap). Each YAFSNP is highlighted by its known rsID and the calculated R² value. The results are saved in a UCSC genome session: http://goo.gl/xrZPD.