A ChIP-seq defined genome-wide map of vitamin D receptor binding: associations with disease and evolution

Abstract

Initially thought to play a restricted role in calcium homeostasis, the pleiotropic actions of vitamin D in biology and their clinical significance are only now becoming apparent. However, the mode of action of vitamin D, through its cognate nuclear vitamin D receptor (VDR), and its contribution to diverse disorders, remain poorly understood. We determined VDR binding throughout the human genome using chromatin immunoprecipitation followed by massively parallel DNA sequencing (ChIP-seq). After calcitriol stimulation, we identified 2776 genomic positions occupied by the VDR and 229 genes with significant changes in expression in response to vitamin D. VDR binding sites were significantly enriched near autoimmune and cancer associated genes identified from genome-wide association (GWA) studies. Notable genes with VDR binding included IRF8, associated with MS, and PTPN2 associated with Crohn's disease and T1D. Furthermore, a number of single nucleotide polymorphism associations from GWA were located directly within VDR binding intervals, for example, rs13385731 associated with SLE and rs947474 associated with T1D. We also observed significant enrichment of VDR intervals within regions of positive selection among individuals of Asian and European descent. ChIP-seq determination of transcription factor binding, in combination with GWA data, provides a powerful approach to further understanding the molecular bases of complex diseases.

Figure 1.

VDR binding intervals and genomic location. (A) Unstimulated cells. (B) Calcitrol-stimulated cells. Enrichment is shown by location with respect to gene structure (y-axis) and binding strength (peak value, maximum number of reads aligned to a genomic position within a ChIP-seq interval) (x–axis). Pie charts summarize the percentage of VDR binding sites based on location. Intergenic regions were defined as at least 5 kb away from the first or last exon of a gene, upstream (promoter) regions defined as within 5 kb of the transcriptional start site, and downstream regions as within 5 kb from the end of the last exon.

Figure 2.

MEME motif analysis for VDR intervals following calcitriol stimulation. The top-scoring motif found by MEME resembles the known VDR element. The nucleotide frequencies of the genomic sequences aligned at the motif are shown in a sequence logo representation (Schneider and Stephens 1990).

Figure 3.

Common traits showing enrichment of VDR binding within intervals identified by GWAS. A total of 47 common diseases and traits were analyzed (see Methods and Supplemental Table 5) and those showing significant enrichment of VDR binding defined by ChIP-seq in two LCLs after calcitriol stimulation with a 1% FDR are shown.

Figure 4.

VDR ChIP-seq analysis for IRF8 and PTPN2. VDR ChIP-seq data shown for two biological replicates of two LCLs, GM10855 and GM10861, either resting or after induction with calcitriol for 36 h. (A) Tracks shown for IRF8 (chr16: 84,467,000–84,537,000) with a novel site of VDR occupancy noted at +3.8 kb relative to the transcriptional start site (TSS), as well as weaker sites at −10.2 kb and +4.7 kb. (B) Validation of VDR binding by ChIP for GM07019, GM07348, and GM10854 analyzed by quantitative real-time PCR. Mean fold difference (±SD) in enrichment of each of the PCR amplicons is expressed relative to input DNA. (C) Tracks also shown for PTPN2 (chr18: 12,750,000–12,900,000) with VDR occupancy in intron 1 (+16.5 kb relative to the TSS of PTPN2) and intron 7 (+80.8 kb) with (D) validation by ChIP for GM07019, GM07348, and GM10854.

Figure 5.

VDR ChIP-seq analysis for RASGRP3, PRKCQ, TNFSF11, and CLPTM1L. VDR ChIP-seq data shown for two biological replicates of two LCLs, GM10855 and GM10861, either resting or after induction with calcitriol for 36 h. (A) RASGRP3 and flanking sequences (chr2: 33,500,000–33,650,000) with rs13385731 (GWAS SNP marker for SLE) located in an intronic VDR binding site; (B) PRKCQ and flanking sequences (chr10: 6,400,000–6,700,000) with rs947474 (GWAS SNP marker for T1D) located in the VDR binding site 78 kb downstream from PRKCQ; (C) TNFSF11 and flanking sequences (chr13: 41,830,000–42,100,000) with rs9594738 (GWAS SNP marker for bone mineral density) 185 kb upstream of TNFSF11; (D) CLPTM1L and flanking sequences (chr5: 1,300,000–1,440,000) with rs4975616 (GWAS SNP marker for lung cancer). Also shown are ChIP-seq and DNase-seq data for GM12878 generated by the ENCODE Project (The ENCODE Project Consortium 2007).

Figure 6.

VDR binding, histone H1 gene cluster, selection, and rs10946808. (A) Enrichment of VDR binding intervals in histone HIST1 gene cluster on chromosome 6p21-22 (chr6: 25,900,000–28,100,000). (B) The region chr6: 26,300,001–26,400,000 reported as showing evidnce of selection includes rs10946808 located 3′ of HIST1H1D. (C) Data from the Human Genome Diversity Project show rs10946808 has a higher minor allele frequency in Asians.