Function Annotation of Hepatic Retinoid x Receptor α Based on Genome-Wide DNA Binding and Transcriptome Profiling

Background Retinoid x receptor α (RXRα) is abundantly expressed in the liver and is essential for the function of other nuclear receptors. Using chromatin immunoprecipitation sequencing and mRNA profiling data generated from wild type and RXRα-null mouse livers, the current study identifies the bona-fide hepatic RXRα targets and biological pathways. In addition, based on binding and motif analysis, the molecular mechanism by which RXRα regulates hepatic genes is elucidated in a high-throughput manner. Principal Findings Close to 80% of hepatic expressed genes were bound by RXRα, while 16% were expressed in an RXRα-dependent manner. Motif analysis predicted direct repeat with a spacer of one nucleotide as the most prevalent RXRα binding site. Many of the 500 strongest binding motifs overlapped with the binding motif of specific protein 1. Biological functional analysis of RXRα-dependent genes revealed that hepatic RXRα deficiency mainly resulted in up-regulation of steroid and cholesterol biosynthesis-related genes and down-regulation of translation- as well as anti-apoptosis-related genes. Furthermore, RXRα bound to many genes that encode nuclear receptors and their cofactors suggesting the central role of RXRα in regulating nuclear receptor-mediated pathways. Conclusions This study establishes the relationship between RXRα DNA binding and hepatic gene expression. RXRα binds extensively to the mouse genome. However, DNA binding does not necessarily affect the basal mRNA level. In addition to metabolism, RXRα dictates the expression of genes that regulate RNA processing, translation, and protein folding illustrating the novel roles of hepatic RXRα in post-transcriptional regulation.


Introduction
Retinoid x receptor (RXR) plays a critical role in metabolism, development, differentiation, proliferation, and cell death by regulating gene expression [1,2]. The expression profile of RXRs and their downstream signaling is altered in a variety of diseases including breast cancer [3] and viral hepatitis [4]. Correspondingly, RXR agonists are implicated in cancer prevention, antiviral therapy, dermatological disease, and metabolic syndromes [2,5,6]. Mediated via RXRs (a, b, and c) and retinoic acid receptors (a, b, and c), retinoic acid exerts its biological effects. In addition to RARs, RXRa is crucial for many other receptors to work. The receptors for fatty acids, oxysterols, bile acids, vitamin D, etc. form dimers with RXRs to regulate gene transcription. RXRa is essential for fetal morphogenesis [7]. Global RXRa knockout is embryonically lethal and the embryos develop myocardiac hypoplasia with reduced liver size at gestational day 12.5 [8]. RXRa is predominantly expressed in the liver [7]. Hepatocyte RXRa knockout mice were produced to study its role in the liver [9]. RXRa-deficient mice have compromised lipid [10], carbohydrate [11], xenobiotic [12], and amino acid homeostasis [13]. In addition, hepatocyte RXRa is also implicated in liver steatosis and inflammation as well as regeneration [14][15][16][17]. Specifically, hepatocyte RXRa-deficient mice are more susceptible than wild type mice to alcohol and non-alcohol-induced steatohepatitis [13][14][15]. Partial hepatectomy-induced liver regeneration is also hampered due to hepatic RXRa deficiency [16]. Liver is the organ that stores and converts retinol to its biological active form i.e. retinoic acid. Liver also produces binding proteins to deliver retinoic acid to target sites. Furthermore, hepatocytes express high level of RXRa. Thus, it is pivotal to understand the biological function of RXRa in the liver.
The current study determines genome-wide RXRa binding in normal mouse livers by chromatin immunoprecipitation using specific anti-RXRa antibody followed by next generation sequencing (ChIP-seq). In addition, microarray was performed to identify genes that are differentially expressed in wild type and RXRa-null mouse livers. Combining the two datasets, we established the relationship between hepatic RXRa-DNA binding and RXRa-dependent gene expression. This global profiling of RXRa binding along with gene expression not only allows us to capture all RXRa downstream targets and pathways, but also helps us to understand the molecular mechanism by which hepatic RXRa regulates gene expression.

Animals and Tissues
The generation of hepatocyte RXRa-deficient mice was described in previous publications [9,10]. The LoxP sequences were inserted into introns flanking the fourth exon of the RXRa gene and DNA binding domain of the gene is deleted after crossing the floxed RXRa allele against a transgenic line in which the cre recombinase is expressed under the control of the albumin promoter. The mutant mice express a truncated protein that lacks the DNA binding domain. However, the ligand binding domain remains intact. Mouse livers were collected at 10 AM from 12week-old wild type (WT) and hepatocyte RXRa knockout (KO) mice, which have C57BL/6 genetic background. Tissues were snap frozen in liquid nitrogen and kept at 280uC for future use. Animal protocols and procedures were approved by the Institutional Animal Care and Use Committee (IACUC) at the University of Kansas Medical Center and the University of California, Davis.

RNA Preparation and Microarray
Total RNA was extracted using TRIzol Reagent (Invitrogen Co., CA) and purified with the RNeasy Mini Kit (Qiagen Inc., CA). The quantity and quality of the total RNA were assessed by Bioanalyzer 2100 (Agilent Technologies, CA). Complementary DNA was made using High Capacity RNA-to-cDNA Kit (Applied Biosystems, CA). Affymetrix chips (MOE 430A 2.0) that covered about 14,000 mouse genes were used. Microarray (n = 3 per group) and data processing as well as the methods used for data validation were described in our previous publication [18]. By Student t-test, genes differentially expressed in wild type and RXRa-knockout mouse livers with p value less than 0.05 were identified as RXRa-dependent genes.
Chromatin immunoprecipitation sequencing (ChIP-Seq) and ChIP-qPCR Frozen livers were fixed in 1% formaldehyde (pH = 7) for 15 minutes before being quenched with 0.125 M glycine. Following cell lysis, the nuclear fraction was extracted and sonicated to produce 300-500 base-pair (bp) DNA fragments. Genomic DNA (Input) was prepared by treating aliquots of chromatin with RNase, proteinase K and heated for de-crosslinking, followed by ethanol precipitation. Chromatin (30 mg) was precleared by Dynase beads (Invitrogen Co., CA) before incubation with a ChIP-quality anti-RXRa antibody (Santa Cruz, CA). An antibody to IgG (Santa Cruz, CA) and RNA Pol II (Millipore, MA) was used as negative and positive control, respectively. Samples were incubated with prepared Dynase beads at 4uC overnight, followed by de-crosslinking and purification.

Data Analysis
Primary image analysis and base calling were performed using Genome Analyzer Pipeline Software (Illumina Inc., CA). All sequenced reads were aligned to mm9 mouse reference genome using bowtie version 0.12.7 [19]. Only uniquely mapped reads were included. Regions with reads enrichment were detected using the Model-based Analysis of ChIP-Seq (MACS v 1.4.1) method [20]. Non-specific peaks with false discovery ratio (FDR) greater than 0.1 were eliminated by comparing them with the IgG background. Peaks were further split by Mali Salmon's Peak Splitter (http://www.ebi.ac.uk/bertone/software.html) and filtered by p-value of Poisson distribution lower than 10 25 . The peak annotation was based on the UCSC genome NCBI/mm9 database. Relative to the transcription start site (TSS) or 39 end of the gene, the peak locations were defined as ingene (summit located within the coding region), promoter (22 kb,0 bp to the TSS), downstream (0,100 kb to the 39 end), upstream (2100 kb , 22 kb to the TSS), or intergenic (more than 100 kb away from TSS or 39 end).

Motif and Pathway Analysis
The sequences that were 100 bp up and downstream from the summits of the top 500 peaks, which had the highest peak scores, Figure 1. Global analysis of ChIP-seq and microarray data. (A) Venn diagrams of RXRa-bound genes (III+IV+V), hepatic genes (I+II+III+IV), and RXRa-dependent genes (II+III). RXRa binding location for (B) RXRa-dependent, (C) RXRa-independent hepatic expressed, and (D) non-hepatic expressed genes. In gene: peak summit located within the coding region; Promoter: peak summit located within 22 kb ,0 bp to the transcription start site (TSS); Upstream: peak summit located within -100 kb,22 kb to the TSS; Downstream: peak summit located within 0 bp ,100 kb to 39 end of the gene; Intergenic: peak summit located outside the above mentioned regions. Red circle: RXRa-dependent genes; Blue circle: hepatic genes; Yellow circle: RXRa-bound genes. I: RXRa-independent genes and lack of RXRa binding; II: RXRa-dependent genes but lack of RXRa binding; III: RXRadependent genes that have RXRa binding; IV: RXRa-independent genes that have RXRa binding; V: non-hepatic genes that have RXRa binding. doi:10.1371/journal.pone.0050013.g001 were subjected to motif analysis using by MEME-ChIP (Multiple EM for Motif Elicitation) [21]. Furthermore, a Hidden Markov Model (HMM) was established based on nuclear receptor binding sites from published work [22] and JASPAR CORE (http:// jaspar.cgb.ki.se/cgi-bin/jaspar db.pl) to predict the specific motif information for all of the RXRa binding sites. All of the biological function and pathway analyses were performed using the Functional Annotation Tool in DAVID (DAVID; http://www. david.niaid.nih.gov). Functional pathways or process with p,0.05 and Bonferroni value ,0.1 were accepted.

Statistical Analysis
For ChIP-qPCR data and microarray data, the difference between two groups was analyzed by Student's t test. P,0.05 was considered statistically significant.

Global RXRa binding and RXRa-dependent gene expression in mouse livers
Microarray data showed that there were 9068 genes with detectable signals (mRNA levels), which are referred to as ''hepatic genes''. Among these, 768 were significantly up-regulated (UIK: up-regulated in KO) and 696 genes were down-regulated (DIK: down-regulated in KO) in the hepatocyte RXRa KO mouse livers. Thus, 16.14% of hepatic genes (1464 out of 9068) were expressed in an RXRa-dependent manner.
There were 109971 confident RXRa peaks detected in 14816 genes in the mouse liver. The average is 7.4 peaks per gene. Remarkably, 79.1% (7174) of the hepatic genes were bound by RXRa ( Figure 1A), which includes 657 (85.5%) UIK and 605 (86.9%) DIK genes. However, 7642 RXRa-bound genes did not have significant signals by microarray in either WT or RXRa KO mouse livers and we refer to those genes as ''non-hepatic genes''.
The frequencies of promoter peaks are similar for RXRadependent (12.2%) and -independent (12.6%) genes ( Figure 1B and C). RXRa binding sites occurred more frequently in promoter (12.6% vs. 8.6%) and coding (49.4% vs. 39.4%) regions for the hepatic than non-hepatic genes ( Figure 1C and D). Conversely, RXRa bound more frequently to the intergenic, up-and downstream regions for the non-hepatic than the hepatic genes. Thus, it seems that the binding location of RXRa could be a determinant factor for the level of hepatic gene expression.
In terms of chromosomal distribution, RXRa binding occurred more frequently on chromosome 5 and 11 (956 and 1115 RXRa binding sites, respectively) than others ( Figure 2). Notably, gene function analysis demonstrated that 32% of RXRa-bound genes on chromosome 5 are associated with alternative splicing, while 52% of RXRa-bound genes on chromosome 11 encode phosphorproteins suggesting the role of RXRa in post-transcriptional modification. In addition, seven RXRa binding sites (peak scores ranging170-470) were mapped on mitochondrial DNA in the areas where genes were not identified by the UCSD Gene Browser. 171 peaks were matched on X chromosome, while only 1 peak was noted on Y chromosome.

Motif Analysis of RXRa binding sites
Since RXRa dimerizes with multiple nuclear receptors, its binding sites are composed of diverse binding motifs. The result predicted by the Hidden Markov Model identified DR1 (12.5%) as the most common motif for RXRa binding in mouse liver followed by DR4, inverted repeat (IR) 1, DR3, and IR3. Furthermore, spacers 1, 3, and 4 were relatively more prevalent than other spacers ( Figure 3A). The results of the MEME-ChIP analysis for the 500 strongest bindings demonstrated that one of the most common motifs contains three half nuclear receptor binding sites, which may form two overlapped DR1s (p = 7.52E-6) ( Figure 3B). The other identified sequence contains a GC box that matches to the motif of specific protein 1 (Sp1) binding site (p = 2.28E-5) ( Figure 3C) [43].

Functional analysis of RXRa based on genetic profiling of hepatocyte RXRa-deficient mouse livers
Biological functional analysis of both UIK and DIK genes demonstrated that RXRa-dependent genes participated predominantly in oxidation/reduction, lipid metabolism, generation of precursor metabolites and energy, cofactor metabolism, carboxylic acid biosynthesis, organic acid biosynthesis, coenzyme metabolism, protein folding, electron transport chain, translation and apoptosis. Among these, an average of 84% of the genes was bound by RXRa ( Table 2). The ratio of UIK and DIK genes in each pathway revealed a distinct role for hepatic RXRa in regulating lipid metabolism and apoptosis. For example, all of the sterol and cholesterol biosynthetic-related genes and most of the steroid and lipid biosynthetic-related genes were up-regulated in KO mouse livers, which indicate the role of hepatic RXRa in the metabolism of sterol and cholesterol. In the apoptosis pathway, 14 out of 20 anti-apoptotic genes were significantly down-regulated in KO livers, whereas only 6 were up-regulated. This finding suggests the pro-apoptotic role of hepatocyte RXRa. In addition to serving as a transcriptional factor, RXRa has a role in regulating translation as many translation-related genes were down-regulated due to RXRa deficiency.
There were four UIK genes, including spleen tyrosine kinase (Sykb), 3-hydroxy-3-methylglutaryl-Coenzyme A synthase 1 (Hmgcs1), isopentenyl-diphosphate delta isomerase (Idi1), and NADP dependent steroid dehydrogenase-like (Nsdhl), which were not bound by RXRa and are also involved in steroid and cholesterol metabolism, suggesting RXRa could directly and indirectly regulate steroid and cholesterol homeostasis ( Figure 4A). The analysis also revealed that the pyrin-related genes including Mediterranean fever (Medv), myeloid cell nuclear differentiation antigen (Mnda), and interferon activated genes 204 and 205 (Ifi204, Ifi205) were significantly induced due to lack of hepatic RXRa, but they do not have RXRa binding sites ( Figure 4A). In addition, there were 25 DIK genes without RXRa binding sites, which were functionally related to intracellular non-membranebounded organelle, ribonucleoprotein complex, ribosome, translation or structural molecule activity ( Figure 4B). Additionally, biological functions of RXRa-bound genes (Table 3) and RXRabound non-hepatic genes (Table 4) were analyzed. The RXRabound genes were implicated extensively in metabolic process, gene expression, protein processing, cell death, and response to stress, which encompassed most of the functional pathways of RXRa-dependent genes. However, the analysis of RXRa-bound non-hepatic genes revealed the potential role of RXRa in other tissues, such as neuron differentiation and development, cell adhesion, motion, and morphogenesis as well as signal transduction.
The interaction between RXRa and other nuclear receptors as well as cofactor genes Surprisingly, many steroid and orphan nuclear receptor genes were bound by RXRa in mouse livers (Table 5). Among these, the mRNA levels of the retinoid acid receptor beta (Rarb), nuclear receptor 1d1 (Nr1d1, Rev-erba), and Nr5a2 was significantly changed due to hepatic RXRa deficiency. These findings were validated by real-time PCR in livers ( Figure 5). Notably, some nuclear receptor cofactors, such as peroxisome proliferative activated receptor c coactivator 1 a (Ppargc1a), b (Ppargc1b), and nuclear receptor co-repressor 2 (Ncor2), were also bound by RXRa and Ppargc1a and Ncor2 was also expressed in an RXRadependent manner ( Figure 5).

Discussion
It has been more than two decades since the cloning of RXRa. Accumulated literature clearly indicates the importance of this nuclear receptor in regulating liver disease processes such as metabolic syndrome, alcoholic liver disease, chronic hepatitis C, and liver cancer, which comprise some of the most serious worldwide health issues today [2,4,5]. Thus, a global profiling of the bona fide RXRa genomic binding sites has become essential in order to illustrate the function and underlying regulatory mechanism mediated by RXRa. Using high throughput genomic methods and the knockout mouse model, the current study not only identifies the in vivo interaction of RXRa and the mouse genome, but also establishes the relationship between binding and hepatic gene expression as well as biological pathways. The hepatocyte RXRa-deficient mice were generated by deleting the DNA binding domain of the RXRa [9]. Thus, the differential gene expression between WT and KO mouse livers is due to lack of direct DNA binding of RXRa. Our data indicated that RXRa bound 78.7% of hepatic expressed genes and more than 80% of RXRa-dependent genes revealing the direct extensive role of RXRa in regulating liver gene expression and liver function. Additionally, 202 RXRa-binding free genes also displayed an RXRa-dependent expression pattern. These genes could be regulated by RXRa indirectly or post-transcriptionally. Moreover, RXRa binding sites were located in all of the mouse chromosomes as well as the mitochondrial DNA. However, the enriched bindings on chromosome 5 and 11 and rare bindings on Y chromosome suggest specific biological functions of RXRa. The current study also presented evidence for the potential role of RXRa in regulating mitochondrial gene function. Since mitochondria DNA is maternally derived, the finding is congruent with the specific binding of RXRa to the X chromosome. These observations suggest a gender-specific effect of RXRa. Taken together; the data indicate the extensive biological action of endogenous RXRa ligands such as retinoic acid and fatty acids in the liver. RXR has three isoforms. Among them, RXRa has the highest expression level in the liver. It would be important to study the differential and redundant role of RXRb and c using the same approach.
Promoter and enhancer regions are regarded as the most important regulatory regions in the control of transcription [44]. Comparing the distribution of RXRa binding sites in non-hepatic, hepatic, and RXRa-dependent genes, RXRa bindings occurred more frequently within the coding and promoter regions for hepatic than non-hepatic genes. However, the distribution profile of RXRa binding was similar between RXRa-dependent andindependent genes. Thus, it seems like direct binding of RXRa does not necessarily affect baseline mRNA levels. There are several potential reasons for this. Firstly, the effect of RXRs might be redundant and the presence of RXRb and c in hepatocyte RXRa-deficient mouse liver may be sufficient to maintain the basal transcriptional machinery of the hepatic RXR target genes. Secondly, albumin-cre recombinase was used to produce the RXRa knockout mice. Since albumin is predominantly expressed in the hepatocyte, it is also possible that some of those RXRabound genes are predominantly expressed in other types of liver cells rather than hepatocytes. Thus, RXRa deficiency does not alter the hepatic mRNA level. Thirdly, RXRa binding can be ''silent''. Active gene transcription occurs upon ligand binding, which leads to recruitment of many other factors and Pol II. It is possible that the biological level of RXRa ligands is not sufficient to induce gene transactivation and thus knockout RXRa has no impact on basal mRNA level. Additional experiments that use pharmacological ligands such as retinoic acids and polyunsaturated fatty acids to activate RXRa are needed to identify exogenous ligand/RXRa-dependent genes. The expression of gene is tissue specific. Figure 1A shows that there are 7642 genes have RXRa binding, but mRNA levels of those genes cannot be detected in the liver. Those genes might be expressed in the non-hepatic tissues such as neuron and skeletal muscle. We speculate that RXRa might have the most extensive binding in the liver genome since it has the highest expression level in the liver. Those binding sites may also exist in the DNA derived from other tissues. As a master nuclear receptor, RXRa plays an essential role in effecting the biological actions of other adopted orphan nuclear  receptors. Based on the ''1-2-3-4-5 rule'', RXR and its partners preferentially bind to specific motifs composed of hexamers (A/ GGGTCA) separated by a various number of nucleotides [45]. For example, DR1 is a preferred motif for peroxisome proliferator activator receptors (PPAR)/RXR heterodimer and RXR homodimer. DR4 is the preferred binding site for liver x receptor and thyroid hormone receptor [1,46], and IR1 is the preferred motif for farnesoid x receptor [47]. DR3, DR4, and ER6 are usually recognized by pregnane x receptor [48]. DR2 and 5 are the main binding motifs for retinoic acid receptor [1,45]. Therefore, the specific heterodimeric partners of RXRa that are involved in the regulatory process can be predicted based on the motif. The present data showed that DR1 was the most common motif, followed by DR4, IR1, and DR3. According to the ''1-2-3-4-5 rule'' mentioned above, PPAR/RXRa heterodimer and RXRa homodimer may be the dominant regulators followed by LXR, FXR, and PXR in the liver. In agreement with this finding, a recent study showed that 67.7% and 25.7% of RXRa binding sites overlap with those of PPARa and LXR, respectively, in normal mouse livers [49]. Furthermore, the finding is consistent with the phenotype of hepatocyte RXRa-deficient mice, which have elevated serum cholesterol and triglyceride levels [16]. In addition, lipid, bile acid, and xenobiotic metabolism are altered when the hepatocyte RXRa-deficient mice are challenged with exogenous ligands for those receptors [10,12]. Notably, analysis of common sequences of the 500 strongest RXRa peaks demonstrated that GC box, the motif of Sp1, may be one of the most popular motifs for RXRa binding. Consistently, it has been shown that there is a thyroid hormone response element overlapping with the GC box in the promoter of the epidermal growth factor receptor, where it is bound competitively by RXR/T3R and Sp1 [42]. Furthermore, a physical interaction between RAR/RXR and Sp1 has also been reported and such interaction synergistically enhances the expression of RA-induced genes in vitro [50]. Our findings suggest the presence of extensive crosstalk between Sp1 and RXRa and their collaborative or competitive regulatory role in regulating liver gene transcription. Recent studies have reported the mouse liver genome-wide binding profile of FXR [26], PXR [33], LXR, and PPARa [49]. Comparing the RXRa ChIP-seq data from this study with others' findings, common bindings were aligned between RXRa and its partners (data not shown). Accordingly, biological functional analysis of RXRa-bound genes (Table 3) encompassed almost all its partners' biological functions. In addition, we identified novel pathways that are specific for RXRa, which include oxidation/ reduction, protein localization, intracellular signaling cascade, regulation of transcription, cofactor metabolic process, cellular response to stress, and cell death. Furthermore, novel potential of RXRa was also unveiled by the analysis of RXRa-bound nonhepatic genes. Based on binding and expression profiling generated from wild type and knockout mice, this paper is the first to establish the relationship between binding and expression in an RXRa-dependent manner.
It is important to note that our results did not show any significant correlation between peak score and the fold change in mRNA level caused by RXRa deficiency. In addition, the binding characteristic (location of the peak) is similar between UIK and DIK genes. This finding suggests that physical interaction of RXRa is essential, but not sufficient for predicting the subsequent transcriptional effect.
The role of retinoids in other organs such as the eye and skin is well known. However, retinoids are stored, processed, and metabolized in the liver. In addition, the liver is also a retinoid target organ and yet the action of retinoids in the liver has been Figure 5. Expression of genes in wild type and RXRa-null livers. RNA extracted from wild type and RXRa-null livers (n = 3-4) were subjected to real time PCR to determine the expression level of the studied genes. Data were normalized to glyceraldehyde 3-phosphate dehydrogenase (GAPDH) mRNA level.**: p,0.05 in comparisons between two groups. Microarray experiment showed that fold change in Rarb, Nr1d1, Nr5a2, Ncor2, and Ppargc1a mRNA levels due to hepatic RXRa deficiency was 1.8, 0.7, 0.5, 1.3, and 0.3, respectively (n = 3, p,0.05). doi:10.1371/journal.pone.0050013.g005 overlooked. This study, for the first time, demonstrates the potential biological effect of endogenous ligands of RXRa in the liver. Taken together, as an active partner of many nuclear receptors, our reported data showed that RXRa and its endogenous ligands control liver metabolism and function in general.