Molecular mechanisms remain unknown for most type 2 diabetes genome-wide association study identified loci. Variants associated with type 2 diabetes and fasting glucose levels reside in introns of ADCY5, a gene that encodes adenylate cyclase 5. Adenylate cyclase 5 catalyzes the production of cyclic AMP, which is a second messenger molecule involved in cell signaling and pancreatic β-cell insulin secretion. We demonstrated that type 2 diabetes risk alleles are associated with decreased ADCY5 expression in human islets and examined candidate variants for regulatory function. rs11708067 overlaps a predicted enhancer region in pancreatic islets. The type 2 diabetes risk rs11708067-A allele showed fewer H3K27ac ChIP-seq reads in human islets, lower transcriptional activity in reporter assays in rodent β-cells (rat 832/13 and mouse MIN6), and increased nuclear protein binding compared with the rs11708067-G allele. Homozygous deletion of the orthologous enhancer region in 832/13 cells resulted in a 64% reduction in expression level of Adcy5, but not adjacent gene Sec22a, and a 39% reduction in insulin secretion. Together, these data suggest that rs11708067-A risk allele contributes to type 2 diabetes by disrupting an islet enhancer, which results in reduced ADCY5 expression and impaired insulin secretion.
Introduction
Genome-wide association studies (GWAS) have identified more than 80 loci associated with type 2 diabetes (1,2); however, the underlying biological and molecular mechanisms responsible for most of these associations remain unknown. One type 2 diabetes association signal overlaps introns 1–3 of the ADCY5 gene. The rs11717195-T and the rs11708067-A alleles were associated (P < 5 × 10−8) with type 2 diabetes risk in a GWAS meta-analysis of individuals of European ancestry (3), and rs11717195-T was associated (P = 2.2 × 10−8) with type 2 diabetes risk in a trans-ancestry GWAS meta-analysis (4). The rs11708067-A allele was associated (P = 9.1 × 10−4) with type 2 diabetes risk in a candidate variant study in individuals from South Asia (5) and in a candidate gene study (P = 4.7 × 10−3) in individuals of African American ancestry (6). rs11717195 and rs11708067 exhibit strong linkage disequilibrium (LD) with each other (r2 ≥0.938, 1000 Genomes phase 3) in individuals of European, East Asian, and South Asian ancestry (7,8). In individuals of sub-Saharan African ancestry, rs143882978 (LD r2 = 0.005 with rs11708067, 1000 Genomes phase 3 African) was also associated (P = 0.02) with type 2 diabetes (9).
Variants at the ADCY5 locus are also associated with other glycemic traits (Supplementary Table 1) that implicate them as variants affecting islet and β-cell function. rs11708067-A was associated (P = 1.7 × 10−14) with higher fasting glucose levels in individuals of European ancestry (10) and higher fasting glucose levels (P = 6.3 × 10−6) in individuals of African ancestry (11). An additional variant, rs2877716, in strong LD (r2 = 0.837, 1000 Genomes phase 3 European [EUR]) (7) with rs11708067, was associated (P = 4.2 × 10−16) with 2-h glucose level (12). Associations with fasting glucose, 2-h glucose, insulinogenic index, and 2-h insulin (rs9883204, n = 2,151, P ≤ 0.05) have also been reported in Asian Indians (13). ADCY5 variant alleles associated with higher glucose levels were also associated with lower birth weight (rs9883204, LD r2 = 0.826 with rs11708067, 1000 Genomes phase 3 EUR, P < 5 × 10−8) (14). The type 2 diabetes–, fasting glucose–, and 2-h glucose–associated variant, rs11708067, was also associated (P = 2.0 × 10−3) with impaired proinsulin-to-insulin conversion as individuals homozygous for the type 2 diabetes risk A allele showed higher proinsulin levels and a higher proinsulin/insulin ratio after an oral glucose tolerance test (15). Transethnic fine-mapping studies (11) and credible set analyses suggest rs11708067 is a functional variant at this locus.
ADCY5 encodes adenylate cyclase 5, which catalyzes generation of cyclic AMP (cAMP) from ATP (16). cAMP, ATP, and calcium (Ca2+) regulate insulin secretion in the pancreatic islet β-cell (17). Previous studies have established ADCY5 as a compelling candidate target gene influencing glucose metabolism, showing association of the type 2 diabetes risk alleles with decreased ADCY5 mRNA expression (18–20). Furthermore, ADCY5 knockdown in human islets reduced glucose-stimulated cAMP levels, inhibiting glucose metabolism toward ATP and glucose-stimulated insulin secretion (18).
Identification of regulatory variants that exhibit allelic differences at GWAS loci can help elucidate the molecular mechanism(s) responsible for the disease associations. Open chromatin data in human pancreatic islets, from DNase I hypersensitivity sequencing (DNase-seq) (21), formaldehyde-assisted isolation of regulatory elements sequencing (FAIRE-seq) (22), and assay for transposase-accessible chromatin with high-throughput sequencing (ATAC-seq) (20) experiments, can help to identify regulatory elements of genes that play a role in type 2 diabetes risk.
In this study, we aimed to identify the variant(s) and molecular mechanism(s) responsible for the type 2 diabetes and glycemic trait associations at the ADCY5 locus. We demonstrate that 1 of 18 candidate variants at the locus, rs11708067, exhibits allelic differences in reporter assays of transcriptional activity in two insulinoma cell lines, consistent with the direction of association between the type 2 diabetes risk alleles and expression of ADCY5. We present evidence that the rs11708067-A type 2 diabetes risk allele exhibits reduced enhancer potential (reduced representation in islet histone 3 lysine 27 acetylation [H3K27ac] chromatin immunoprecipitation sequencing [ChIP-seq] reads) in vivo. Additionally, rs11708067-A is bound by different transcription factor complexes in vitro. We also show that the orthologous enhancer region in 832/13 cells reduces Adcy5 expression and insulin secretion. Together, these data provide a candidate molecular mechanism for the regulatory function of genetic variation at this type 2 diabetes locus.
Research Design and Methods
Association With Human Islet Gene Expression and Additional Glycemic Traits
We evaluated evidence of association between type 2 diabetes and fasting glucose GWAS variants and gene expression levels in human islets by analyzing existing RNA-seq–based expression quantitative trait loci (eQTL) data from 112 individuals (20,23). We looked up the evidence of association between rs11708067 and gene expression level for genes with a transcription start site located within 1 Mb (Supplementary Fig. 1A and B). After rs11708067, which is the variant most strongly associated with ADCY5 expression level, no additional variant was significant (false discovery rate [FDR] <5%) in a stepwise forward selection approach, as previously described (20). The Genotype-Tissue Expression (GTEx) Project portal was accessed to look for variant association with ADCY5 expression (www.gtexportal.org).
We also evaluated associations between rs11708067 and 21 glucose- and insulin-related traits from ∼8,500 individuals without diabetes from the Metabolic Syndrome in Men (METSIM) study (Supplementary Table 2) (24,25). Trait values containing skewed distributions were log transformed, and we winsorized all traits at 5 SDs from the mean. We adjusted all trait measurements for BMI, age, and age squared and inverse-normalized the residuals prior to association analyses (26). We assumed an additive genetic model when testing for trait–variant associations and used a linear mixed model to correct for relatedness, as applied in Efficient Mixed-Model Association eXpedited (EMMAX) (27).
Identification of a Potential Transcriptional Enhancer Region Overlapping rs11708067
We used LDlink (8) to identify variants in strong LD (r2 >0.8) with rs11708067 in the 1000 Genomes Project phase 3 EUR (7) data set. We used FAIRE-seq (22), DNase-seq (21), chromatin states (28), and ATAC-seq data (20) from human pancreatic islets to identify potential transcriptional enhancer regions at ADCY5.
H3K27ac ChIP-seq in a Primary Human Islet Sample
We obtained human pancreatic islets through the Integrated Islet Distribution Program (IIDP) in accordance with regulations of human subject research. We performed H3K27ac ChIP-seq experiments and read mapping as previously described (21,28). We aligned sequence reads to the UCSC Genome Browser build hg19 considering variants from the 1000 Genomes Project using GSNAP (29) with AA-ALIGNER (30). When sequences were the same, we aligned the five sequences with the highest overall quality. We determined significance of allelic imbalance at rs11708067 using an exact binomial test with a null expectation of 0.5, based on the number of reads containing the reference allele and the total number of reads at the heterozygous site.
Cell Culture
MIN6 mouse insulinoma β-cells (31) were grown in DMEM (Sigma-Aldrich, St. Louis, MO) with 10% FBS, 1 mmol/L sodium pyruvate, and 0.1 mmol/L β-mercaptoethanol. INS-1–derived 832/13 rat insulinoma β-cells (a gift from C. Newgard, Duke University, Durham, NC) were grown in RPMI-1640 medium (Corning, NY) supplemented with 10% FBS, 10 mmol/L HEPES, 2 mmol/L L-glutamine, 1 mmol/L sodium pyruvate, and 0.05 mmol/L β-mercaptoethanol. We maintained cell lines at 37°C and 5% CO2. We seeded MIN6 or 832/13 cells (200,000 cells per well) in a 24-well plate 1 day before transfections.
Transcriptional Reporter Assays
We PCR-amplified 231-base pair (bp) segments (centered around rs11708067) from the DNA of individuals homozygous for each rs11708067 allele. Primers are listed in Supplementary Table 3. We cloned amplicons into the multiple cloning region of firefly luciferase reporter vector pGL4.23 (Promega, Fitchburg, WI) in both orientations with respect to the minimal promoter as described previously (32,33). We transfected three to five independent plasmids containing each allele in duplicate into MIN6 mouse or 832/13 rat insulinoma cells. For transfection of MIN6 cells, we added 250 ng of plasmid DNA to each well with Lipofectamine LTX (Thermo Fisher Scientific, Waltham, MA) and Opti-MEM (Gibco, Waltham, MA). For transfection of 832/13 cells, we added 770 ng of plasmid to each well using Lipofectamine 2000 (Thermo Fisher Scientific) and Opti-MEM (Gibco, Thermo Fisher Scientific). To control for transfection efficiency, we cotransfected phRL-TK Renilla luciferase reporter vector (Promega) with the pGL4.23 plasmids containing test sequences into MIN6 and 832/13 cells (80 ng each well). Cells were incubated at 37°C with 5% CO2 for 48 h. We performed luciferase assays with cell lysates using the Dual-Luciferase Reporter Assay System (Promega). We normalized firefly luciferase readings to Renilla luciferase readings and compared the readings to pGL4.23 minimal promoter empty vector control activity using unpaired two-sided t tests.
Electrophoretic Mobility Shift Assays
Complementary oligonucleotide probes (21 bp) centered on rs11708067 alleles (Supplementary Table 3) were synthesized by Integrated DNA Technologies, Inc. (Coralville, IA) and labeled with biotin on the 5′ end. We incubated single-stranded oligonucleotides in 10 mmol/L Tris, 50 mmol/L NaCl, 1 mmol/L EDTA, pH 8.0 at 95°C for 5 min followed by 70 1-min cycles with temperature decreasing by 1°C. We prepared nuclear lysates from MIN6 cells and 832/13 cells using Pierce NE-PER Nuclear and Cytoplasmic Extraction Reagents (Thermo Fisher Scientific). We performed electrophoretic mobility shift assay (EMSA) reactions using the LightShift Chemiluminescent EMSA Kit (Thermo Fisher Scientific). We incubated binding buffer and Poly(dI-dC) with 2–4 μg of 832/13 or MIN6 nuclear lysate at room temperature for 10 min, and then added biotin-labeled DNA probes. For competition reactions, we added 40- to 240-fold excess unlabeled double-stranded probe and incubated the reactions for 10 min before adding the biotin-labeled DNA. We then performed gel electrophoresis and transfer, wash, and detection steps as described previously (32,33). For supershift assays, we added 4–8 μg of the following antibodies to the nuclear lysate for 20 min before adding labeled DNA probes: CEBPA (sc-61×), CEBPB (sc-150× and sc-7962×), NKX2.2 (sc-15015×), FOXA2 (sc-9187×), PDX1 (sc-14662×), p300 (sc-585×), NKX6.1 (sc-15030×), MAFB (sc-22830×), YY1 (sc-281×), RXRA (sc-553×) antibodies (all from Santa Cruz Biotechnology, Inc., Dallas, TX).
CRISPR-Cas9–Mediated Deletion of Adcy5 Enhancer Element
We used BLAT (http://genome.ucsc.edu) to identify the human genome sequence of the enhancer region containing rs11708067 (hg19) in the rat genome (rn5) and designed guide RNA oligonucleotides to target the conserved enhancer region using CRISPR Design (http://crispr.mit.edu/). Three guide RNAs upstream of rs11708067 were cloned into a pSpCas9(BB)-2A-GFP (34) vector (PX458, a gift from F. Zhang [Broad Institute of MIT and Harvard, Cambridge, MA], Addgene plasmid #48138), and three guide RNAs downstream of rs11708067 were cloned into a modified PX458 vector containing an mCherry marker in place of GFP. Target DNA sequences for guide RNAs are listed in Supplementary Table 3 and were validated in each plasmid by Sanger sequencing. We transfected three combinations of guide RNA pair containing GFP and mCherry vectors into 832/13 cells using Lipofectamine 2000 (Thermo Fisher Scientific), generating three deletions of the enhancer region in the rat genome. We performed similar transfections using GFP and mCherry vectors lacking guide sequences to generate GFP/mCherry empty vector control clones. After incubating the 832/13 cells at 37°C with 5% CO2 for 72 h, we selected single 832/13 cells positive for both GFP and mCherry fluorescence by fluorescence-activated cell sorting performed by the University of North Carolina at Chapel Hill (UNC) Flow Cytometry Core Facility using the MoFlo XDP sorter. We plated single GFP- and mCherry-positive cells in each well and allowed the cells to grow for 21 days. To validate homozygous deletions of the enhancer region, we used QuickExtract (Epicentre, Madison, WI) to extract DNA from clones, performed diagnostic PCR with primers spanning the targeted genomic enhancer region (listed in Supplementary Table 3), and performed Sanger sequencing.
Evaluation of Adcy5 and Sec22a Gene Expression
We extracted RNA from CRISPR-Cas9–edited clones, GFP and mCherry empty vector control clones, or untransfected 832/13 cells using the RNeasy Plus Mini kit (Qiagen, Hilden, Germany) from two independent experiments. We plated 200,000 cells in 24-well plates, extracted RNA after 24 h, and synthesized cDNA using the SuperScript III First-Strand Synthesis System (Thermo Fisher Scientific). We measured gene expression of Adcy5, Sec22a, and B2m control by quantitative PCR using TaqMan assays (Rn00575059_m1, Rn00589273_m1, Rn00560865_m1; Applied Biosystems, Thermo Fisher Scientific) and the standard curve method. We generated standard curves using 1:3 serial dilutions of cDNA, measured in triplicate. We normalized Adcy5 and Sec22a gene expression level to the mean expression level of the control B2m gene. We plotted gene expression for 832/13 cells, homozygous deletion clones, and GFP/mCherry empty vector control clones as a percentage of the average gene expression measured for the control clones. We compared gene expression for intact control clones with enhancer homozygous deletion clones using unpaired t tests.
Insulin Secretion Assays
We plated 250,000 CRISPR-Cas9–edited, GFP/mCherry empty vector control, and untransfected 832/13 cells into 24-well plates. We measured insulin secretion at 3 mmol/L glucose and 18 mmol/L glucose concentrations 48 h after plating or when cells were 90% confluent. We performed insulin secretion assays on three separate occasions; on each occasion, we measured insulin secretion for each sample in duplicate or triplicate using Rat Insulin ELISA (Mercodia, Winston-Salem, NC). We compared insulin secretion for intact control clones with enhancer homozygous deletion clones using unpaired t tests.
Results
rs11708067 Is Associated With ADCY5 Expression Level in Primary Human Pancreatic Islets and With Glycemic Traits in the METSIM Study
We examined evidence of association between the type 2 diabetes risk allele rs11708067-A and expression level of nearby (<1 Mb) genes by evaluating islet eQTL data from 112 human pancreatic islet samples (20,23). rs11708067-A allele was associated with decreased ADCY5 mRNA level (P = 3.09 × 10−9) (Supplementary Fig. 1A and B), consistent with previous reports (18,19). The same lead variants and LD proxies most strongly associated with type 2 diabetes and glucose-related traits are most strongly associated with ADCY5 mRNA level, strongly suggesting that these signals are coincident. rs11708067 is not associated with mRNA level (FDR >5%) of any other gene within 1 Mb. No eQTL was observed in GTEx V6p at FDR ≤5% for other tissues.
We evaluated the association between rs11708067 and 21 glycemic traits in ∼8,500 Finnish males without diabetes from the METSIM study. Among these traits, we observed the strongest evidence of associations of the rs11708067-A allele with three measures of plasma glucose concentration after an oral glucose challenge: glucose area under the curve between 30 min and 2 h (N = 8,591, β = 0.086, P = 4.3 × 10−5), glucose concentration after 2 h (N = 8,620, β = 0.084, P = 5.8 × 10−5), and glucose area under the curve between 0 min and 2 h (N = 8,591, β = 0.083, P = 8.4 × 10−5) (Supplementary Table 2).
rs11708067 Overlaps Evidence of Open Chromatin and a Predicted Enhancer Region in Human Pancreatic Islets
We next aimed to determine whether the reported GWAS lead variants or linked variants (proxies based on LD) showed evidence of allelic differences in regulatory activity. In 1000 Genomes phase 3 EUR data (7), 17 variants exhibited strong LD (r2 >0.8) with any of the lead type 2 diabetes or glycemic trait lead GWAS variants (rs11708067, rs11717195, rs2877716, rs9883204). Each of the 18 candidates, including rs11708067, is located within ADCY5 introns (Fig. 1). We used DNase-, FAIRE-, and ATAC-seq data in primary human pancreatic islets to investigate potential regulatory regions overlapping variants at this locus (20–22). Four of the 18 candidate variants overlap predicted enhancer regions (28) in human pancreatic islets (Supplementary Table 4). Of these, rs11708067 overlaps an ATAC-seq peak in pancreatic islets (20,35) and is the only variant of the 18 to overlap both a DNase and FAIRE peak in human pancreatic islets (Fig. 1 and Supplementary Table 4).
rs11708067 overlaps pancreatic islet open chromatin and a strong enhancer chromatin state in ADCY5 intron 3. The 50-kb region shown includes all 17 variants in strong LD (r2 >0.8, 1000 Genomes phase 3 EUR) with the type 2 diabetes–associated single nucleotide polymorphism rs11708067 (in boldface type) (18 total noncoding candidate variants located within ADCY5 introns 1–3). Selected Encyclopedia of DNA Elements (ENCODE) open chromatin (DNase, FAIRE) (21,22), chromatin accessibility (ATAC-seq) (20), and RNA-seq and chromatin state tracks (28) are shown. The black vertical line above the pancreatic islet open chromatin tracks represents the 231-bp segment containing rs11708067 that was tested in luciferase reporter assays. GM12878, lymphoblastoid cells from B cells; H1 ES, embryonic stem cells; HepG2, human hepatocellular carcinoma cells; HMEC, human mammary epithelial cells; HSMM, human skeletal muscle myoblasts; HUVEC, human umbilical vein endothelial cells; islets, pancreatic islets; K562, leukemia cell line; NHEK, epidermal keratinocytes; NHLF, lung fibroblasts.
rs11708067 overlaps pancreatic islet open chromatin and a strong enhancer chromatin state in ADCY5 intron 3. The 50-kb region shown includes all 17 variants in strong LD (r2 >0.8, 1000 Genomes phase 3 EUR) with the type 2 diabetes–associated single nucleotide polymorphism rs11708067 (in boldface type) (18 total noncoding candidate variants located within ADCY5 introns 1–3). Selected Encyclopedia of DNA Elements (ENCODE) open chromatin (DNase, FAIRE) (21,22), chromatin accessibility (ATAC-seq) (20), and RNA-seq and chromatin state tracks (28) are shown. The black vertical line above the pancreatic islet open chromatin tracks represents the 231-bp segment containing rs11708067 that was tested in luciferase reporter assays. GM12878, lymphoblastoid cells from B cells; H1 ES, embryonic stem cells; HepG2, human hepatocellular carcinoma cells; HMEC, human mammary epithelial cells; HSMM, human skeletal muscle myoblasts; HUVEC, human umbilical vein endothelial cells; islets, pancreatic islets; K562, leukemia cell line; NHEK, epidermal keratinocytes; NHLF, lung fibroblasts.
Allelic Imbalance of H3K27ac Reads Overlapping rs11708067 in Human Pancreatic Islets
Acetylation of H3K27ac marks active enhancer regions (36,37). Thus, we evaluated H3K27ac ChIP-seq data in a primary human islet sample heterozygous at rs11708067 and investigated the alleles in sequencing reads spanning rs11708067. The H3K27ac ChIP-seq reads showed evidence of allelic imbalance, with 9 reads containing rs11708067-A reference allele compared with 69 reads containing the rs11708067-G allele (Fig. 2) (exact binomial P = 1.2 × 10−14). These data suggest that in the islet cellular chromatin context, rs11708067 is located within an enhancer region marked by H3K27ac.
Evidence of allelic imbalance in H3K27ac ChIP-seq reads at rs11708067. H3K27ac ChIP-seq reads in a primary human islet sample heterozygous at rs11708067. Blue indicates reads that contain the G allele of rs11708067, red indicates reads that contain the A allele of rs11708067, and gray indicates reads in the peak that do not overlap rs11708067.
Evidence of allelic imbalance in H3K27ac ChIP-seq reads at rs11708067. H3K27ac ChIP-seq reads in a primary human islet sample heterozygous at rs11708067. Blue indicates reads that contain the G allele of rs11708067, red indicates reads that contain the A allele of rs11708067, and gray indicates reads in the peak that do not overlap rs11708067.
rs11708067 Shows Allelic Differences in Enhancer Activity
To investigate whether rs11708067 shows allelic differences in transcriptional regulatory activity, we tested a 231-bp DNA segment containing rs11708067 using luciferase reporter assays in 832/13 rat insulinoma and MIN6 mouse cell lines. We examined each of the rs11708067 alleles in the sequence in both orientations with respect to a minimal promoter. In 832/13 cells, we observed 1.6-fold greater enhancer activity for the A allele compared with empty vector control and 3.3-fold greater enhancer activity for the G allele compared with an empty vector control in the forward orientation, with significant differences between the alleles (Fig. 3A, P = 0.0007). In the reverse orientation, we observed 0.84-fold luciferase activity (A allele) compared with empty vector control and 1.3-fold enhancer activity (G allele) compared with empty vector control, with significant differences between the alleles (Fig. 3A, P = 0.003). In MIN6 cells, we observed 2.6-fold greater enhancer activity for the A allele compared with empty vector control and 8.3-fold greater enhancer activity for the G allele compared with empty vector control in the forward orientation, with significant differences between the alleles (Supplementary Fig. 2) (P = 1.9 × 10−3). In the reverse orientation, we observed 2.3-fold (A allele) and 13.9-fold (G allele) greater enhancer activity compared with empty vector control with significant differences between the alleles (Supplementary Fig. 2) (P = 4.9 × 10−5). The lower transcriptional activity of rs11708067-A is consistent with the lower expression of rs11708067-A in the eQTL association with ADCY5 in human islets. In both cell types and both orientations, the rs11708067-G allele showed higher luciferase activity, suggesting that rs11708067 can have an allelic effect on transcriptional activity.
rs11708067 exhibits allelic differences in transcriptional activity and protein binding in 832/13 rat insulinoma cells. A: The 231-bp segments containing allele A or G of rs11708067 were cloned into a pGL4.23 luciferase reporter vector upstream of the minimal promoter in both orientations (forward and reverse). The relative luciferase activities of plasmids transfected into 832/13 cells normalized to an empty vector control are shown on the y-axis. Error bars represent SD of three to five independent clones per allele (t tests). Black triangles, empty vector; black circles, rs11708067-A allele; white circles, rs11708067-G allele. B: EMSAs with biotin-labeled probes containing either the A or G allele of rs11708067. Probes were incubated with 2.2 µg 832/13 rat insulinoma nuclear lysate. The arrows indicate differential protein binding to the A allele, which is competed away by 240-fold excess competitor DNA containing the A allele (lane 3). EMSAs were performed six separate times with consistent results. Antibodies to CEBPB and NKX2.2 (8 µg) were tested for supershifts in lanes 4 and 5, respectively, and lanes 9 and 10, respectively.
rs11708067 exhibits allelic differences in transcriptional activity and protein binding in 832/13 rat insulinoma cells. A: The 231-bp segments containing allele A or G of rs11708067 were cloned into a pGL4.23 luciferase reporter vector upstream of the minimal promoter in both orientations (forward and reverse). The relative luciferase activities of plasmids transfected into 832/13 cells normalized to an empty vector control are shown on the y-axis. Error bars represent SD of three to five independent clones per allele (t tests). Black triangles, empty vector; black circles, rs11708067-A allele; white circles, rs11708067-G allele. B: EMSAs with biotin-labeled probes containing either the A or G allele of rs11708067. Probes were incubated with 2.2 µg 832/13 rat insulinoma nuclear lysate. The arrows indicate differential protein binding to the A allele, which is competed away by 240-fold excess competitor DNA containing the A allele (lane 3). EMSAs were performed six separate times with consistent results. Antibodies to CEBPB and NKX2.2 (8 µg) were tested for supershifts in lanes 4 and 5, respectively, and lanes 9 and 10, respectively.
Alleles of rs11708067 Exhibit Differential Protein Binding
To determine whether transcription factors show allelic differences in binding to DNA at rs11708067, we performed EMSAs with 21-bp biotinylated probes containing either the rs11708067-A or rs11708067-G alleles using nuclear lysate from 832/13 or MIN6 cells. We observed four bands with greater protein binding for the rs11708067-A allele compared with the rs11708067-G allele (Fig. 3B, lane 2 vs. lane 7, and Supplementary Fig. 3B, lane 2 vs. lane 7). Addition of unlabeled probe containing the A allele of rs11708067 competed away the stronger A-allele bands (Fig. 3B, lane 3, and Supplementary Fig. 3A and B, lanes 3). Using bioinformatic prediction tools JASPAR (38), TRANSFAC (39,40), and The Islet Regulome Browser (41) to search for transcription factor motif matches, we found that the DNA region spanning rs11708067 overlaps a NKX2.2 transcription factor ChIP-seq peak (41) and a predicted binding motif for CEBP. However, incubation of the EMSA reaction with antibodies to CEBPB or NKX2.2 did not result in supershifts (Fig. 3B, lanes 4 and 5, respectively, and lanes 9 and 10, respectively; Supplementary Fig. 3B, lanes 5 and 10). Additional antibodies we tested in EMSAs, including FOXA2, PDX1, p300, NKX6.1, MAFB, YY1, and RXRA, also did not demonstrate evidence of supershifts (Supplementary Fig. 4). In summary, multiple undefined proteins from both mouse and rat β-cell nuclear lysates bind specifically to the A allele of rs11708067.
Deletion of the Enhancer Element Reduced Adcy5 Expression and Insulin Secretion
To investigate whether the enhancer element affects the function of Adcy5 in 832/13 cells, we generated deletions of the enhancer element using CRISPR-Cas9 technology. We successfully generated three clones with homozygous deletions of the enhancer element (Fig. 4A and Supplementary Figs. 5–8) and measured expression of Adcy5 and nearby gene Sec22a. Compared with five control clones with an intact enhancer, the three clones harboring homozygous deletions of the enhancer element exhibited no effect on expression of Sec22a (Fig. 4B) (P = 0.52) but did show a mean 64% reduction in Adcy5 expression (Fig. 4C) (P = 0.002). We then investigated whether deletion of the enhancer element led to effects on insulin secretion. The three homozygous deletion-containing clones showed a mean 39% reduction of insulin secretion at basal 3 mmol/L glucose concentration compared with the five control clones (Fig. 4D, P = 0.0006; Supplementary Fig. 9). Effects on insulin secretion after treatment with 18 mmol/L glucose were variable (Supplementary Fig. 10). These results show that the enhancer element containing rs11708067 affects Adcy5 expression and furthermore can alter the role of Adcy5 in insulin secretion.
Deletion of the orthologous enhancer element in 832/13 cells leads to reduced Adcy5 expression and reduced insulin secretion. A: Three homozygous deletions (indicated by the genomic regions in between the Xs) of the orthologous Adcy5 enhancer region were generated using CRISPR-Cas9 technology in 832/13 rat insulinoma cells. UCSC Genome Browser chromosome 11 coordinates are shown for the rat genome rn5. PhyloP Cons indicates basewise conservation with 13 species (rat, mouse, guinea pig, human, chimp, rhesus, cow, dog, panda, opossum, chicken, turkey, and zebrafish) based on the phyloP track from the UCSC Genome Browser. Gene expression normalized to B2m control (shown on the y-axes for B and C) or insulin secretion (shown on the y-axis for D) was measured in 832/13 rat β-cells. B: Normalized Sec22a nearby control gene expression (displayed as a percent of average expression of intact control clones). C: Normalized Adcy5 expression (displayed as a percent of average expression of intact control clones). D: Insulin secretion at 3 mmol/L glucose concentration (displayed as pg insulin/µg total protein). B and C: For gene expression, error bars indicate the SD of three to five clones per sample for two independent experiments. D: Insulin secretion assays were performed in duplicate or triplicate, and error bars indicate the SD of two to five clones per sample. Unpaired t tests were performed to compare gene expression or insulin secretion of enhancer homozygous deletion clones with intact control clones. Black circles, wild-type (WT) 832/13 clones; black triangles, intact control clones; white circles, enhancer homozygous deletion clones.
Deletion of the orthologous enhancer element in 832/13 cells leads to reduced Adcy5 expression and reduced insulin secretion. A: Three homozygous deletions (indicated by the genomic regions in between the Xs) of the orthologous Adcy5 enhancer region were generated using CRISPR-Cas9 technology in 832/13 rat insulinoma cells. UCSC Genome Browser chromosome 11 coordinates are shown for the rat genome rn5. PhyloP Cons indicates basewise conservation with 13 species (rat, mouse, guinea pig, human, chimp, rhesus, cow, dog, panda, opossum, chicken, turkey, and zebrafish) based on the phyloP track from the UCSC Genome Browser. Gene expression normalized to B2m control (shown on the y-axes for B and C) or insulin secretion (shown on the y-axis for D) was measured in 832/13 rat β-cells. B: Normalized Sec22a nearby control gene expression (displayed as a percent of average expression of intact control clones). C: Normalized Adcy5 expression (displayed as a percent of average expression of intact control clones). D: Insulin secretion at 3 mmol/L glucose concentration (displayed as pg insulin/µg total protein). B and C: For gene expression, error bars indicate the SD of three to five clones per sample for two independent experiments. D: Insulin secretion assays were performed in duplicate or triplicate, and error bars indicate the SD of two to five clones per sample. Unpaired t tests were performed to compare gene expression or insulin secretion of enhancer homozygous deletion clones with intact control clones. Black circles, wild-type (WT) 832/13 clones; black triangles, intact control clones; white circles, enhancer homozygous deletion clones.
Discussion
Identification of the functional regulatory variants at GWAS loci and the genes they regulate can help elucidate the mechanisms by which variants alter gene function. In this study, we showed evidence of a functional regulatory variant at the ADCY5 GWAS for type 2 diabetes (3,4) and glycemic traits (10–12,42). We confirmed that the type 2 diabetes risk alleles are associated with lower expression of ADCY5 in human pancreatic islets (18–20) and showed that candidate variant rs11708067 exhibits allelic differences in regulatory activity. Specifically, the rs11708067-A allele showed lower transcriptional activity, stronger binding of nuclear proteins, and lower enhancer activity in the native chromatin context of human islets based on H3K27ac islet ChIP-seq reads. Additionally, CRISPR-Cas9–mediated homozygous deletion of the enhancer element in rat insulinoma cells diminished Adcy5 gene expression and reduced insulin secretion. We propose that the rs11708067-A allele binds to one or more transcriptional repressor proteins and reduces ADCY5 transcriptional activity, leading to reduced insulin secretion and increased type 2 diabetes risk.
We provide evidence that rs11708067 is a regulatory variant at the ADCY5 locus. Our data are consistent with the fact that rs11708067 was the strongest candidate for an effect on regulatory function based on its location in a putative islet enhancer element. rs11708067 also showed the highest posterior probability (highest log-Bayes factor 12.39) of driving the association signal from a credible set analysis (11). Considering that some GWAS loci have multiple functional regulatory variants responsible for the association (43), other variants at this locus may also contribute to allelic differences in transcriptional activity.
Despite evidence of a NKX2.2 ChIP-seq peak in human pancreatic islets and a predicted binding motif for CEBP overlapping rs11708067, we did not observe evidence of supershifts using antibodies to NKX2.2 and CEBPB. These factors may not bind to the sequence spanning rs11708067, or overexpression of these proteins may be necessary to elicit a supershift. The absence of a supershift may also be attributed to taking the protein and DNA interactions out of their genomic DNA and chromatin context.
ADCY5 is the strongest candidate gene underlying the metabolic trait associations at this locus. ADCY5 is strongly expressed in human pancreatic islets (18,44) and eQTL data showed a strong coincident association between the GWAS variant rs11708067 and ADCY5 expression in islets (19,20). ADCY5 is also highly expressed in H1 human embryonic stem cells (Fig. 1) (28), which may suggest a role for ADCY5 expression in development, consistent with an association between an ADCY5 variant and birth weight (14,45). The fetal insulin hypothesis supports the reported associations between the type 2 diabetes risk alleles and low birth weight (46). This hypothesis suggests that genetic variants contribute to decreased insulin secretion by the fetus and/or decreased insulin sensitivity in fetal tissues. Considering that insulin is a growth factor, these reductions in insulin secretion and/or insulin sensitivity ultimately lead to reduced fetal growth and birth weight (47). ADCY5 may also act in other tissues that could contribute to biological mechanisms of type 2 diabetes; for example, alterations in TCF7L2 type 2 diabetes–associated gene expression alter glucose and insulin metabolism in liver and pancreatic islets (48–50).
Adenylate cyclase 5 has been shown to couple glucose signals to cAMP production, increase Ca2+, and influence β-cell connectivity (18). Our METSIM study trait association data support a role for rs11708067 in glucose homeostasis after a glucose stimulus (Supplementary Table 2). However, in our experiments, deletion of the enhancer element decreased Adcy5 expression and showed lower insulin secretion at basal (3 mmol/L) glucose concentration (Fig. 4D and Supplementary Fig. 9), but we did not observe consistent alterations in insulin secretion at 18 mmol/L glucose (Supplementary Fig. 10). These contrasting results may be attributed to our deletion of only a single regulatory element, the variable responsiveness of 832/13 cells to glucose stimulation, or that the deletion may not perfectly represent the human regulatory element containing rs11708067. Nonetheless, the effect of deleting the regulatory element on gene expression supports the role of this noncoding sequence in the cellular genomic context.
Taken together, our data identify a functional variant at the ADCY5 locus, rs11708067, and provide evidence of a potential mechanism by which this variant affects expression of ADCY5, leading to altered fasting glucose levels and risk of type 2 diabetes (Supplementary Fig. 11). The type 2 diabetes risk allele rs11708067-A is likely to bind multiple repressor proteins that weaken enhancer activity and reduce ADCY5 expression. These experimental data provide consistent evidence supporting the role of a noncoding variant in the molecular basis of type 2 diabetes.
Article Information
Acknowledgments. The authors thank Lori Bonnycastle and Amy Swift (National Human Genome Research Institute) for providing helpful advice regarding CRISPR-Cas9 experiments, Peter Chines (National Human Genome Research Institute) for performing transcription factor binding predictions, Jennifer Kulzer (University of North Carolina at Chapel Hill) and Kyle Gaulton (University of California, San Diego) for helpful discussions, the UNC Flow Cytometry Core Facility for cell sorting in CRISPR-Cas9 experiments, the GTEx project for gene expression and quantitative trait loci data, and the Encyclopedia of DNA Elements (ENCODE) and National Institutes of Health Roadmap Epigenomics Consortia, data analysis and coordination centers, and production laboratories that generated the chromatin and ChIP sequencing data used for variant annotation.
Funding. This study was supported by National Institutes of Health grants T32HL069768 (T.S.R.), T32GM007092 (T.S.R.), F31HL127984 (M.E.C.), R01DK072193 (K.L.M.), U01DK105561 (K.L.M.), R01DK093757 (K.L.M.), T32GM067553 (M.L.B.), R00DK092251 (M.L.S.), and U01DK062370 (M.B.); National Human Genome Research Institute Division of Intramural Research project number Z01HG000024 (F.S.C.); American Heart Association Predoctoral Fellowship 13PRE16930025 (M.L.B.); Academy of Finland grants 77299 and 124243 (M.L.); the Finnish Heart Foundation (M.L.); the Finnish Diabetes Foundation (M.L.); Finnish Funding Agency for Technology and Innovation (TEKES) contract 1510/31/06 (M.L.); and the Commission of the European Community HEALTH-F2-2007-201681 (M.L.). The UNC Flow Cytometry Core Facility is supported in part by the Lineberger Comprehensive Cancer Center, University of North Carolina (P30 CA016086 Cancer Center Core support grant). Research reported in this article was supported in part by the North Carolina Biotech Center institutional support grant 2005-IDG-1016. The GTEx project was supported by the Common Fund of the Office of the Director of the National Institutes of Health, National Cancer Institute, National Human Genome Research Institute, National Heart, Lung, and Blood Institute, National Institute on Drug Abuse, National Institute of Mental Health, and National Institute of Neurological Disorders and Stroke. The data used for the analyses described in this article were obtained from the GTEx Portal on 23 May 2017.
Duality of Interest. No potential conflicts of interest relevant to this article were reported.
Author Contributions. T.S.R. and K.L.M. wrote and edited the manuscript and contributed to experimental design. T.S.R. performed luciferase reporter and EMSA and gene expression experiments. M.E.C. contributed to experimental design and edited the manuscript. S.V. performed luciferase reporter, CRISPR-Cas9, and insulin secretion assays and edited the manuscript. M.L.B., B.N.W., and A.V. performed regulatory data analyses. R.P.W. performed eQTL analyses. M.A.M. performed H3K27ac ChIP experiments. G.J.K. contributed to experimental design. R.K. edited the manuscript. Y.W. and A.U.J. performed association analyses. National Institutes of Health Intramural Sequencing Center (NISC) Comparative Sequencing Program contributed sequencing for H3K27ac ChIP-seq experiments. M.R.E. performed human islet collection and contributed experimental data. J.K. and M.L. contributed METSIM association data. L.J.S. and M.B. contributed islet regulatory data and edited the manuscript. F.S.C. contributed islet regulatory data. S.C.J.P. contributed islet regulatory data and discussion and edited the manuscript. M.L.S. performed H3K27ac ChIP experiments, contributed islet regulatory data and discussion, and edited the manuscript. K.L.M. is the guarantor of this work and, as such, had full access to all the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.
Prior Presentation. An abstract describing this work was presented as a poster presentation at the annual meeting of The American Society of Human Genetics, Baltimore, MD, 6–10 October 2015.