Identification of polymorphic antioxidant response elements in the human genome

Abstract

Single nucleotide polymorphisms (SNPs) in transcription factor binding sites (TFBSs) may affect the binding of transcription factors, lead to differences in gene expression and phenotypes and therefore affect susceptibility to environmental exposure. We developed an integrated computational system for discovering functional SNPs in TFBSs in the human genome and predicting their impact on the expression of target genes. In this system, we (i) construct a position weight matrix (PWM) from a collection of experimentally discovered TFBSs; (ii) predict TFBSs in SNP sequences using the PWM and map SNPs to the upstream regions of genes; (iii) examine the evolutionary conservation of putative TFBSs by phylogenetic footprinting; (iv) prioritize candidate SNPs based on microarray expression profiles from tissues in which the transcription factor of interest is either deleted or over-expressed and (v) finally, analyze association of SNP genotypes with gene expression phenotypes. The application of our system has been tested to identify functional polymorphisms in the antioxidant response element (ARE), a cis-acting enhancer sequence found in the promoter region of many genes that encode antioxidant and Phase II detoxification enzymes/proteins. In response to oxidative stress, the transcription factor NRF2 (nuclear factor erythroid-derived 2-like 2) binds to AREs, mediating transcriptional activation of its responsive genes and modulating in vivo defense mechanisms against oxidative damage. Using our novel computational tools, we have identified a set of polymorphic AREs with functional evidence, showing the utility of our system to direct further experimental validation of genomic sequence variations that could be useful for identifying high-risk individuals.

Figure 1

The workflow for in silico identification of polymorphic antioxidant response elements (AREs). Cylinders represent databases used, and rectangles contain computational applications. First, a position weight matrix (PWM) model is constructed from functional AREs and used to search for novel AREs in the input SNP sequences; Next, SNPs selected by PWM scoring are mapped into the upstream regions of genes; Then, a phylogentic footpring approach is used to examine evolutionary conservation in order to reduce excess false-positives, and microarray expression profiles are used to prioritize candidate SNPs that map to oxidative stress induced genes; Last, the genotype-expression association is computed to test SNP's function in silico.

Figure 2

Modeling AREs (NRF2 binding sites) based on 57 functional ARE sequences discovered by experimental methods. These 21-base sequences were used as input to our system to obtain (2A) position frequency matrix, (2B) position weight matrix, used for direct scoring of DNA sequences, and (2C) sequence logo, a visualization showing contributions of individual positions and individual nucleotides to the overall ARE motif, generated by the WebLogo program (http://weblogo.berkeley.edu). In the logo, the letters are sorted so the most common one is on top, the height of each letter is made proportional to its frequency, and the height of the entire stack is adjusted to signify the ‘information content’ (measured in bits) of the sequences at that position. The logo displays both significant residues and subtle sequence patterns.

Figure 3

The ΔPWM (difference in PWM scores of an allele-pair) scores of 4967 putative ARE SNPs grouped by the position in ARE motif. The mean (circle and solid line) and the standard deviation (Y error bars) of ΔPWM scores at each position were plotted, as well as the number of ARE SNPs (square and dash line) at each position.

Figure 4

A phylogenetic footprinting approach used to detect conserved response elements in multiple species. (A) Survey of AREs in upstream 5-kb regions of all genes in human, mouse, and rat genomes using position weight matrix, and found 2388 orthologous genes have AREs in upstream; (B) Identification of conserved AREs by examining the multiple sequence alignment of upstream 5-kb regions of orthologous genes. An example of multiple sequence alignment shows a conserved, functional ARE, and three non-conserved putative AREs in thioredoxin (TXN) gene upstream 5-kb region.

Figure 5

An example of association of SNP genotype with gene expression phenotype. The first three tracks, obtained from the HapMap Project genome browser, display chromosomal region, human gene, and genotyped SNPs and their relative allelic frequencies (little pie charts). The bottom track shows the linear regression plots of gene expression level against SNP genotype. Among 4 association tests, the expression levels of hemoglobin epsilon 1 (HBE1) gene are associated with genotypes of the putative ARE SNP (rs7130110) in upstream (4.7 kb). The presence of the “C” allele correlated with higher basal levels of HBE1 gene expression.