Concerted effects of heterogeneous nuclear ribonucleoprotein C1/C2 to control vitamin D-directed gene transcription and RNA splicing in human bone cells.

Traditionally recognized as an RNA splicing regulator, heterogeneous nuclear ribonucleoprotein C1/C2 (hnRNPC1/C2) can also bind to double-stranded DNA and function in trans as a vitamin D response element (VDRE)-binding protein. As such, hnRNPC1/C2 may couple transcription induced by the active form of vitamin D, 1,25-dihydroxyvitamin D (1,25(OH)2D) with subsequent RNA splicing. In MG63 osteoblastic cells, increased expression of the 1,25(OH)2D target gene CYP24A1 involved immunoprecipitation of hnRNPC1/C2 with CYP24A1 chromatin and RNA. Knockdown of hnRNPC1/C2 suppressed expression of CYP24A1, but also increased expression of an exon 10-skipped CYP24A1 splice variant; in a minigene model the latter was attenuated by a functional VDRE in the CYP24A1 promoter. In genome-wide analyses, knockdown of hnRNPC1/C2 resulted in 3500 differentially expressed genes and 2232 differentially spliced genes, with significant commonality between groups. 1,25(OH)2D induced 324 differentially expressed genes, with 187 also observed following hnRNPC1/C2 knockdown, and a further 168 unique to hnRNPC1/C2 knockdown. However, 1,25(OH)2D induced only 10 differentially spliced genes, with no overlap with differentially expressed genes. These data indicate that hnRNPC1/C2 binds to both DNA and RNA and influences both gene expression and RNA splicing, but these actions do not appear to be linked through 1,25(OH)2D-mediated induction of transcription.

INTRODUCTION

Steroid hormones regulate gene transcription via nuclear receptors (NR) that bind to specific DNA response elements around and within target genes. Transcriptional regulation by NRs involves recruitment of coregulator proteins that coordinate interactions with the basic transcriptional machinery (1). NR coregulatory proteins have structural similarities to RNA splicing proteins, suggesting that some of these coregulators could act as coupling factors linking transcription and RNA splicing (2). This dual role for NR and their coregulators appears to be part of the process that maintains both the magnitude and correctness of the gene products generated in response to steroid hormone signals (3). Given the importance of this facet of gene expression, a reciprocal arrangement by which components of the splicing machinery act as NR transcriptional coregulators seems likely.

A link between RNA splicing and the regulation of transcription has been described for the stable small nuclear ribonucleoproteins that make up the core of the spliceosome (4,5) and the splicing factors that facilitate formation of the spliceosome and RNA splicing itself (6,7). Prominent amongst this latter group is the heterogeneous nuclear ribonucleoprotein (hnRNP) family of RNPs that associate with pre-mRNA (heterogeneous nuclear RNA) and have multiple functions in RNA processing (8). Owing to their ability bind specifically to both single- (9–11) and double-stranded (11–14) DNA elements, hnRNPs may also regulate the products of transcription by alternative means (15). In previous studies we have shown that hnRNPs can modulate transcriptional responses to the active form of vitamin D, 1,25-dihydroxyvitamin D (1,25(OH)2D), by competing with the nuclear vitamin D receptor (VDR) for binding to a vitamin D response elements (VDRE) that control target gene expression (16). Studies using cells from subhuman primates with a physiologically relevant state of relative insensitivity to the actions of 1,25(OH)2D (17) identified an hnRNP-like VDRE-binding protein (VDRE-BP) (18,19) that attenuated 1,25(OH)2D-mediated transcription by direct competition with the 1,25(OH)2D-activated VDR for a VDRE (18).

Subsequent studies indicated that the human VDRE-BP was hnRNPC1/C2 (20), with hnRNPC1/C2-induced vitamin D-insensitivity postulated as a cause of hereditary vitamin D resistant rickets (HVDRR) in a patient in the absence of a mutation in the VDR gene (19). In this particular case, dysregulation of 1,25(OH)2D-induced transcription appears to be due to disruption of the normal cyclical patterns of VDRE occupancy by direct competition between VDR and hnRNPC1/C2 for DNA binding (20). To date functional analysis of hnRNPC1/C2-VDR interactions has focused on feedback suppression of expression of the CYP24A1, the gene encoding the 24-hydroxylase enzyme which promotes the catabolism of 1,25(OH)2D (20). CYP24A1 is known to be alternatively spliced, to generate in-frame mRNA variants that include skipping of exons 1 and 2 (CYP24A1-SV) (21), and exon 10 (CYP24A1 variant 2) (22). However, studies using cells from the HVDRR patient described previously suggest that the transcriptional interaction between VDR and hnRNPC1/C2 involves many more genes (23). In the current study we have utilized CYP24A1,as a model target gene in combination with genome wide analyses to assess the role of hnRNPC1/C2 in coordinating transcriptional and RNA splicing responses to 1,25(OH)2D.

MATERIALS AND METHODS

Reagents and cell culture

Human MG63 osteosarcoma cells (American Type Culture Collection, Manassas, VA, USA) were used as a model for vitamin D signaling. Previous studies have shown that these cells demonstrate vitamin D-induced gene expression responses common to normal human osteoblastic cells such as osteocalcin and alkaline phosphatase (24,25). MG63 cells were cultured in Dulbecco's modified Eagle medium (Invitrogen) with 10% Fetal Bovine Serum at 37°C with 5% CO2. Crystalline 1,25-dihydroxyvitamin D3 (1,25(OH)2D3, Enzo Life Sciences, Farmingdale, NY, USA) was reconstituted in ethanol. Ethanol (0.1%) was used as vehicle treatment.

Analysis of 25OHD metabolism by MG63 cells

Activity of the 24-hydroxylase enzyme in MG63 cells was quantified using tritiated 25-hydroxyvitamin D (3H-25OHD, PerkinElmer, Waltham, MA, USA) substrate as described previously (26). 3H-25OHD and 3H-24,25-dihydroxyvitamin D (3H-24,25(OH)2D) were quantified using High Performance Liquid Chromatography. Data were reported as fmoles of vitamin D metabolites produced/hr/μg cellular protein from n = 3 separate cell preparations.

Western blot analysis

Cell lysates were prepared in RIPA buffer (150 mM sodium chloride, 1.0% Triton X-100, 0.5% sodium deoxycholate, 0.1% sodium dodecyl sulphate (SDS), 50 mM Tris, pH 8.0) with 1x protease inhibitor cocktail (Thermo, Rockford, lL, USA). Protein samples were separated by SDS-PAGE. Antibodies included: hnRNP C1/C2 (Santa Cruz Biotechnology, # sc-32308); CYP24A1 (Santa Cruz Biotechnology, # sc-32166); β-actin (Santa Cruz Biotechnology, # sc-81178).

Quantitative RT-PCR (qRT-PCR) analyses

Total RNA was extracted from MG63 cells and reverse transcribed using RNeasy mini kit (Qiagen, Valencia, CA, USA) and SuperScript III Reverse Transcriptase (Invitrogen) as described previously (23). Polymerase chain reaction (PCR) amplification of cDNA was carried out using TaqMan system reagents (Applied Biosystems) as described previously (23), under following conditions: (i) 95°C for 10 min; 40 cycles of [95°C for 15 s, 60°C for 1 min] or (ii) 95°C for 10 min, 40 cycles of [95°C for 15 s, 55°C for 30 s, 72°C for 30 s] and then finally 72°C for 2 min using Maxima SYBR Green qPCR Master Mixes (Thermo). Target genes were normalized to 18S (Primers from Applied biosystems, Carlsbad, CA, USA) or GAPDH expression. The primers sequences used are described in the Supplementary Table S1. All reactions were performed in triplicate for each biological experiment.

Chromatin immunoprecipitation-qPCR (ChIP-qPCR)

Chromatin immunoprecipitation (ChIP) was performed as described previously (20) with several modifications. Briefly, following treatment with 10 nM 1,25(OH)2D or vehicle for 3 h, cells were washed twice with phosphate-buffered saline (PBS) and crosslinked with 1% formaldehyde at room temperature for 10 min. After quenching of cross-linking with glycine at final concentration of 125 mM, cell monolayers were washed twice with cold PBS and then scraped into 5 ml cold PBS. Collected cells were pelleted by centrifugation for 5 min at 1000 g and the resulting pellets resuspended in 1 ml of ChIP lysis buffer (5 mM Pipes, pH 8.0, 85 mM KCl, 0.5% Nonidet P-40, 1 mM dithiothreitol, 0.25 mM phenylmethylsulfonyl fluoride and 1 g/ml each of pepstatin, leupeptin and aprotinin) with protease inhibitor cocktail (Thermo). The resulting chromatin lysates were sonicated to yield sheared DNA fragments of sizes between 300 and 1000 bp, and then centrifuged for 10 min, 8000 g at 4°C. Supernatants were removed and aliquots used as ‘Input’ sample. For each immunoprecipitation, sheared chromatin was diluted 1:10 with ChIP dilution buffer (1% Triton X-100, 2 mM EDTA, 150 mM NaCl and 16.7 mM Tris-HCl, pH 8.1) with 1x protease inhibitor cocktail (Thermo). The chromatin was incubated at 4°C overnight with 50 μl of magnetic Dynabeads protein A (Life Technologies AS, Oslo, Norway) pre-conjugated with antibodies. The bead–antibody–antigen complexes were then subjected to serial washes of 10 min each with 2 × 1 ml ChIP-RIPA (1 mM EDTA, 0.7% sodium deoxycholate, 1% NP-40, 0.5 M LiCl, 50 mM HEPES, pH 7.6), 2 × 1 ml ChIP-RIPA buffer and 2 × 1 ml TE buffer (1 mM Tris, pH 7.6, 1 mM EDTA). The resulting immune complexes were then resuspended in 100 μl reversal dilution buffer (1% SDS, 0.1 M NaHCO3), and vortexed for 30 s every 5 min for 30 min. Formaldehyde cross-linking was reversed by heating at 70°C 2 h. DNA was extracted using QIAquick Gel Extraction Kit (Qiagen) and analyzed by qRT-PCR. ChIP grade antibodies against VDR (Santa Cruz Biotechnology, Santa Cruz, CA, USA; sc-13133), hnRNPC1/2 (Santa Cruz Biotechnology, # sc-32308) were used. Normal mouse or rabbit IgG (Santa Cruz Biotechnology, # sc-2025, sc-2027) was used as the negative control. qRT-PCR was performed as described above with primer set flanking enhancer VDREs in the proximal CYP24A1 (VDRE1 [−150 bp] CGCCCT- [CAC]TCACCT [−136 bp] and VDRE2 is [−258 bp] CGCACC[CGC]TGAACC [−244 bp]). Primers were designed to amplify a region of DNA from −304 to −63 bp relative to the start site of transcription (Supplementary Table S1). A separate primer set spanning a region not harboring a VDRE in calponin 1 (CNN1) was used as a negative control (Supplementary Table S1). The amplified products were normalized to Input DNA content and presented as percent of Input. ChIP experiments were performed in three biological independent replicates.

RNA immunoprecipitation (RIP) from cell lysates

RNA immunoprecipitation (RIP) was performed using MILLIPORE Magna RIP™ RNA-Binding Protein Immunoprecipitation Kit (Millipore, # 17–700) according to manufacturer's instructions. Antibodies to hnRNPC1/C2 (Santa Cruz Biotechnology, # sc-32308), VDR (Santa Cruz Biotechnology, # sc-13133), normal mouse or rabbit IgG (Santa Cruz Biotechnology, # sc-2025, sc-2027) were used. Anti-SNRNP70 (Millipore, Temecula, CA, USA) was a positive control antibody. RNA extraction and qPCR were as described above.

Small interfering RNA (siRNA) knockdown of hnRNPC1/C2

Two different human hnRNPC1/C2 stealth select small interfering RNA (siRNA) (Invitrogen, HSS179304 and HSS17930) and a stealth siRNA negative control low GC (Invitrogen) were used. Knockdown of hnRNPC1/C2 was performed as described previously (27). Briefly, MG63 cells were plated at 105 cells/well in a 6-well cluster plate in 2 ml growth medium and 0.03 × 106 cells/well in a 24-well plate in 0.5 ml of growth medium. Cells were transfected at 20–30% confluency with hnRNPC or control siRNA at a final concentration of 5 nM, using Lipofectamine RNAiMAX (Invitrogen) transfection reagent according to the manufacturer's protocol. Ninety six hours after transfection, cells were harvested and knockdown efficiency assessed by qRT-PCR and Western blot. All transfection experiments were performed in triplicate. Transfection efficiency of siRNAs in MG63 cells was determined by using Block-iT Alexa Fluor Red Fluorescent Control Oligo (Invitrogen, cat# 14750-100). Optimization of siRNA knockdown of hnRNPC1/C2 is shown in Supplementary Figure S1.

Fluorescently labeled PCR

Fluorescently labeled PCR to visualize endogenous, alternatively spliced CYP24A1 exon 10 products was modified from previously described methods (28). Briefly, a 22 nt universal tag sequence (designated ‘GFPN’), 5′-CGTCGCCGTCCAGCTCGACCAG-3′ derived from GFP N-terminal region, was added to the 5′ end of the forward primer during oligo synthesis, while the reverse primer remained untagged. A fluorescently labeled universal primer 5′-FAM-CGTCGCCGTCCAGCTCGACCAG-3′ was added as a third primer in the PCR reaction, under the following amplification conditions: (i) 98°C for 10 s, 32 cycles of [98°C for 30 s, 50°C for 30 s, 72°C for 1 min] and then (ii) 72°C for an additional 10 min. Primer sets are shown in Supplementary Table S1. Reaction products were resolved on 5% TBE-urea polyacrylamide gels and visualized using a BIO-RAD ChemiDocTM MP Gel Imaging system (Bio-Rad, 170–8280). Bands were quantified using Image lab software (Bio-Rad). Exon inclusion levels were calculated as the intensity of the exon inclusion band divided by the total intensity of the exon inclusion and skipping bands.

Analysis of alternative RNA splicing using CYP24A1 minigene constructs

CYP24A1 mini-gene splicing cassettes were cloned using genomic DNA isolated from MG63 cells with the Qiagen DNeasy kit. A DNA fragment containing the 3′ portion of exon 9, intron 9, exon 10, intron 10 and the 5′ portion of exon 11 of CYP24A1 was synthesized with Q5 Hot Start High-Fidelity 2x master mix using a forward primer containing Not I, 5′ ATATATGCGGCCGCGGCAACAGTTCTGGGTGAAT 3′ and a reverse primer containing Apa I 5′ ATATATGGGCCCATTGTCTGTGGCCTGGATGT 3′. A pcDNA 3.1(+) and PCR fragment was sequentially digested with Not I and Apa I and purified with the Qiagen DNA gel extraction kit for isolation by agarose gel electrophoresis. Ligation of DNA fragments was accomplished with NEB T7 ligase, and library efficiency DH5a competent cells (Thermo) were transformed and screened.

Genomic DNA from MG63 was also used as template to generate PCR fragments containing specific proximal promoters with Q5 Hot Start High-Fidelity 2x master mix. The PCR primer sets used to generate CYP24A1 (−1197/+87 bp), interleukin-6 (IL-6) (−1183/+31 bp), β-actin (ACTB) (−472/+49 bp), Glyceraldehyde-3-Phosphate Dehydrogenase (GAPDH) (−487/+20) and osteocalcin (BGLAP) (−834/+46 bp) gene fragments are shown in Supplementary Table S2. The CMV promotor of pcDNA3.1(+) was removed by sequential digestion with Nru I and Kpn I and replaced by ligation of specific proximal promoters described above with forward primers containing the Nru I site and reverse primers containing the Kpn I site. An additional CYP24A1 promoter fragment was generated containing mutations in the VDRE introduced using the Agilent QuikChange Site-Directed Mutagenesis Kit with primers designed via Agilent QuikChange primer design program (http://www.genomics.agilent.com/primerDesignProgram.jsp). VDRE-1: AGGTGAGCGAGGGCG (wild-type), TCCTGAGCGTCCGCG (mutated). VDRE-2: AGTTCACCGGGTGTG (wild-type), TCCTCACCGTCCGTG (mutated). Half-sites are shown in bold and underlined sequences show introduced mutations (29).

For transfection of minigene constructs, MG63 cells were cultured in 12-well plates and transfected at 30% confluency with Lipofectamine 3000 (Invitrogen) according to manufacturer protocols. Forty eight hours post-transfection, RNA was isolated by Trizol followed by in-column DNAse treatment according to Qiagen RNeasy protocol. cDNA was prepared from 500 ng RNA by SuperScript III reverse transcriptase. Minigene splice products were generated by PCR using primer sets shown in Supplementary Table S3 with Q5 Hot Start High-Fidelity 2x master mix using the following amplification conditions: 98°C for 30 s; 25 to 36 cycles of 98°C for 10 s, 60°C for 20 s, 72°C for 1 min. PCR products were analyzed by agarose electrophoresis with ethidium bromide; band intensity was determined by ImageJ verision 1.49 and exon skipping percentage was calculated.

High-throughput RNA sequencing (RNAseq) and genomic mapping

RNA extracted was extracted from MG63 cells using RNAeasy mini kit (Qiagen, Valencia, CA, USA) with on-column DNase treatment to remove genomic DNA. cDNA libraries were prepared using the Illumina TruSeq RNA Sample Preparation Kit (illumina). High-throughput sequencing was performed using an Illumina HiSeq2500 (paired-end, non-strand-specific 107-bp read length). Knockdown and control samples were sequenced together in two flow cells on four lanes. Reads from each of four lanes were combined for all analyses. A total of 525 million 107-bp paired-end reads were generated by sequencing (Supplementary Table S4). Reads were mapped to human transcripts (Ensembl, release 72) and genome (hg19) using Tophat software (v1.4.1) (30) allowing up to 3 bp mismatches per read and up to 2 bp mismatches per 25 bp seed. More than 60% of all RNAseq reads mapped uniquely to the human genome (Supplementary Table S4). RNAseq data were deposited with the NCBI Gene Expression Omnibus <>.

http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?token=qtyzgcaifzkxvev&acc=GSE63086.

Analyses of differential gene expression/alternative splicing

Cuffdiff (v.2.2.0) (31) was used to calculate gene expression levels with the FPKM metric (fragments per kilobase of exon per million fragments mapped) and differential gene expression (DEG) between two conditions assessed at an false detection rate (FDR) < 1% and minFPKM > 0.1. rMATS (v3.0.8)] (30–33) (http://rnaseq-mats.sourceforge.net/) was used to identify differential alternative splicing (AS) events corresponding to five major types of AS patterns. AS events were identified as those with significant difference in inclusion levels (|ΔPSI| ≥ 0.05) between two groups, FDR < 0.05. Stringency for AS events was defined for varying PSI and compared with different stringency for DEG (Supplementary Figure S2) to provide comparable DEG and differentially splice genes for each cell condition.

RT-PCR validation of differentially spliced exons identified by RNA-seq

RNA extraction, reverse transcription and PCR amplification were as described previously (23). Sequences for primer pairs used are listed in Supplementary Table S1.

Data normalization and statistical analysis

Data for ChIP-qPCR assays were reported as % input, and normalized as described previously (23). Briefly, each ChIP DNA fraction Ct value was normalized to the input DNA fraction to account for chromatin preparation differences: △Ct[normalized ChIP] = {Ct[ChIP] - [Ct[input] - log2 (input dilution factor)]}, where the dilution factor was 100. Triplicate normalized ChIP Ct values were meaned. To calculate % Input for each ChIP fraction, linear conversion of normalized ChIP △Ct: % Input = 2 (-ΔCt [normalized ChIP]. Data for RIP-qPCR were presented and normalized as for ChIP-qPCR. Data were compared statistically using an unpaired Student's t-test and one way analysis of variance (ANOVA) with Tukey's multiple-comparison post-hoc test as indicated.

RESULTS

1,25(OH)2D-induced 24-hydroxylase activity is associated with increased expression of CYP24A1 and CYP24A1-splice variants

MG63 cells treated with 1,25(OH)2D showed increased 24-hydroxylase activity (conversion of 3H-25OHD to 3H-24,25-dihydroxyvitamin D (3H-24,25(OH)2D); Figure 1A). This was associated with increased expression of the full length CYP24A1-24-hydroxylase protein (Figure 1C). By contrast, 1,25(OH)2D had no effect on expression of mRNA (Figure 1B) or protein (Figure 1C) for hnRNPC1/C2. Time course qRT-PCR analyses showed that 1,25(OH)2D-induced 24-hydroxylase activity was associated with increased expression of full length mRNA for CYP24A1 (Figure 1D) and two splice variants of the CYP24A1 gene, CYP24A1-variant 2 (Figure 1E) and CYP24A1-splice variant (CYP24A1-SV) (Figure 1F).

Figure 1.

Effect of 1,25(OH)2D on 24-hydroxylase and hnRNPC1/C2 in MG63 cells. (A) (Left) 24-hydroxylase enzyme activity (conversion of 3H-25OHD to 3H-24,25(OH)2D) in MG63 cells 24 h after addition of vehicle or 1,25(OH)2D (10 nM). Data are mean ± SD fmol 3H-24,25(OH)2D/hr/μg cell protein (n = 3, **P < 0.01, Student's t-test). (Right) Representative high performance liquid chromatography analyses for vehicle (top) and 1,25(OH)2D (bottom) treatment. (B) qRT-PCR analysis of hnRNPC1/C2 mRNA expression after treatment with vehicle or 1,25(OH)2D (10 nM) for 12 h. Data are mean ± SD, n = 3. n.s., non-significant. (C) Western blot analysis of 24-hydroxylase and hnRNPC1/C2 with vehicle or 1,25(OH)2D (10nM) 12 h; N.T., No Treatment. (D) Expression of mRNA for full length CYP24A1, E. CYP24A1-variant 2 (exon 10 skipped), F. CYP24A1-SV (exon 1 and 2 skipped) following treatment with 1,25(OH)2D (10 nM) up to 12 h. Data are mean ± SD (n = 3, P < 0.05, 1-way ANOVA). Lower panels D and F show schematic illustrations of the full-length, wild-type CYP24A1 and the exon organization of the splice variants described above.

hnRNPC1/C2 binds to DNA and pre-mRNA for CYP24A1

ChIP-qPCR analyses showed that in the absence of 1,25(OH)2D, there was low level binding of both VDR and hnRNPC1/C2 to a region of the CYP24A1 gene promoter (−304 to −63 bp upstream of the transcriptional start site for CYP24A1) containing a VDRE (Figure 2A). However, in contrast to hnRNPC1/C2, binding of VDR to the CYP24A1 proximal promoter was increased following treatment with 1,25(OH)2D (Figure 2A). Parallel RNA immunoprecipitation (RIP)-qPCR was carried out using primers targeted at pre-mRNA, intronless mature mRNA or both pre- and mature mRNA for CYP24A1 (Figure 2B–D). Pre-mRNA for CYP24A1 showed significantly higher binding of hnRNPC1/C2 following treatment with 1,25(OH)2D. For the mature CYP24A1 mRNA, hnRNPC1/C2 showed low level association that was unaffected by 1,25(OH)2D (Figure 2C). RIP with antibodies for VDR showed no association of the receptor with pre- or mature RNA for CYP24A1 (data not shown).

Figure 2.

CYP24A1 DNA and RNA binding of hnRNPC1/C2. (A) Chromatin immunoprecipitation (ChIP) using antibodies to hnRNPC1/C2 (C1/C2) or vitamin D receptor (VDR). Precipitated DNA was PCR amplified with primers to a VDRE-containing region of CYP24A1 promoter and a VDRE-negative region of calponin 1 (CNN1) (negative control). Data are % of input DNA. (B–D) RNA immunoprecipitation (RIP) analyses using antibodies to hnRNPC1/C2 (IgG as a control). Precipitated RNA was PCR amplified with primers encompassing CYP24A1 pre-mRNA, mature CYP24A1 RNA or both as indicated. Data are % of input RNA. Cells were treated with vehicle or 1,25(OH)2D (10 nM) for 3 h. Data are mean ± SD (n = 3, *P < 0.05, ***P < 0.001, Student's t-test, #, undetectable). n.s. = non-significant.

Effects of hnRNPC1/C2 knockdown on expression of CYP24A1 and splice variants

To assess the functional impact of hnRNPC1/C2 on vitamin D-mediated transcription and pre-mRNA splicing, siRNA was used to knockdown hnRNPC1/C2 in MG63 cells (Figure 3A). RT-PCR analyses showed divergent effects of hnRNPC1/C2 knockdown on full length CYP24A1 mRNA and splice variants (Figure 3B–D). In vehicle-treated cells, CYP24A1 mRNA was decreased and CYP24A1-variant 2 increased with hnRNPC1/C2 knockdown (Figure 3B and C). Knockdown of hnRNPC1/C2 had no effect on CYP24A1 expression when performed in conjunction with exposure to 1,25(OH)2D. Expression of CYP24A1-variant 2 was induced by 1,25(OH)2D in a fashion similar to that for full length CYP24A1, but this effect was significantly enhanced with hnRNPC1/C2 knockdown (Figure 3C); induction of CYP24A1-variant 2 was 62-fold ± 3, and 72 ± 8 in KD1 and KD2 cells versus 45 ± 3 in control cells. CYP24A1-SV was unaffected by hnRNPC1/C2 knockdown in either vehicle or 1,25(OH)2D-treated cells (Figure 3D).

Figure 3.

Effect of HnRNPC1/C2 knockdown on CYP24A1 transcription and pre-mRNA splicing. MG63 cells were transfected with two hnRNPC1/C2 siRNAs (KD1 and KD2) or scrambled control siRNA (Ctrl). (A) Western blot analysis of CYP24A1 24-hydroxylase and hnRNPC1/C2. Cells were treated with vehicle or 1,25(OH)2D (10 nM) for 12 h. (B–D) qRT-PCR analysis of mRNA for: (B) wildtype CYP24A1; (C) CYP24A1-variant 2; (D) CYP24A1-SV in KD1, KD2 and Ctrl cells. (E) CYP24A1 exon 10 inclusion by fluorescently labeled RT-PCR. (Left) Inclusion (gray) and exclusion (white) of exon 10. (Right) Gel visualization of PCR products with exon 10 inclusion or skipping. For D–G, cells were treated with vehicle or 1,25(OH)2D (10 nM) for 12 h. Data are mean ± SD (n = 3, *P < 0.05, **P < 0.01, Student's t-test).

To further characterize generation of the CYP24A1-variant 2 in MG63 cells, CYP24A1 exon 10-inclusion was quantified using fluorescently tagged RT-PCR (Figure 3E). Vehicle-treated cells showed 96.9% inclusion of exon 10, which was reduced to 93.9% and 89.0% after knockdown of hnRNPC1/C2. In 1,25(OH)2D-treated cells, inclusion of exon 10 was similar to vehicle (97.8%). However, in 1,25(OH)2D-treated cells with hnRNPC1/C2 knockdown, inclusion of exon 10 decreased to 73.7% and 75.7%. An alternatively spliced product (lacking exon 10) was present in each of the hnRNPC1/C2 knockdown treatments, but the PCR band corresponding to this product was more pronounced in 1,25(OH)2D-treated cells (Figure 3E).

VDRE expression and mRNA splicing in a CYP24A1 mini-gene model

To assess a possible link between VDRE and mRNA splicing for 1,25D-VDR-target genes, a CYP24A1 mini-gene expression construct was used to mimic alternative splicing associated with CYP24A1-variant 2 (exons 9–11). In the presence of a conventional CMV promoter or other promoters lacking a VDRE (IL-6, ACTB and GAPDH), exon skipping to generate a transcript lacking exon 10 ranged from 15–22% of correctly spliced transcript (Figure 4). By contrast, minigene constructs containing functional VDRE(s) (CYP24A1 and BGLAP) showed less exon 10 skipping (5.37 ± 0.67% and 7.00 ± 1.91%, respectively). The enhanced splicing fidelity associated with the CYP24A1 bearing functional VDREs was diminished by mutation of specific nucleotides associated with VDR-VDRE interaction (10.03 ± 2.50%).

Figure 4.

Alternative splicing of CYP24A1 is modulated by VDRE. Splicing of minigene constructs transfected into MG63 cells containing DNA for exons 9–11 of CYP24A1 linked to promoter DNA from: CMV; CYP24A1 containing two wild-type VDREs; CYP24A1 bearing mutated, non-VDR binding VDREs; VDRE-negative IL-6, ACTB and GAPDH promoters; and VDRE-positive BGLAP. (A) Example of a gel showing bands for correctly spliced minigene transcripts containing exons 9, 10 and 11 (upper band) and alternatively spliced transcripts containing only exons 9 and 11 (skipped exon 10, lower band). (B) Percent exon 10 skipping for each minigene-promoter construct (mean ± S.D., n = 3 separate transfections). * = statistically different from CYP24A1, P < 0.05. Constructs with functional VDRE are highlighted as boxes.

Genome-wide analysis of the effects of hnRNPC1/C2 knockdown on transcription and mRNA splicing

To define the genome-wide impact of hnRNPC1/C2 on 1,25(OH)2D-induced transcription and RNA splicing, high throughput RNA-seq analysis of MG63 cells was carried out. Knockdown of hnRNPC1/C2 resulted in 3500 differentially expressed genes relative to vehicle-treated control cells (Supplementary Data S1), whilst treatment with 1,25(OH)2D resulted in 324 differentially expressed genes relative to vehicle control cells (Supplementary Data S2), Comparing these two data sets, 153 differentially expressed genes, including CYP24A1, were common to hnRNPC1/C2 knockdown and 1,25(OH2D-treatment (Figure 5A). Of the 171 differentially expressed genes specific to treatment with 1,25(OH)2D, 100 were also observed when comparing hnRNPC1/C2 knockdown with combined 1,25(OH)2D and hnRNPC1/C2 knockdown, whilst a further 87 differentially expressed genes were only induced by 1,25(OH)2D in hnRNPC1/C2 knockdown cells (Figure 5A) (Supplementary Data S3). Cells with combined hnRNPC1/C2 knockdown and 1,25(OH)2D treatment also showed a high number of differentially expressed genes (3844) relative to vehicle control cells (Supplementary Data S4). Of these, 2932 differentially expressed genes (76.3%) were common to hnRNPC1/C2 knockdown alone, with 912 being specific to treatment with 1,25(OH)2D in combination with hnRNPC1/C2 knockdown (Figure 5B). The distinct nature of these differentially expressed genes is also illustrated by heatmap representation of the differentially expressed genes within the various sectors of the Venn diagram comparisons shown in Figure 5A. The heatmap comparison for hnRNPC1/C2 knockdown versus 1,25(OH)2D and hnRNPC1/C2 knockdown (KD versus 1,25+KD) shows the strongest contrast between these treatments (see Supplementary Figure S3).

Figure 5.

Differential gene expression by 1,25(OH)2D or hnRNPC1/2 knockdown. Weighted Venn diagrams depicting overlap of differentially expressed genes: (A) Control versus hnRNPC1/C2 knockdown (Ctrl versus KD; maize), control versus 1,25 (10 nM, 6 h) (Ctrl versus 1,25; grey) or hnRNPC1/C2 knockdown vs 1,25(OH)2D combined with hnRNPC1/C2 knockdown (KD versus 1,25+KD; blue). (B) Control versus hnRNPC1/C2 knockdown (Ctrl versus KD, maize) or combined 1,25(OH)2D and hnRNPC1/C2 knockdown (Ctrl versus 1,25+KD, green). (C) Up- and down-regulated differentially expressed genes: (left) Control versus 1,25(OH)2D (Ctrl versus 1,25) or hnRNPC1/C2 knockdown versus 1,25(OH)2D combined hnRNPC1/C2 knockdown (KD versus 1,25+KD); (right) Control versus hnRNPC1/C2 knockdown (Ctrl versus KD). (D) Weighted Venn diagrams depicting overlap of differentially expressed genes between Ctrl versus 1,25(OH)2D and KD versus 1,25(OH)2D and hnRNPC1/C2 knockdown. Numbers of differentially expressed genes are shown within each segment. Heatmaps show patterns of down-regulated (blue) differentially expressed genes (upper heatmaps) and up-regulated (red) differentially expressed genes (lower heatmaps) for each treatment/knockdown comparison group and the overlapping differentially expressed genes common to those comparisons. * includes CYP24A1.

Differentially expressed genes induced by knockdown of hnRNPC1/C2 alone were 57% down-regulated (including CYP24A1) and 43% up-regulated (Figure 5C). By contrast, treatment with 1,25(OH)2D alone favored up-regulation of gene expression (70% of differentially expressed genes) (Figure 5C). This included induction of the CYP24A1 gene (Supplementary Figure S4), a response that was also observed when comparing hnRNPC1/C2 knockdown with the combination of hnRNPC1/C2 knockdown and 1,25(OH)2D (Figure 5C). In total 322 genes were upregulated following treatment with 1,25(OH)2D. Of these, 25% (82) were induced by 1,25(OH)2D alone, whilst 30% (95) were only induced by a combination of 1,25(OH)2D treatment and hnRNPC1/C2 knockdown (Figure 5D). The remaining 45% of genes (145/322) were common to treatment with 1,25(OH)2D or 1,25(OH)2D with hnRNPC1/C2 knockdown; importantly, this group included CYP24A1 (Figure 5D). Of the 171 genes down-regulated by 1,25(OH)2D, 33% (56) were only suppressed by 1,25(OH)2D alone, while 43% (74) were only suppressed by the combination of 1,25(OH)2D treatment with hnRNPC1/C2 knockdown (Figure 5D); the remaining 24% (41/171) were common to treatment with 1,25(OH)2D or combination 1,25(OH)2D treatment with hnRNPC1/C2 knockdown (Figure 5D).

Alternative splicing was assessed using a statistical algorithm for RNAseq, multivariate analysis of transcript splicing (MATs) (32,34), that detects differentially spliced genes corresponding to the five major modes of types of AS described in metazoan organisms (Figure 6A and Supplementary Data S5). Using this program, AS was analyzed following knockdown of hnRNPC1/C2, with or without 1,25(OH)2D. In MG63 cells 3921 differentially spliced genes (FDR < 0.05, |ΔPSI| ≥ 0.05) were observed following hnRNPC1/C2 knockdown, with the majority of these (2713) being skipped exons (Figure 6B). Differentially spliced genes induced by hnRNPC1/C2 knockdown in MG63 cells were compared with previously reported differentially spliced genes in hnRNPC1/C2 knockdown HeLa cells (27); 542 differential splicing events being common to MG63 and HeLa cells (Supplementary Figure S5). RT-PCR confirmed differentially spliced genes following hnRNPC1/C2 knockdown in 12 transcripts identified by MATs (Supplementary Figure S6).

Figure 6.

Effect of 1,25(OH)2D and hnRNPC1/C2 knockdown on alternative splicing. (A) Five major types of alternative splicing patterns examined: SE, Skipped Exon; MXE, Mutually Exclusive Exon; A5SS, Alternative 5′ Splice Site; A3SS, Alternative 3′ Splice Site; RI, Retained Intron. (B) Summary of the different differentially expressed genes for each type of alternative splicing event identified by RNA-seq data in MG63 cells with various combinations of no treatment (Ctrl), 1,25(OH)2D (10 nM 1,25, 6 h) exposure and/or hnRNPC1/C2 knockdown (KD). Run from MATS3.0.8 with significant events at FDR < 0.05, |ΔPSI| ≥ 0.05. * includes CYP24A1.

In the absence of hnRNPC1/C2 knockdown, compared to vehicle alone-treated cells, treatment with 1,25(OH)2D resulted in only 10 differentially spliced genes (FDR < 0.05, |ΔPSI| ≥ 0.05) (Figure 6B). Relative to treatment with 1,25(OH)2D only, exposure to combined hnRNPC1/C2 knockdown and 1,25(OH)2D resulted in 4525 differentially spliced genes (FDR < 0.05, |ΔPSI| ≥ 0.05), with one of these being exon 10 skipping in CYP24A1 to generate CYP24A1-variant 2 (ΔPSI: −0.103, Supplementary Figure S4); this CYP24A1 alternative splicing was also observed when comparing hnRNPC1/C2 knockdown alone with 1,25(OH)2D-treatment alone.

Exon skipping (including the skipped exon 10 of CYP24A1) was far more evident when comparing cells treated with 1,25(OH)2D in the presence of hnRNPC1/C2 knockdown versus 1,25(OH)2D alone (Figure 6). In this setting, the dominant determinant of skipped exons continued to be hnRNPC1/C2 knockdown, with 2284 skipped exons being common to 1,25(OH)2D versus 1,25(OH)2D and hnRNPC1/C2 knockdown when compared to vehicle control versus 1,25(OH)2D and hnRNPC1/C2 knockdown (Supplementary Figure S7). However, this comparison also showed that 889 skipped exon events (including exon 10 skipping in CYP24A1) were only observed when comparing 1,25(OH)2D with combined 1,25(OH)2D and hnRNPC1/C2 knockdown (Supplementary Figure S7). These data indicated that a discreet set of AS are associated with 1,25(OH)2D treatment in the presence of hnRNPC1/C2 knockdown, with these events being distinct from AS due to hnRNPC1/C2 knockdown alone.

Of the 3500 differentially expressed genes induced by hnRNPC1/C2 knockdown (see Figure 5), 584 (16.7%) were also shown to be differentially spliced genes (Figure 7A). By contrast none of the 324 differentially expressed genes induced by exposure to 1,25(OH)2D alone were also differentially spliced under these conditions (Figure 7B), and only 1 of 3 differentially spliced genes induced by hnRNPC1/C2 knockdown and 1,25(OH)2D treatment, compared to hnRNPC1/C2 knockdown alone, was also a differentially expressed gene (Figure 7C). When compared to cells only receiving 1,25(OH)2D, cells receiving 1,25(OH)2D in the presence of hnRNPC1/2 knockdown showed 3572 differentially expressed genes and 2510 differentially spliced genes, with 706 (19.8%/28.1%) being both differentially expressed and spliced (Figure 7D). Comparison of vehicle-treated cells with those exposed to 1,25(OH)2D in combination with hnRNPC1/C2 knockdown showed 3844 differentially expressed genes and 2535 differentially spliced genes, with 696 (18.1%/27.5%) being both differentially expressed and spliced (Figure 7E).

Figure 7.

Overlap of differentially expressed genes and differentially spliced genes. Venn diagrams depicting overlap of DEG and DSG induced by: (A) hnRNPC1/C2 knockdown (Ctrl versus KD); (B) 1,25(OH)2D exposure (Ctrl versus 1,25; (10 nM 1,25, 6 h); (C) hnRNPC1/C2 knockdown alone (KD) compared to 1,25(OH)2D exposure combined with hnRNPC1/C2 knockdown (KD versus 1,25+KD); (D) 1,25(OH)2D-treated cells compared to 1,25(OH)2D treatment combined with hnRNPC1/C2 knockdown ((1,25 versus 1,25+KD); and (E) combination of 1,25(OH)2D treatment and hnRNPC1/C2 knockdown relative to control cells (Ctrl versus 1,25+KD). Absolute numbers of DEG and DSG shown within each segment.

DISCUSSION

The over-arching aim of the current study was to investigate a potential link between the established RNA splicing effects of hnRNPC1/C2 and its alternative role in transcriptional regulation as a DNA-binding protein. In particular we postulated that effects of the active form of vitamin D, 1,25(OH)2D, on gene expression may involve coordination of transcriptional regulation and RNA splicing, with hnRNPC1/C2 linking these two effects. Knockdown of hnRNPC1/C2 has been used previously to show how hnRNPC1/C2 acts as a splicing enhancer by enforcing exon inclusion and preventing exon skipping (8). In HeLa cells hnRNPC1/C2 competes with the core splicing factor U2 auxillary factor 65 (U2AF65) to prevent the introduction of cryptic splice sites (27). To date, studies of DNA binding of hnRNPC1/C2 and associated transcriptional effects have focused on over-expression of hnRNPC1/C2 (20). The current study therefore utilized siRNA knockdown of hnRNPC1/C2 to assess both the transcriptional and splicing impact of hnRNPC1/C2 using both candidate gene (CYP24A1) and genome-wide strategies.

In mammalian biology, enhanced expression of CYP24A1 is the most well recognized and robust response in vitro to 1,25(OH)2D exposure (35). In MG63 osteoblastic cells this effect (Figure 1), was associated with binding of hnRNPC1/C2 to both DNA and RNA for CYP24A1 (Figure 2). As reported previously (20), chromatin association of hnRNPC1/C2 is reciprocal with VDR attesting to its competition with the VDR for VDRE binding. By contrast, 1,25(OH)2D-enhanced binding of hnRNPC1/C2 to CYP24A1 pre-mRNA was independent of VDR, suggesting a complementary but mechanistically distinct response to 1,25(OH)2D. The CYP24A1 gene is also characterized by established splice variants that encode enzymatically inactive alternative forms of the 24-hydroxylase protein (21,22). Expression of the heme-binding domain-deficient CYP24A1 variant 2 in vehicle-treated MG63 cells was potently induced by hnRNPC1/C2 knockdown, indicating a significant shift toward exon-skipping in these cells under these conditions. This effect was increased further when 1,25(OH)2D was present to drive transcription of CYP24A1. These data underline the potential role of hnRNPC1/C2 in coordinating regulation of both transcription and splicing of the CYP24A1 gene. The data also highlight different functional responses to knockdown of hnRNPC1/C2 relative to previously reported studies of hnRNPC1/C2 over-expression where the predominant effect was suppression of 1,25(OH)2D-mediated transcription of CYP24A1 (20).

Genome-wide mRNA expression studies for hnRNPs have been limited, especially under conditions of targeted provocation (36). As a consequence, the goal of work presented here was to determine if the concerted actions of hnRNPC1/C2 and 1,25(OH)2D to alter VDR-VDRE-directed expression and alternative splicing of the CYP24A1 gene are reflected in other gene products in the MG63 human osteoblast-like transcriptome. The number of differentially expressed genes associated with knockdown of hnRNPC1/C2 in MG63 cells (3500) was similar to the number of differentially spliced genes (3921), underlining the importance of hnRNPC1/C2 beyond its traditional roles in RNA processing and transport. These data also confirmed observations in Figure 3 showing induction of CYP24A1 by 1,25(OH)2D, and suppression of CYP24A1 with hnRNPC1/C2 knockdown. Even though hnRNPC1/C2 plays a role in lowering the baseline expression of CYP24A1, 1,25(OH)2D was still able to induce expression of CYP24A1 in the presence of hnRNPC1/C2 knockdown (Figures 3 and 5). This indicates that while hnRNPC1/C2 is not essential for 1,25(OH)2D-induced transcription of this gene, it is able to modulate this response.

Further, genome-wide analysis of differentially expressed genes showed that the combinatorial actions of 1,25(OH)2D and altered hnRNPC1/C2 expression are not restricted to CYP24A1. Genome-wide strategies have been used to characterize chromatin-binding of 1,25(OH)2D and associated effects on gene expression (37). The number of differentially expressed genes identified following treatment with 1,25(OH)2D (324) is similar to previous DNA array analysis of primary cultures of human osteoblasts (336) (38). That being said, the impact of hnRNPC1/C2 on 1,25(OH)2D-induced differentially expressed genes is complicated by several factors. First is the strong effect of hnRNPC1/C2 on overall gene expression (Figure 5A and B). Second is the fact that knockdown of hnRNPC1/C2 affected both up- and down-regulation of gene expression by 1,25(OH)2D (Figure 5D). Third is the observation that the impact of hnRNPC1/C2 knockdown on 1,25(OH)2D-induced gene expression was not universal, as many differentially expressed genes specifically induced by 1,25(OH)2D were unaffected by hnRNPC1/C2 knockdown (Figure 5A). Nevertheless, the facts that (i) combined hnRNPC1/C2 knockdown and exposure to 1,25(OH)2D resulted in 912 differentially expressed genes that were distinct from those differentially expressed with hnRNPC1/C2 knockdown alone, and (ii) 87 differentially expressed genes were only induced by 1,25(OH)2D under conditions of hnRNPC1/C2 knockdown (Figure 5B and Supplementary Data S3) endorses a role for hnRNPC1/C2 in the normal regulation of gene expression by 1,25(OH)2D.

Previous studies have characterized transcriptome binding patterns of hnRNPs (39) and transcriptome-wide changes in alternative splicing (27). Data presented here showing genome-wide changes in differentially spliced genes after hnRNPC1/C2 knockdown are consistent with its role in maintaining splicing fidelity (27,40). In the absence or presence of 1,25(OH)2D the majority of differentially spliced genes associated with hnRNPC1/C2 knockdown were skipped exons, endorsing the previously described role of hnRNPC1/C2 in preventing this form of alternative splicing (8). However, despite well documented actions on the transcriptome (41), exposure to 1,25(OH)2D alone had minimal effects on alternative splicing of RNA (Figure 6), not-withstanding the observations in Figure 4 showing enhanced expression of the exon 10-skipped CYP24A1-variant 2 in 1,25(OH)2D-treated cells. This may in part reflect the more stringent statistical analysis based on FDR; analysis of alternative splicing based on P-values showed 139 skipped exons in 1,25(OH)2D-treated cells, including skipping of exon 10 of CYP24A1 (data not shown).

RNA-seq measures steady-state mRNA levels. The altered mRNA levels of certain differentially spliced genes may result from the degradation of specific, alternatively spliced transcripts by the mRNA nonsense-mediated decay pathway. A second (and not mutually exclusive) explanation is that some of the transcriptional actions of hnRNPC1/C2 are linked to its effects on RNA splicing. The latter hypothesis was supported by data from mini-gene constructs, which showed that a CYP24A1 exon 10 skipping event was less likely to occur in the presence of a functional VDRE (Figure 4). Transcriptomic data also showed overlap between the transcription and splicing effects of hnRNPC1/C2 (Figure 7D). However, this only occurred in 20% of differentially spliced genes for 1,25(OH)2D versus combined 1,25(OH)2D and hnRNPC1/C2 knockdown (Figure 7D), indicating that any putative link between gene expression and splicing effects of hnRNPC1/C2 is not universal. Likewise, although comparison of hnRNPC1/C2 knockdown with knockdown plus 1,25(OH)2D resulted in 355 differentially expressed genes, only one of these was also alternatively spliced under the same conditions. These data suggest that although hnRNPC1/C2 plays a dual role in regulating transcription and RNA splicing, these actions are not necessarily linked directly through 1,25(OH)2D-VDR chromatin signaling. However, it is important to recognize that hnRNPC1/C2 also plays a role in the export of transcripts from the nucleus. Specifically hnRNPC1/C2 delineates the length of transcripts for export as either mRNA (>200–300 nucleotides) or snRNA (<200–300 nucleotides) by competing for RNA binding with RNA adapter proteins (42). It is therefore possible that hnRNPC1/C2 links 1,25(OH)2D-VDR-driven transcription with RNA export rather than splicing. In future studies, it will be important to assess the impact of hnRNPC1/C2 on the size and sub-cellular localization of transcripts in bone cells induced by active vitamin D metabolites.