TY - JOUR AU - Meza-Zepeda, Leonardo A. AB - Abstract Differentiation of osteoblasts from mesenchymal stem cells (MSCs) is an integral part of bone development and homeostasis, and may when improperly regulated cause disease such as bone cancer or osteoporosis. Using unbiased high-throughput methods we here characterize the landscape of global changes in gene expression, histone modifications, and DNA methylation upon differentiation of human MSCs to the osteogenic lineage. Furthermore, we provide a first genome-wide characterization of DNA binding sites of the bone master regulatory transcription factor Runt-related transcription factor 2 (RUNX2) in human osteoblasts, revealing target genes associated with regulation of proliferation, migration, apoptosis, and with a significant overlap with p53 regulated genes. These findings expand on emerging evidence of a role for RUNX2 in cancer, including bone metastases, and the p53 regulatory network. We further demonstrate that RUNX2 binds to distant regulatory elements, promoters, and with high frequency to gene 3′ ends. Finally, we identify TEAD2 and GTF2I as novel regulators of osteogenesis. Stem Cells 2014;32:2780–2793 Mesenchymal stem cells, Differentiation, Osteogenesis, Epigenomics, RUNX2, Transcription factors Introduction Mesenchymal stem or stromal cells (MSCs) have multilineage potential, being capable of differentiation towards the osteogenic, adipogenic, chondrogenic as well as other lineages [1]. MSCs are present in bone marrow and other mesenchymal tissues and may migrate to sites of injury or inflammation to contribute to tissue repair [2]. MSCs can be purified and differentiated into osteoblasts in vitro providing a convenient tool for molecular investigation of osteogenesis [1]. More knowledge about this process is important for basic research purposes and clinical application of MSCs in skeletal regenerative medicine [3], as well as to understand its dysregulation in bone disease, such as cancer. In this regard, it is interesting to note that osteosarcoma, the most common primary cancer of bone, may develop due to genetic and epigenetic changes that interrupt osteogenic differentiation from MSCs [4]. Stem cell differentiation is driven by highly regulated changes in gene expression, which at the level of transcription is associated with DNA or chromatin binding of transcriptional regulators and local changes in chromatin landscape [5-7]. The development of high-throughput methods, such as RNA sequencing (RNA-Seq) and chromatin immunoprecipitation and sequencing (ChIP-Seq), has enabled genome-wide characterization of such alterations [8, 9]. For example, RNA-Seq has allowed systematic discovery of alternate promoter and exon usage [10], whereas ChIP-Seq has revealed regions of specific histone modifications, and how these change with altered gene activity [5, 8]. ChIP-Seq has also allowed for detection of genomic binding sites and target genes for a wide range of transcription factors (TFs) [11, 12]. Moreover, integrative genomic analysis has enabled the identification of functional regions such as promoters and enhancers throughout the genome, which subsequently can be investigated for enrichment of transcription factor binding site (TFBS) profiles to aid the discovery of novel regulators of biological processes [5]. Runt-related transcription factor 2 (RUNX2, also known as CBFA1) is indispensable for bone development, demonstrated by skeletal defects in humans with RUNX2 haploinsufficiency [13] or absence of ossified bone in Runx2 nullizygous mice [14, 15]. In osteogenesis, RUNX2 is involved in regulation of proliferation, migration, and commitment and differentiation to the osteogenic lineage [16-19]. RUNX2 regulates gene expression in response to upstream cues, including from the TGFB, BMP, and Wnt signaling pathways [20-23]; however, downstream target genes in osteoblasts are less well characterized. In addition to its role in development, dysregulation of RUNX2 is also known to be involved in tumorigenesis; RUNX2 promotes neoplastic development in synergy with the MYC oncogene in the hematopoietic lineage [24], and recent evidence implicates RUNX2 in tumor growth and metastatic properties of breast and prostate cancer cells [25-27]. In osteosarcoma, RUNX2 is frequently overexpressed, which is linked to poor prognosis [28]. Interestingly, evidence for interplay between RUNX2 and p53 is beginning to emerge, including the suppression of p53-mediated gene expression by RUNX2 in response to DNA damage [29], and evidence for p53-dependent regulation of RUNX2 protein levels [30]. We here describe transcriptomic and epigenomic signatures of osteogenic differentiation from a nontumorigenic immortalized MSC cell line. Furthermore, we report for the first time the genome-wide binding sites of the bone master regulator RUNX2 in human osteoblasts, and reveal more than 2,000 target genes which shed new light on the transcriptional regulation by RUNX2 in bone biology. We also used the epigenomic signatures generated by histone ChIP-Seq to identify promoters with differentiation-induced chromatin alterations and predicted TEAD2 and GTF2I as novel regulators of these promoters. Finally, by silencing TEAD2 and GTF2I expression, we confirmed a role for these TFs in regulation of osteogenic differentiation. Materials and Methods Cell Culture The iMSC#3 cells were culture in minimum essential medium alpha (MEM Alpha) with Glutamax (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) supplemented with 10% fetal calf serum (FCS; Integro, Zaandam, The Netherlands, http://www.integro.nl) (growth medium; GM). For osteogenic differentiation, GM was supplemented with 67 µM ascorbic acid 2-phosphate, 3.5 mM β-glycerol phosphate, and 10 nM dexamethasone (all from Sigma-Aldrich, St. Louis, MO, http://www.sigmaaldrich.com) (differentiation medium; DM) and cells were cultured for 28 days. NHOst cell cultures were expanded and subsequently differentiated for 28 days in media and supplements provided with the Clonetics Normal Human Osteoblast Cell System (Lonza, Walkersville, MD, http://www.lonza.com) as recommended by Lonza. Primary human MSCs were isolated from bone marrow that was collected from the hip of a healthy male and female donor after written informed consent. Subsequently, mononuclear cells were isolated by standard density centrifugation, transferred to a tissue culture flask, and incubated overnight. Nonadherent cells were washed away and the adherent cells were cultured in MEM Alpha with Glutamax (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) supplemented with 20% FCS (Integro, Integro, Zaandam, The Netherlands, http://www.integro.nl). Alizarin Red S Staining Cells were rinsed in phosphate buffered saline (PBS; Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) and fixed in 70% cold ethanol for 1 hour at 4°C and rinsed in water. Cells were then exposed to 40 mM Alizarin Red S (Sigma-Aldrich, St. Louis, MO, http://www.sigmaaldrich.com) for 15 minutes at room temperature (RT), rinsed in water and PBS, and dehydrated in 70% and 100% ethanol. Dried stained cells were imaged. For quantification, Ca2+-bound Alizarin Red S was extracted from cells by addition of 10% acetic acid (vol/vol) and incubated with agitation at RT for 30 minutes. Cells were scraped off, transferred to a tube, vortexed, incubated for 10 minutes at 85°C, and put on ice. Samples were centrifuged for 15 minutes at 1,300 rpm, the supernatant transferred to fresh tube, and 10% NH4OH was added to a final concentration of 1%. The solution was transferred to a 96-well plate and absorbance was read at 405 nm. RNA Isolation, cDNA Synthesis, and qPCR RNA was prepared from cells cultured in GM or DM by TriReagent (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) according to the manufacturer’s recommendations. Isolated RNA was dissolved in water and used as template for complementary DNA (cDNA) synthesis, microarray hybridization, or RNA-Seq. cDNA was made from 1 µg total RNA using the High Capacity RNA-to-cDNA Master Mix (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) according to the manufacturer’s recommendations. Real-time PCR was performed using TaqMan assays detecting RUNX2 and reference transcripts B2M and TBP (Supporting Information Table S1) on an ABI 7900HT Fast Real-Time PCR machine (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com). For quantification of GTF2I, GTF2IRD1, purine-rich element binding protein A (PURA), TEAD2, and ZNF148 transcript levels after siRNA knockdown, cDNA was prepared using the TaqMan Gene Expression Cells-to-Ct kit (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) (early time points) or TriReagent (late time points) and PCR was run with TaqMan assays for each TF and housekeeping genes (TBP and B2M) (Supporting Information Table S1). Relative changes in transcript levels were determined using the method from Pfaffl [31]. SDS-PAGE and Western Blotting For Western blotting, cells were lysed in 2% SDS and 10 mM Tris-HCl, pH 7.6 and 40 µg lysate was separated on a NuPage 4%–12% bis Tris gel (Life Technologies, Carlsbad, CA, http://www.lifetechnologies.com) and transferred to a 0.45 µm PVDF membrane (Merck KGaA, Darmstadt, Germany, http://www.merckgroup.com). Hybridization was performed with anti-RUNX2 or anti-Histone H4 antibody and HRP-conjugated antibodies (Supporting Information Table S2) in Tris-buffered saline with 5% milk, followed by incubation with PierceSuperSignal West Dura Chemiluminescent Substrate (Thermo Scientific, Fremont, CA, http://www.thermoscientific.com) and detection of signal on film. Microarray Hybridization Expression profiling of iMSC#3 and primary human MSCs was performed at the Oslo University Hospital Genomics Core Facility using the Illumina HumanHT12 v4 Expression BeadChip (provided by Illumina, San Diego, CA, http://www.illumina.com) following manufacturer’s protocol. Data extraction and initial quality control of the bead summary raw data were performed essentially as previously described [32]. Expression data were quantile normalized and visualized in GenomeStudio v 2011.1 (Illumina, San Diego, CA, http://www.illumina.com). DNA Methylation Profiling Genomic DNA was purified from cell pellets using the Wizard Genomic DNA purification Kit (Promega, Madison, WI, http://www.promega.com), bisulfite converted using the Zymo EZ DNA Methylation Kit (Zymo Research, Irvine, CA, http://www.zymoresearch.com), and hybridized to Infinium HumanMethylation450 BeadChIP arrays (Illumina, San Diego, CA, http://www.illumina.com) following manufacturer’s protocol. Data extraction and quality control were done in GenomeStudio v2011.1 and the Methylation module v1.9 (both provided by Illumina, San Diego, CA, http://www.illumina.com). For each sample, normalized average beta values for every probe were calculated using beadarray internal controls with no background subtraction. Visualization was performed in GenomeStudio. Chromatin Immunoprecipitation ChIP was performed according to a previously described protocol [33] with some modifications. Briefly, 2–6 million cells per ChIP were cross-linked in 1% formaldehyde at RT followed by addition of glycine to a 0.125 M final concentration. Cells were washed in cold PBS, pelleted, and a 6× pellet volume of lysis buffer (50 mM Tris-HCl, pH 8.0, 10 mM EDTA, 1% [wt/vol] SDS) was added. Chromatin was fractionated using a Covaris S220 instrument and the lysate spun at 15,000g in a table top centrifuge to remove insoluble components. For immunoprecipitation, the lysate was diluted 1:10 in RIPA buffer (10 mM Tris-HCl, pH 7.5, 140 mM NaCl, 1 mM EDTA, 0.5 mM EGTA, 1% [vol/vol] Triton X-100, 0.1% [wt/vol] SDS, and 0.1% [wt/vol] Na-deoxycholate) and added to protein A beads precomplexed with 10 µg antibody (Supporting Information Table S2) and incubated for 2 hours at 4°C. Beads were washed 3× in RIPA buffer and 1× in TE pH8 followed by addition of elution buffer (20 mM Tris-HCl, pH 7.5, 5 mM EDTA, 50 mM NaCl, 1% [wt/vol] SDS) containing DNase-free RNase (Roche, Indianapolis, IN, http://lifescience.roche.com) and incubation at 37°C for 30 minutes. Proteinase K was added and the samples incubated at 68°C for 5 hours for reversal of crosslinks. DNA was extracted once with phenol/chloroform/isoamylalcohol (25:24:1) and once with chloroform/isoamylalcohol (24:1) followed by ethanol precipitation. Precipitated DNA was resuspended in water and used for high-throughput sequencing. For histone and RUNX2 ChIPs from iMSC#3 cells, enrichment was determined relative to a control ChIP with anti-H3 antibody. For RUNX2 ChIP from NHOst cells, enrichment was determined relative to input DNA. Library Preparation and Sequencing RNA sequencing libraries were made using the Illumina TruSeq RNA Sample Preparation Kit v2 (RS-122–2001, Illumina, San Diego, CA, http://www.illumina.com). Total RNA was quality assessed on an Agilent Technologies (Santa Clara, CA, http://www.genomics.agilent.com) 2100 Bioanalyzer and processed according to manufacturer’s instructions. For each sample 1 µg total RNA was used as starting material. For ChIP-Seq the DNA Sample Prep Kit (IP-102–1001; Illumina, San Diego, CA, http://www.illumina.com) was used to build DNA libraries of 40 ng precipitated DNA per sample. The samples were processed as recommended except for size selection of adapter ligated DNA which was done after instead of before amplification. The resulting libraries were sequenced on the Genome Analyzer IIx using TruSeq SBS Kit v5—GA (FC-104–5001; Illumina, San Diego, CA, http://www.illumina.com), performing a paired-end run with 75 bp sequence length (2 × 75) for the RNA libraries and a single read run with 36 bp sequence length for the ChIP libraries. Real-time analysis and base calling were done by Illumina’s software packages SCS2.9/RTA1.9 and Off-line Basecaller v1.9. Sequence Analysis Aligned reads were extended to 200 bp to reflect the ChIP fragments. For histone ChIP-Seq, 500 bp bins were generated genome-wide for all experiments, and for each bin the coverage was calculated for both ChIP and control which were used to calculate fold change values and False Discovery Rate corrected (FDR-corrected) Poisson p-values. Thresholds for calling enriched bins were a minimum of 100 reads, fivefold enrichment against control and p-value ≤ .05. Threshold for calling differentially modified bins was an absolute log2-fold change ≥2. Peak detection was performed for all RUNX2 experiments, followed by fold change ≥4 and significance estimation (p ≤ .05). Enforcing equal read counts was used as a normalization method for all experiments. All ChIP-Seq analysis was performed using custom R/Bioconductor scripts. Gene Association For promoter marks, Ensembl transcripts were used to identify transcripts with a 5′ end within 2 kb of bins. The transcript with the closest 5′ end was selected. Microarray expression data were integrated using Refseq IDs. Motif Enrichment Analysis For each histone modification, the top 5,000 increasing and decreasing bins were selected. Bins were merged if they were next to each other. The regions were scanned with motifs from Jaspar, Transfac, and Uniprobe with a score threshold of 80%. Fisher test was used to compare the number of sequences with at least one hit and the number of sequences with no hits between the two sets. Wilcoxon signed-rank test was used to compare the rank of total number of hits between the two sets. Motifs with a p-value < .05 from both tests were defined as significant. Ingenuity Pathway Analysis All analysis were performed with species filter set to “human and mouse” and confidence set to “experimentally confirmed.” RNA-Seq Analysis Tophat [34] was used as an exon junction aware aligner, and Cufflinks [35] was used for quantification of expression levels. Tracks were generated for visualization in the UCSC Genome Browser. Small Interfering RNA Knockdown During Osteogenic Differentiation Cells were plated in 250 µl GM in 48-well plates at 3,000 cells per cubic centimeter 24 hours prior to transfection. Silencer select siRNAs purchased from Life Technologies (Supporting Information Table S3) were introduced to cells using Lipofectamine 2000 according to the manufacturer’s recommendations. Briefly, a transfection cocktail of 0.25 µl 50 nm siRNA and 0.25 µl Lipofectamine 2000 in 50 µl OptiMEM was prepared and added to each well. Cells were exposed to siRNA for 20 hours before removal and addition of DM. The procedure was repeated on Day 7, Day 14, and Day 21 where prior to transfection DM was replaced with GM for the duration of siRNA exposure. Cells were harvested 72 hours after knockdown and the level of mRNA from the targeted transcripts was determined by quantitative PCR. Accession of Data All sequencing data have been made public by submission to the European Bioinformatics Institute (ENA) with ENA Study Accession numbers ERP003789 (RNA-Seq data) and ERP003787 (ChIP-Seq data). The DNA methylation data have been deposited in the Gene Expression Omnibus repository under the accession number GSE49585. Results Differentiation of iMSC#3 Cells to the Osteogenic Lineage Elicits Global Gene Expression Changes To investigate osteogenic differentiation from MSCs we used a nontumorigenic cell line (iMSC#3) generated in our laboratory which is immortalized by ectopic expression of telomerase (Skårn et al., manuscript submitted for publication). This cell line can be manipulated by transfection and retains multilineage differentiation potential after many population doublings. Moreover, the iMSC#3 cells strongly resemble primary human MSCs at the gene expression level, indicating that immortalization did not drastically change their MSC identity (Supporting Information Fig. S1). For osteogenesis, iMSC#3 cells were cultured in DM for 28 days and the formation of a mineralized matrix was indicated by robust Alizarin Red S staining (Fig. 1A). To further examine the extent of differentiation, global gene expression from undifferentiated (Day 0) and differentiated (Day 28) cells was measured by RNA-Seq analysis and compared. Transcripts corresponding to a total of 1,462 and 1,695 genes were found to be upregulated or downregulated, respectively, at least twofold or more upon differentiation (Fig. 1B). Replicate experiments showed high-level of correlation (Supporting Information Fig. S2). Using ingenuity pathway analysis (IPA) we find that upregulated genes were associated with proliferation, cell death, skeletal development, and cellular movement, whereas downregulated genes were connected to cell cycle regulation (Fig. 1C). Known osteogenic markers, including ALPL (alkaline phosphatase), RUNX2, ATF4 (activating transcription factor 4), OMD (osteomodulin), and FOXO1 (forkhead box O1) were found to be upregulated, whereas MSC markers including MCAM (melanoma cell adhesion molecule), ENG (endoglin), CLCF1 (cardiotrophin-like cytokine factor 1), TPM1 (tropomyosin 1 alpha), and CD44 (CD44 molecule [Indian blood group]) were downregulated (Fig. 1D). Open in new tabDownload slide Osteogenic differentiation of immortalized mesenchymal stem cell (iMSC)#3 cells induces gene expression changes and mineralization. (A): Alizarin red S staining of iMSC#3 cells before (Day 0) and after 28 days of culture in differentiation medium. (B): Scatter plot of mean transcript levels from two biological replicate differentiation experiments, as determined by RNA-Seq. The mean log10 fragments per kilobase of exon per million fragments mapped per gene values are plotted for iMSC#3 Day 0 and Day 28 cells. Red spots indicate a twofold or more change in transcript level between Day 0 and Day 28 cells, with the number of upregulated and downregulated genes shown in red. (C): Ingenuity pathway analysis (IPA) of molecular and cellular functions for genes with more than or equal to twofold upregulation or downregulation in transcript level, with the significance of association to each IPA category indicated (p-value interval). (D): Ratio of transcript levels for osteogenic markers ALPL, RUNX2, ATF4, OMD, and FOXO1, and for MSC markers MCAM, ENG, CLCF1, TPM1, and CD44 from two biological replicates of RNA-Seq. Open in new tabDownload slide Osteogenic differentiation of immortalized mesenchymal stem cell (iMSC)#3 cells induces gene expression changes and mineralization. (A): Alizarin red S staining of iMSC#3 cells before (Day 0) and after 28 days of culture in differentiation medium. (B): Scatter plot of mean transcript levels from two biological replicate differentiation experiments, as determined by RNA-Seq. The mean log10 fragments per kilobase of exon per million fragments mapped per gene values are plotted for iMSC#3 Day 0 and Day 28 cells. Red spots indicate a twofold or more change in transcript level between Day 0 and Day 28 cells, with the number of upregulated and downregulated genes shown in red. (C): Ingenuity pathway analysis (IPA) of molecular and cellular functions for genes with more than or equal to twofold upregulation or downregulation in transcript level, with the significance of association to each IPA category indicated (p-value interval). (D): Ratio of transcript levels for osteogenic markers ALPL, RUNX2, ATF4, OMD, and FOXO1, and for MSC markers MCAM, ENG, CLCF1, TPM1, and CD44 from two biological replicates of RNA-Seq. Epigenomic Profiling of Differentiating iMSC#3 Cells To investigate genome wide epigenetic changes in osteogenesis, we profiled the distribution of histone H3 modifications representing different functional states of chromatin by ChIP-Seq at Day 0 and Day 28. We examined the distribution of H3K4m3, enriched at promoters [8, 36], H3K9ac and H3K27ac, enriched at active regulatory elements [37, 38], H3K36m3, found in bodies of actively transcribed genes [36], and the repressive modification H3K27m3 which typically show a broad distribution pattern covering promoter, coding, and intergenic regions [8, 39]. ChIP-Seq and RNA-Seq enrichment profiles for the ALPL locus at Day 0 and Day 28 are illustrated (Fig. 2A). The profiles show upregulated transcript levels at Day 28 coinciding with increased enrichment of H3K4m3 and H3K36m3 in promoter and coding regions, respectively. Multiple H3K27ac and H3K9ac peaks are indicative of active regulatory elements, whereas repressive H3K27m3 methylation is found upstream of the gene (Fig. 2A). Open in new tabDownload slide Osteogenic differentiation elicits global chromatin changes. (A): Snapshot from UCSC Genome Browser of the ALPL gene showing RNA-Seq enrichment profiles and ChIP-Seq enrichment profiles of H3K4m3, H3K9ac, H3K27ac, and H3K27m3 in iMSC#3 cells at Day 0 and Day 28. Scale indicates intensity of enrichment. RefSeq annotations of three COL1A1 transcript variants are depicted. (B): Clustered H3 modification dynamics genome wide. For plotting, the log2-fold change of the ratio (Day 28/Day 0) of the number of reads pr 500 bp region was clustered by principal component. (C): Violin plots of the log2 ratios (Day 28/Day 0) of the number of reads pr 500 bp-region genome wide or in regions covering ±2 kb of ENSEMBL annotated 5′ gene ends (±2 kb of 5′) for H3K9ac, H3K27ac, H3K4m3, and H3K27m3. Only 500 bp-regions with a minimum of 100 reads in at least one of the conditions (Day 0 or Day 28) were included. Total number of regions for each modification is shown above each plot. (D): Height of bars indicates number of promoter regions with increase or decrease (log2-fold ≥ 2 or −2, respectively) of H3K4m3, acetylation (H3K9ac and/or H3K27ac), and H3K27m3 at Day 28 relative to Day 0. Red or blue color indicates the number of promoters associated with genes that show a corresponding increase or decrease in transcript level, respectively. (E): Ingenuity pathway analysis (IPA) of bio functions, including molecular and cellular functions, diseases and disorders, and physiological system development and function, for genes with increased or decreased (log2-fold ≥ 2 or −2, respectively) enrichment of acetylation (H3K9ac and/or H3K27ac), H3K4m3, or H3K27m3 in their promoter regions at Day 28 relative to Day 0, irrespective of change in gene expression. The significance of association to each IPA category is indicated (p-value interval). Abbreviation: ALPL, alkaline phosphatase. Open in new tabDownload slide Osteogenic differentiation elicits global chromatin changes. (A): Snapshot from UCSC Genome Browser of the ALPL gene showing RNA-Seq enrichment profiles and ChIP-Seq enrichment profiles of H3K4m3, H3K9ac, H3K27ac, and H3K27m3 in iMSC#3 cells at Day 0 and Day 28. Scale indicates intensity of enrichment. RefSeq annotations of three COL1A1 transcript variants are depicted. (B): Clustered H3 modification dynamics genome wide. For plotting, the log2-fold change of the ratio (Day 28/Day 0) of the number of reads pr 500 bp region was clustered by principal component. (C): Violin plots of the log2 ratios (Day 28/Day 0) of the number of reads pr 500 bp-region genome wide or in regions covering ±2 kb of ENSEMBL annotated 5′ gene ends (±2 kb of 5′) for H3K9ac, H3K27ac, H3K4m3, and H3K27m3. Only 500 bp-regions with a minimum of 100 reads in at least one of the conditions (Day 0 or Day 28) were included. Total number of regions for each modification is shown above each plot. (D): Height of bars indicates number of promoter regions with increase or decrease (log2-fold ≥ 2 or −2, respectively) of H3K4m3, acetylation (H3K9ac and/or H3K27ac), and H3K27m3 at Day 28 relative to Day 0. Red or blue color indicates the number of promoters associated with genes that show a corresponding increase or decrease in transcript level, respectively. (E): Ingenuity pathway analysis (IPA) of bio functions, including molecular and cellular functions, diseases and disorders, and physiological system development and function, for genes with increased or decreased (log2-fold ≥ 2 or −2, respectively) enrichment of acetylation (H3K9ac and/or H3K27ac), H3K4m3, or H3K27m3 in their promoter regions at Day 28 relative to Day 0, irrespective of change in gene expression. The significance of association to each IPA category is indicated (p-value interval). Abbreviation: ALPL, alkaline phosphatase. Differentiation Is Associated With Changes of the Chromatin Landscape For global analysis of histone modification changes, the fold change (Day 28/Day 0) of reads aligning to 500 bp regions throughout the genome was calculated. In general, we find that differentiation-induced changes of H3K4m3, H3K9ac, H3K27ac, and H3K36me3 enrichment show positive correlation, whereas a change of repressive methylation largely anticorrelate with changes of activating marks (Fig. 2B). Further we find that induction of differentiation was associated with an overall genome-wide loss of acetylation at H3K9 and H3K27, whereas total levels of H3K4m3 and H3K27m3 levels remained relatively constant (Fig. 2C; all regions). Restricting the analysis to regions surrounding gene promoters we found an overall loss of acetylation, however less pronounced than genome wide, meanwhile H3K27m3 enrichment was overall increased at promoters (Fig. 2C; ±2 kb of 5′). Enrichment ratios for H3K4m3 remained unchanged at promoters (Fig. 2C; ±2 kb of 5′). A large number of promoter regions showed increased or decreased histone modification enrichment without a change in transcript levels from the corresponding gene (Fig. 2D). In particular, reduction of acetylation occurred at several thousand promoters, with a great extent of overlap between H3K9 and H3K27 (data not shown), affecting genes involved in gene regulation, cellular survival, growth, and proliferation (Fig. 2D, 2E). Interestingly, genes with increased acetylation or decreased H3K27m3 levels at their promoters were strongly associated with oncogenesis. Increased H3K27m3 enrichment was strongly associated with nervous system development (Fig. 2E). For H3K4m3, increased enrichment was associated with developmental processes whereas decreased enrichment was associated with cell cycle (Fig. 2E). In general, a positive correlation was found between upregulated expression and increased H3K4m3, H3K9ac, and H3K27ac enrichment (Supporting Information Fig. S3). For H3K27m3 changes, little correlation was found with regard to expression on a global scale (Supporting Information Fig. S3) and only 133 and 81 upregulated or downregulated genes showed significant changes of enrichment in the promoter region, respectively (data not shown). IPA analysis found no significant enrichments of functions for these genes. In addition to histone modifications, we investigated global changes in DNA methylation status at 450,000 CpG sites in Day 0 and Day 28 cells using the HumanMethylation450 BeadChIP assay. The assay covers 99% of RefSeq genes, with an average coverage of 17 CpG sites per gene. We found altered methylation status at only 1,273 CpG sites, with decreased methylation at 1,121 sites and increased methylation 152 sites. Altered methylation status was associated with 710 genes, which is an average of 1.6 CpG sites per gene and with only 27 genes affected at more than two CpG sites (Supporting Information Fig. S4; Supporting Information Table S4). Thus, DNA methylation remains relatively constant and does not appear to be an important regulatory factor governing osteogenic differentiation from MSCs. Alternative Transcript Variants of the RUNX2 Gene Given the role of RUNX2 as a master regulator of osteogenic differentiation, we investigated expression of the gene in more detail. In Day 0 cells, mainly transcripts originating from the P1 promoter were detected by RNA-Seq (Fig. 3A). Upon induction of differentiation, the upstream P2 promoter was induced as shown by an increase of RNA-Seq reads from the two first exons at Day 28 (Fig. 3A). In addition, increased enrichment of H3K4m3 at the P2 promoter and H3K36m3 in the coding region between P2 and P1 promoter at Day 28 supports a differentiation-induced activation of the P2 transcriptional start site (Fig. 3A). TaqMan qPCR using primers covering two different exon junctions, as illustrated in Figure 3A, 3B, validated a strong increase in transcripts from the P2 promoter and with only a slight increase from the P1 promoter (Fig. 3C); however, no increase in protein level was detected (Fig. 3D). Interestingly, expression of a non-RefSeq but ENCODE-annotated exon located between exons 2 and 3 of the RUNX2 gene was induced by differentiation (Fig. 3A, 3B; asterisk). Transcripts containing this exon do not read through to any of the downstream exons as evidenced by RNA-Seq alignment reads (Supporting Information Fig. S5) and constitute a significant part of the upregulated transcripts from the RUNX2 locus. This short RNA corresponds to an ENCODE-annotated transcript (Fig. 3A) with an unknown function, but predicted to be protein coding. Open in new tabDownload slide Expression of RUNX2 in osteogenic differentiation. (A): RNA-Seq profiles and ChIP-Seq profiles of H3K4m3 and H3K36m3 of the RUNX2 locus in iMSC#3 cells at Day 0 and Day 28 of osteogenic differentiation. Transcript variants from the P2 promoter (RefSeq 1 and 2) and from the P1 promoter (RefSeq 3) are depicted. I and II indicate the amplified regions of two different TaqMan assays. Arrows indicate gene 5′ ends (RefSeq) and direction of transcription. * indicates transcripts from an ENCODE, but not RefSeq, annotated exon. The area within the dashed box is shown in more detail in (B). The 5′ region of the SUPT3H (suppressor of Ty 3 homolog [Saccharomyces cerevisiae]) gene encoded on the minus strand is depicted. (B): More detailed view of the area within the dashed box in (A). RNA-Seq profiles show upregulation of transcripts from the two first exons of the RUNX2 gene together with the third exon of the non-RefSeq ENCODE-annotated transcript at Day 28 but not at Day 0. (C): Fold change (Day 28/Day 0) in RUNX2 transcript levels as determined by TaqMan quantitative PCR using assay I and II as shown in (A). (D): Western blot showing bands specific for RUNX2 in Day 0 and Day 28 cells, and histone H4 as a loading control, with protein size in kilodalton as determined by SDS-PAGE (kDa). Open in new tabDownload slide Expression of RUNX2 in osteogenic differentiation. (A): RNA-Seq profiles and ChIP-Seq profiles of H3K4m3 and H3K36m3 of the RUNX2 locus in iMSC#3 cells at Day 0 and Day 28 of osteogenic differentiation. Transcript variants from the P2 promoter (RefSeq 1 and 2) and from the P1 promoter (RefSeq 3) are depicted. I and II indicate the amplified regions of two different TaqMan assays. Arrows indicate gene 5′ ends (RefSeq) and direction of transcription. * indicates transcripts from an ENCODE, but not RefSeq, annotated exon. The area within the dashed box is shown in more detail in (B). The 5′ region of the SUPT3H (suppressor of Ty 3 homolog [Saccharomyces cerevisiae]) gene encoded on the minus strand is depicted. (B): More detailed view of the area within the dashed box in (A). RNA-Seq profiles show upregulation of transcripts from the two first exons of the RUNX2 gene together with the third exon of the non-RefSeq ENCODE-annotated transcript at Day 28 but not at Day 0. (C): Fold change (Day 28/Day 0) in RUNX2 transcript levels as determined by TaqMan quantitative PCR using assay I and II as shown in (A). (D): Western blot showing bands specific for RUNX2 in Day 0 and Day 28 cells, and histone H4 as a loading control, with protein size in kilodalton as determined by SDS-PAGE (kDa). Genome-Wide Characterization of Binding Sites for the Osteogenic Master Regulator RUNX2 In order to identify RUNX2 target genes in osteoblasts we performed ChIP-Seq using a monoclonal antibody against RUNX2 to detect binding sites genome-wide in differentiated iMSC#3 (iMSC#3 Day 28) cells. For validation, RUNX2 ChIP-Seq was also performed on primary NHOst osteoblasts. In iMSC#3 cells, a total of 9,549 RUNX2 bound regions were identified by annotation to the closest 5′ gene end within a 50 kb region (Supporting Information Table S5). The binding sites were distributed close to genes and at distal sites (Fig. 4A). Ninety-two percent of the 500 most enriched peaks in iMSC#3 Day 28 cells and 76% of the 4,000 peaks most enriched in Day 28 cells were also enriched in NHOst cells (Fig. 4B). De novo motif discovery at RUNX2 peaks produced a motif identical with the RUNX family motif in 99% of analyzed sequences, further supporting specific enrichment for RUNX2 (Fig. 4C). Open in new tabDownload slide Genome-wide analysis of RUNX2 binding sites by ChIP-Seq in iMSC#3 Day 28 cells. (A): Pie chart showing the number of RUNX2 enriched regions at intergenic or intragenic sites (more than 2 kb from ENSEMBL-annotated gene 5′ or 3′ ends), or in regions covering ±2 kb of gene 5′ or 3′ ends. (B): Overlap of regions with RUNX2 binding between iMSC#3 Day 28 and NHOst cells for the 500 or 4,000 most strongly enriched sites in iMSC#3 Day 28 cells. (C): Motif returned by de novo motif discovery from DNA sequences in regions bound by RUNX2 in iMSC#3 Day 28 cells and the previously described RUNX1 motif from the Jaspar database. (D): Ingenuity pathway analysis (IPA) analysis of molecular and cellular functions or diseases and disorders for 2,167 RUNX2 target genes defined by binding sites within ±2 kb of 5′ or 3′ gene end, with the significance of association (p-value interval) and the number of target genes (# Genes) associated with each IPA category indicated. (E): IPA analysis of upstream regulators for RUNX2 target genes as defined in (D), with the significance of association (p-value) and the number of target genes (# Genes) associated with each upstream regulator indicated. (F): Plot showing the density of RUNX2 binding distributed around the 5′ and 3′ ends of genes. (G): Graph shows the percentage of RUNX2 enriched regions overlapping with H3K4m3, H3K27ac, or H3K4m3 and H3K27ac (Both) enrichment, or with neither modification (None). Overlap is shown for all RUNX2 enriched regions (All), regions within ±2 kb of 5′ gene ends (5′) or within ±2 kb of 3′ gene ends (3′). (H): ChIP-Seq enrichment profiles of RUNX2, H3K4m3, H3K27ac, and H3 at the COL1A1 locus in iMSC#3 Day 28 cells. Scale indicates intensity of enrichment. (I): IPA analysis of molecular and cellular functions genes with 5′ only, 3′ only, or 5′ and 3’ binding of RUNX2 in iMSC#3 Day 28 cells (as defined in (D)), with the significance of association to each category indicated (p-value interval). The number of genes associated with 5′ only, 3′ only, or 5′ and 3′ enrichment of RUNX2 is shown in parenthesis. (J): IPA analysis of predicted upstream regulators of gene categories described in (I), with p-values indicating the significance of association. Abbreviation: iMSC, immortalized mesenchymal stem cell. Open in new tabDownload slide Genome-wide analysis of RUNX2 binding sites by ChIP-Seq in iMSC#3 Day 28 cells. (A): Pie chart showing the number of RUNX2 enriched regions at intergenic or intragenic sites (more than 2 kb from ENSEMBL-annotated gene 5′ or 3′ ends), or in regions covering ±2 kb of gene 5′ or 3′ ends. (B): Overlap of regions with RUNX2 binding between iMSC#3 Day 28 and NHOst cells for the 500 or 4,000 most strongly enriched sites in iMSC#3 Day 28 cells. (C): Motif returned by de novo motif discovery from DNA sequences in regions bound by RUNX2 in iMSC#3 Day 28 cells and the previously described RUNX1 motif from the Jaspar database. (D): Ingenuity pathway analysis (IPA) analysis of molecular and cellular functions or diseases and disorders for 2,167 RUNX2 target genes defined by binding sites within ±2 kb of 5′ or 3′ gene end, with the significance of association (p-value interval) and the number of target genes (# Genes) associated with each IPA category indicated. (E): IPA analysis of upstream regulators for RUNX2 target genes as defined in (D), with the significance of association (p-value) and the number of target genes (# Genes) associated with each upstream regulator indicated. (F): Plot showing the density of RUNX2 binding distributed around the 5′ and 3′ ends of genes. (G): Graph shows the percentage of RUNX2 enriched regions overlapping with H3K4m3, H3K27ac, or H3K4m3 and H3K27ac (Both) enrichment, or with neither modification (None). Overlap is shown for all RUNX2 enriched regions (All), regions within ±2 kb of 5′ gene ends (5′) or within ±2 kb of 3′ gene ends (3′). (H): ChIP-Seq enrichment profiles of RUNX2, H3K4m3, H3K27ac, and H3 at the COL1A1 locus in iMSC#3 Day 28 cells. Scale indicates intensity of enrichment. (I): IPA analysis of molecular and cellular functions genes with 5′ only, 3′ only, or 5′ and 3’ binding of RUNX2 in iMSC#3 Day 28 cells (as defined in (D)), with the significance of association to each category indicated (p-value interval). The number of genes associated with 5′ only, 3′ only, or 5′ and 3′ enrichment of RUNX2 is shown in parenthesis. (J): IPA analysis of predicted upstream regulators of gene categories described in (I), with p-values indicating the significance of association. Abbreviation: iMSC, immortalized mesenchymal stem cell. We identified 2,167 ENSEMBL genes bound by RUNX2 within ±2 kb of either the 5′ or 3′ gene end and defined those genes as RUNX2 target genes (Supporting Information Table S6). IPA analysis of RUNX2 target genes showed highly significant association to molecular and cellular functions involved in regulation of cellular growth and proliferation, cellular movement, cell death and survival, development, and cellular assembly and organization (Fig. 4D, Supporting Information Table S7A). A strong enrichment for cancer-related functions was found (Fig. 4D, Supporting Information Table S7B); including identification of 154 RUNX2 target genes known to be part of the p53 regulated network according to IPA database (Fig. 4E, Supporting Information Table S8). Moreover, RUNX2 targets genes were associated with platelet-derived growth factor beta-beta (PDGF BB), Fas ligand (FAS), Alpha-catenin, and transforming growth factor beta 1 (TGFB1) signaling pathways (Fig, 4E, Supporting Information Table S8). Interestingly, we found targeting of genes encoding known protein interaction partners of RUNX2, including ATF4, STAT1, and FOS (Supporting Information Table S6, Supporting Information Fig. S9), suggesting regulatory loops between these proteins. To investigate potential functional implications of RUNX2 binding, we examined the expression dynamics of RUNX2 target genes as well as histone modification dynamics of RUNX2 bound sites. First, we found that 372 upregulated and 339 downregulated genes are bound by RUNX2 at Day 28 (Supporting Information Fig. S6A, Supporting Information Table S9), which means that around 30% of all differentially expressed genes in our osteogenic model are bound by RUNX2 at Day 28. Most of these genes have not previously been described as RUNX2 targets. IPA analysis shows that the genes are associated with processes including proliferation, migration, and cell death and development, with a significant fraction of the upregulated genes being linked to skeletal development and function (Supporting Information Fig. S6B). We then went on to investigate whether RUNX2 bound genomic regions associated with differentially expressed genes also showed changes in histone modification enrichment status at Day 28 compared to Day 0; however, less than 1% of the regions were affected (Supporting Information Fig. S6C). This indicates that histone modification change is not a major mediator of RUNX2 binding at gene proximal sites during differentiation, at least not for the modifications investigated here. Finally, we looked at the histone modification status of all RUNX2 bound regions and found that between 100 and 200 regions are affected by H3K4m3, H3K9ac, or H3K27ac change, with more regions showing increased enrichment than decreased (Supporting Information Fig. S6D). Almost no RUNX2 bound sites were affected by H3K27m3 enrichment change. We also find that histone modification changes at RUNX2 binding sites are more frequent at distal sites (>2 kb from the start or end of gene) compared to gene proximal sites (±2 kb within a gene 5′ or 3′ end). Overall, RUNX2 does not seem to bind genomic sites in response to histone enrichment changes induced by differentiation. An outline of the workflow used to identify RUNX2 bound regions is shown in Supporting Information Figure S7. From the total list of RUNX2 peaks, we identified 97 microRNA (miRNA) miRBase-annotated sequences to be associated with RUNX2 binding (Supporting Information Table S5). Sequences from miRBase represent the predicted hairpin portion of a miRNA transcript, and does not encompass information on the full length of the primary transcript. Thus, for these genes we cannot define RUNX2 binding to 5′ or 3′ ends, as for other types of genes. It is however interesting to note that strong RUNX2 enrichment was found in proximity to annotated miRNAs with a known role in bone biology, including the mir23a-27a-24-2 cluster, mir-22, mir-21, mir-143, mir-145, mir-221, and mir-214 [40-47] (Supporting Information Fig. S8). Other miRNA genes associated with strong RUNX2 enrichment include mir-4530, mir-3619, mir-let7a3 and let7b, mir-612, and mir-199a1 (Supporting Information Fig. S8). RUNX2 bound with similar frequency at the 5′ and 3′ ends of genes, and at higher density in intragenic regions than intergenic regions (Fig. 4F). Approximately half of all RUNX2 bound regions overlapped with H3K27ac enrichment, and about one-third of these regions were also bound by H3K4m3 (Fig. 4G). This indicates that RUNX2 binds both active promoters and active enhancers. A small fraction overlapped with H3K4m3 only. A larger overlap with H3K27ac (approximately 75%) and H3K4m3 (50%) enrichment was found at the 5′ end of genes, indicating that RUNX2 promoter binding primarily occurs at transcriptionally active genes. At gene 3′ regions there was less overlap with H3K27ac and H3K4m3 than at promoters (Fig. 4G). RUNX2 and histone ChIP-Seq profiles at the COL1A1 gene are shown as an example of these patterns (Fig. 4H). Among all RUNX2 target genes, 901 had 5′ binding only, 724 had 3′ binding only, whereas 541 genes had both 5′ and 3′ binding (Fig. 4I; Supporting Information Table S10A–S10C). IPA analysis showed that genes with different patterns of 5′ and 3′ RUNX2 binding were associated with different cellular functions (Fig. 4I). Genes with both 5′ and 3′ binding were significantly associated with cellular movement and development, as well as ERK, PDGF BB, and TGFB upstream signaling, whereas 5′ binding was associated with cell-to-cell signaling and 3′ binding with cell death and survival (Fig. 4I, 4J). More examples of RUNX2 ChIP-Seq enrichment profiles at target genes are shown in Supporting Information Fig. S9). Identification of Novel TFs Regulating Osteogenesis To identify other and potentially novel TFs involved in regulation of osteogenic differentiation, we performed a TFBS enrichment analysis in promoters and enhancers with increased H3K4m3, H3K9ac, and H3K27ac at Day 28 relative to Day 0 using the Uniprobe, Jaspar, and Transfac databases. The H3K4m3 regions were the only set of regions with clearly significant differential motif enrichment, thus we focused our analysis on this histone mark. To further narrow down the search, we looked for motifs exclusively in promoters (±2 kb from 5′ gene end) of genes with increased H3K4m3 enrichment in their promoter after differentiation, and which also had upregulated gene expression, as determined by RNA-Seq, at Day 28 relative to Day 0. An outline of the workflow used to identify novel TFs is shown in Supporting Information Figure S10. Motifs for TFs with nondetectable expression in iMSC#3 cells were excluded and the remaining most significant motifs were for TEA domain family member 2 (TEAD2), PURA, GTF2I repeat domain containing 1 (GTF2IRD1), zinc finger protein 148 (ZNF148), and general transcription factor IIi (GTF2I) (Fig. 5A). Of these, the TEAD2 motif was most frequent and there was extensive overlap in promoters harbouring each of the motifs, especially for TEAD2, GTF2IRD1, and PURA (Fig. 5B). IPA analysis of predicted TEAD2 target genes shows association to cell morphology, development (including skeletal development; data not shown), growth and proliferation, cellular movement, and cell death (Fig. 5C). To investigate a potential regulatory role in osteogenesis, expression of each TF was repeatedly knocked down by two separate small interfering RNAs (siRNAs) during differentiation of iMSC#3 cells. Treatment with each of the siRNAs was performed in 7–10 replicate wells and knockdown efficiency was measured by quantitative PCR 72 hours after each knockdown (Supporting Information Fig. S11). At Day 28, the cell cultures were stained by Alizarin Red S, imaged, and staining was taken as a measure of differentiation. For GTF2I and TEAD2, staining was significantly reduced compared to control cells for both siRNAs used (Fig. 5D, 5E), strongly indicating a specific inhibition of differentiation by the knockdown. For GTF2IRD1 and ZNF148, a significant reduction in staining was observed for one of the siRNAs but not for the other (Supporting Information Fig. S12). For PURA knockdown, we did not observe any significant effect on mineralization (Supporting Information Fig. S12). Open in new tabDownload slide Identification of novel TFs regulating osteogenesis. (A): DNA motifs representing the transcription factor binding site (TFBS) for TEAD2, PURA, GTF2IRD, ZNF148, and GTF2I TFs from TRANSFAC with accession number (Acc), rate of enrichment in upregulated versus downregulated promoters (Ratio), and enrichment p-values (Wilcoxon test) (p-value). (B): Figure indicates presence (black field) or absence (white field) of ZNF148, GTF2I, PURA, GTF2IRD1, or TEAD2 TFBS at promoter regions (±2 kb of TSS) of genes with increased (log2-fold ≥1) gene expression and increased (log2-fold ≥1) H3K4m3 enrichment at Day 28 relative to Day 0. Each black or white field represents a promoter region. (C): Ingenuity pathway analysis (IPA) analysis of molecular and cellular functions for genes enriched in TEAD2 motif is shown, with the significance of association (p-value) and the number of genes (# Genes) associated with each IPA category indicated. (D): Alizarin red S staining of iMSC#3 cells at Day 28 of differentiation after siRNA-mediated knockdown of GTF2I or TEAD2 expression at Day 0, Day 7, Day 14, and Day 21 with two separate siRNAs per gene (1 or 2). Control cells were cultured in DM and exposed to a control siRNA (Ctrl) or no siRNA (DM), or cultured in GM with no siRNA (GM). A representative image for each condition is shown. (E): Quantification of Alizarin red S staining extracted from cell cultures for conditions shown in (C). The absorbance at 405 nm for Alizarin red S extracted from 7 to 10 replicate wells per condition is shown with the standard deviation indicated by error bars. Statistical significant reduction in Alizarin red S staining compared to Ctrl siRNA (two-tailed paired Student’s t test) is indicated by *** (p ≤ .0001), ** (p ≤ .001), or * (p ≤ .05). Abbreviations: DM, differentiation medium; GM, growth medium. Open in new tabDownload slide Identification of novel TFs regulating osteogenesis. (A): DNA motifs representing the transcription factor binding site (TFBS) for TEAD2, PURA, GTF2IRD, ZNF148, and GTF2I TFs from TRANSFAC with accession number (Acc), rate of enrichment in upregulated versus downregulated promoters (Ratio), and enrichment p-values (Wilcoxon test) (p-value). (B): Figure indicates presence (black field) or absence (white field) of ZNF148, GTF2I, PURA, GTF2IRD1, or TEAD2 TFBS at promoter regions (±2 kb of TSS) of genes with increased (log2-fold ≥1) gene expression and increased (log2-fold ≥1) H3K4m3 enrichment at Day 28 relative to Day 0. Each black or white field represents a promoter region. (C): Ingenuity pathway analysis (IPA) analysis of molecular and cellular functions for genes enriched in TEAD2 motif is shown, with the significance of association (p-value) and the number of genes (# Genes) associated with each IPA category indicated. (D): Alizarin red S staining of iMSC#3 cells at Day 28 of differentiation after siRNA-mediated knockdown of GTF2I or TEAD2 expression at Day 0, Day 7, Day 14, and Day 21 with two separate siRNAs per gene (1 or 2). Control cells were cultured in DM and exposed to a control siRNA (Ctrl) or no siRNA (DM), or cultured in GM with no siRNA (GM). A representative image for each condition is shown. (E): Quantification of Alizarin red S staining extracted from cell cultures for conditions shown in (C). The absorbance at 405 nm for Alizarin red S extracted from 7 to 10 replicate wells per condition is shown with the standard deviation indicated by error bars. Statistical significant reduction in Alizarin red S staining compared to Ctrl siRNA (two-tailed paired Student’s t test) is indicated by *** (p ≤ .0001), ** (p ≤ .001), or * (p ≤ .05). Abbreviations: DM, differentiation medium; GM, growth medium. Discussion We have used an immortalized nontumorigenic MSC cell line, iMSC#3, to investigate gene regulatory changes during osteogenic differentiation. Differentiation was accompanied by changes in gene expression associated with decrease in proliferation and upregulation of genes involved in osteoblast function as observed by others [48, 49], and in line with the notion that differentiation is associated with arrest of cell cycle progress. In addition, genes related to cell death and survival were upregulated in accordance with previous observations of apoptosis being involved in osteogenic mineralization [50]. Furthermore, upregulation of genes involved in cellular movement indicates a role for migration during in vitro osteogenesis. A fraction of these genes were bound by RUNX2 in Day 28 cells (Supporting Information Fig. S6), supporting previous findings of RUNX2 as a regulator of migration in osteogenic differentiation [19]. Upregulation or downregulation of gene expression correlated with increased or decreased enrichment of activating histone modifications, respectively, as reported by others [5, 8]. Dynamic changes in H3K27m3, however, known to occur at promoters during embryonic stem cell differentiation [8, 39], affected only a few differentially expressed genes, and thus seem to be less significant in regulation of osteoblast differentiation from MSCs. Extensive epigenetic remodeling was however found at a large number of promoters of genes with no detectable change in expression. We found increased H3K27m3 enrichment in promoters of genes involved in neurogenesis, indicating that the neurogenic potential of MSCs becomes restricted upon osteogenic differentiation. Notably, inhibition of Enhancer of Zeste 2, the methyltransferase responsible for H3K27m3 has been shown to enhance neurogenic differentiation of MSCs [51]. Conversely, genes with either a decrease of H3K27m3 or an increase of H3K27ac at their promoter regions were strongly associated with oncogenesis in the IPA analysis, indicating that histone modification changes associated with gene activation occurs at genes that are susceptible for dysregulation in cancer development. Collectively, these changes could be involved in a prepatterning mechanism poising genes for activation or repression before executing the later parts of the program. Moreover, we found that differentiation induced genome-wide reduction of histone H3K9 and H3K27 acetylation levels. Global decrease in H3K9ac at gene promoters during osteogenic differentiation from MSCs was previously described [52]. Moreover, decreased histone acetylation levels were also found in myogenic, adipogenic, astrocyte, and oligodendrocyte differentiation [7, 53, 54]. Together, these findings indicate that deacetylation at promoters and distal regulatory elements could be a common phenomenon in cellular differentiation, perhaps associated with a region-wide condensation of chromatin leaving the genome less flexible with regard to transcription and differentiation potential. The mechanism for acetylation loss could be a shift of balance between histone acetylase and histone deacetylase activity, or potentially histone cleavage [55]. Investigation of RUNX2 expression revealed transcription from both the P1 and P2 promoter at Day 28, indicating that both protein isoforms are present in differentiated cells. The isoforms give rise to two proteins which only differ by 19 amino acids in the N terminus [56]. Evidence from mutant mice suggests that the role of alternative promoter usage is to facilitate spatiotemporal regulation of RUNX2 expression throughout development rather than producing proteins with different functions [57]; however, functional differences between the isoforms have also been reported [56]. Intriguingly, differentiation is associated with a strong increase in RUNX2 transcript levels not reflected at the protein level correlating with the generation of a short transcript terminating at a novel exon upstream of the P1 promoter. This short transcript could perhaps explain the discrepancy between RUNX2 mRNA and protein levels also observed in other studies [58]. The transcript is annotated as protein-coding (ENCODE); however, any functional role remains unknown. Although predicted to be protein coding, one could speculate on the basis of the observed discrepancy between transcription and protein abundance that a potential function is independent of translation and rather involves an RNA-mediated regulatory mechanism that could affect other transcripts from the RUNX2 locus or from other genes. Notably, we also find expression of the short transcript in osteosarcomas (Lorenz et al., unpublished data). RUNX2 is a known regulator of proliferation and differentiation in osteogenesis [16, 17]. We here reveal several hundred novel RUNX2 target genes linked to these processes (Supporting Information Table S7A); our data thus strongly support previous findings and help elucidate gene regulatory networks involved in executing these aspects of osteogenic differentiation. RUNX2 target genes were also related to cellular organization, including regulation of the actin cytoskeleton, microtubules, and focal adhesion as well as migration and apoptosis (Supporting Information Table S7A), indicating that RUNX2 is involved in regulation of these processes during differentiation [19, 59], and could also possibly elucidate the oncogenic effects of RUNX2 dysregulation in cancer, and in particular metastasis [27]. Two previous studies have reported genome-wide profiling of binding sites for RUNX2 in a cancer context. First, ChIP-Seq of overexpressed and tagged RUNX2 in a prostate cancer cell line revealed predominant binding in introns and intergenic regions and less frequent binding at gene promoters, however, with target genes associated with secretory functions [60]. A second study identified 2,339 RUNX2 target genes in the SAOS-2 osteosarcoma cell line by ChIP and promoter microarrays, revealing binding to genes associated with cell adhesion and motility [61]. Although RUNX2 binding was assessed using a different experimental protocol and analysis platform, by comparing our dataset to that reported from SAOS-2 cells we identify only 243 common genes (Supporting Information Table S11). Altogether, this may point to shift in RUNX2 binding specificity due to an oncogenic context, involving loss of regulation of its normal target genes. There is emerging evidence for regulatory interplay between p53 and RUNX2. In osteoprogenitors from p53-null mice RUNX2 expression levels are elevated [62], and RUNX2 protein level was shown to be downregulated by the p53-dependent microRNA mir-34c in U2OS cells [30]. Furthermore, RUNX2 was found to interact with p53 at the BAX gene promoter to downregulate its expression and suppress apoptosis in response to DNA damage [29]. We identify 154 genes that are direct or indirect downstream targets of p53; including genes involved in apoptosis and cell cycle regulation (Supporting Information Table S8), thus indicating a widespread role for RUNX2 in the p53 regulatory network. In our analysis, PDGF BB was another predicted upstream regulator of RUNX2 target genes. PDGF BB is a known regulator of proliferation and migration of MSCs [63], and has been reported to have both a stimulatory and inhibitory effect on osteogenic differentiation [64, 65]. To our knowledge, no reports have connected PDGF signaling to RUNX2 function. Our finding of RUNX2 targeting to a significant number of known PDGF target genes could suggest such a regulation of RUNX2, perhaps associated with proliferative expansion in the early phase of differentiation. Although, binding of RUNX2 to PDGF target genes in mineralized cells is suggestive of a role in later stages of differentiation. Surprisingly, we found frequent binding of RUNX2 at or immediately downstream of 3′ gene ends, with in many cases a broad region of enrichment (Supporting Information Fig. S9). Genes with RUNX2 enrichment at both 5′ and 3′ ends were strongly associated with cellular movement, development, TGFB1, and PDGF BB signaling, whereas 5′ binding was associated with cell-to-cell signaling and 3′ binding with cell death. This could suggest that different upstream regulatory pathways influence the pattern of RUNX2 binding at target genes. Moreover, RUNX2 bound more frequently in the absence of activating histone marks at 3′ ends compared to at 5′ ends, suggesting that 3′ binding occurs in a different context. Binding near or at the 3′ end of genes has been described for other TFs and has been suggested to play a role in regulation of noncoding RNAs transcribed antisense to the corresponding protein coding gene [66]. We have preliminary strand-specific total RNA-Seq data that indicate that antisense transcription is not widely associated with RUNX2 3′ binding in iMSC#3 Day 28 cells (data not shown). Intriguingly, some genes with strong 3′ binding have miRNA genes located in close association with the enriched region, thus RUNX2 binding could be involved regulation of the miRNA rather than the protein coding gene. This was however only the case for a small number of genes. Other possible roles of RUNX2 3′ binding could be in the regulation of transcriptional termination, or involvement in looping to the promoter, as described for RNAPII at highly expressed genes [67]. More experiments are needed to determine the biological role of 3′ RUNX2 binding. We found RUNX2 protein levels to be relatively unaffected by differentiation. Thus, the osteogenic effect of RUNX2 is likely explained by binding to different genomic sites before and after induction of differentiation. However, by investigating potential differentiation-induced histone modification changes at RUNX2 bound sites we found no such correlation. Thus, in our model, binding of RUNX2 predominantly occurs at sites that did not change histone modification status between Day 0 and Day 28. Thus, targeting of RUNX2 to chromatin does not appear to be driven by dynamic change in the histone modifications investigated here. Based on regions with differentiation induced histone modification changes, we predicted five additional TFs to regulate differentiation and demonstrated by knockdown that GTF2I and TEAD2 had a negative effect on mineralization. The gene encoding GTF2I is located in a region of chromosome 7 which is deleted in Williams-Beuren syndrome (WBS) and has been proposed as candidate gene for the syndrome [68]. WBS affects several aspects of development, including cranial features and body height, indicating that bone development is affected [68]. GTF2I knockdown was reported to augment expression of the ALPL and BGLAP genes at early stages of osteoblast differentiation of mouse C2C12 myoblasts [69]. To the contrary, our findings indicate that GTF2I knockdown has a negative effect on differentiation; however, we did not investigate gene expression at the early differentiation stage. Moreover, the use of different cell types and differentiation protocols could explain the inconsistency. TEAD2 is an effector of the Hippo signaling pathway which is crucial for the regulation of organ size by balancing proliferation and differentiation [70]. TAZ, a mediator of Hippo signaling, is a key modulator of mesenchymal differentiation by activating osteogenic genes through interaction with RUNX2 [71]. Moreover, overexpression of TAZ in glioma cells induced osteogenic differentiation in a TEAD-dependent manner [72]. To our knowledge, no studies have addressed TEAD2 specifically in normal bone differentiation. Our results indicate that TEAD2 expression is required to achieve mineralization in later stages of osteogenic differentiation. Conclusion To conclude, we provide new insight into the regulatory landscape of osteoblast differentiation. We describe genome-wide chromatin remodeling, also independently of detectable changes in gene expression, affecting developmental and cancer-related genes. Furthermore, we provide new evidence to explain how RUNX2 regulates processes linked to cellular proliferation, differentiation, and cancer, by revealing genomic bindings sites and target genes in osteoblasts. Finally, based on prediction of TFBS in regions of altered histone modification status we identify and provide functional evidence for a role of two transcription factors, TEAD2 and GTF2I, in regulation of osteogenic differentiation. Acknowledgments We are grateful to The Research Council of Norway for funding this work through the Stem Cell Research Program. We thank the Genomics Core Facility at Oslo University Hospital for running DNA methylation and gene expression arrays, as well as RNA and ChIP sequencing, Meng-Yu Wang at The Department of Tumor Biology, Oslo University Hospital for preparing primary hMSC cell cultures, and Philippe Collas from the University of Oslo for scientific discussions. This work was supported by the Research Council of Norway (Grant Number 193056); the Norwegian Cancer Society (Grant Number PR-2007-0163); the University of Oslo and the charities Defence Against Cancer (“Forsvar mot kreft”), and Support group for children with cancer (“Støtteforeningen for kreftsyke barn). Author Contributions A-M.H.: conception and design, collection and/or assembly of data, data analysis and interpretation of data, and manuscript writing; J.C.B.: data analysis and interpretation and scientific discussions; K.G.H. and S.L.: collection and/or assembly of data and data analysis and interpretation of data; J.P. and T.S.M.: data analysis and interpretation; J.S.: collection and/or assembly of data; O.M.: provision of study material or patients and scientific discussions; L.A.M-Z.: conception and design, financial support, provision of study material, scientific discussions, manuscript writing, and final approval of manuscript. Disclosure of Potential Conflicts of Interest The authors indicate no potential conflicts of interest. References 1 Pittenger MF , Mackay AM, Beck SC et al. Multilineage potential of adult human mesenchymal stem cells . Science 1999 ; 284 : 143 – 147 . Google Scholar Crossref Search ADS PubMed WorldCat 2 da Silva ML , Chagastelles PC, Nardi NB. Mesenchymal stem cells reside in virtually all post-natal organs and tissues . J Cell Sci 2006 ; 119 ( Pt 11 ): 2204 – 2213 . Google Scholar PubMed OpenURL Placeholder Text WorldCat 3 Schmitt A , van GM , Imhoff AB et al. Application of stem cells in orthopedics . Stem Cells Int 2012 ; 2012 : 394962 . Google Scholar Crossref Search ADS PubMed WorldCat 4 Wagner ER , Luther G, Zhu G et al. Defective osteogenic differentiation in the development of osteosarcoma . Sarcoma 2011 ; 2011 : 325238 . Google Scholar PubMed OpenURL Placeholder Text WorldCat 5 Mikkelsen TS , Xu Z, Zhang X et al. Comparative epigenomic analysis of murine and human adipogenesis . Cell 2010 ; 143 : 156 – 169 . Google Scholar Crossref Search ADS PubMed WorldCat 6 Boyer LA , Plath K, Zeitlinger J et al. Polycomb complexes repress developmental regulators in murine embryonic stem cells . Nature 2006 ; 441 : 349 – 353 . Google Scholar Crossref Search ADS PubMed WorldCat 7 Asp P , Blum R, Vethantham V et al. Genome-wide remodeling of the epigenetic landscape during myogenic differentiation . Proc Natl Acad Sci USA 2011 ; 108 : E149 – E158 . Google Scholar Crossref Search ADS PubMed WorldCat 8 Mikkelsen TS , Ku M, Jaffe DB et al. Genome-wide maps of chromatin state in pluripotent and lineage-committed cells . Nature 2007 ; 448 : 553 – 560 . Google Scholar Crossref Search ADS PubMed WorldCat 9 Nagalakshmi U , Wang Z, Waern K et al. The transcriptional landscape of the yeast genome defined by RNA sequencing . Science 2008 ; 320 : 1344 – 1349 . Google Scholar Crossref Search ADS PubMed WorldCat 10 Wilhelm BT , Marguerat S, Watt S et al. Dynamic repertoire of a eukaryotic transcriptome surveyed at single-nucleotide resolution . Nature 2008 ; 453 : 1239 – 1243 . Google Scholar Crossref Search ADS PubMed WorldCat 11 Johnson DS , Mortazavi A, Myers RM et al. Genome-wide mapping of in vivo protein-DNA interactions . Science 2007 ; 316 : 1497 – 1502 . Google Scholar Crossref Search ADS PubMed WorldCat 12 MacQuarrie KL , Fong AP, Morse RH et al. Genome-wide transcription factor binding: Beyond direct target regulation . Trends Genet 2011 ; 27 : 141 – 148 . Google Scholar Crossref Search ADS PubMed WorldCat 13 Mundlos S , Otto F, Mundlos C et al. Mutations involving the transcription factor CBFA1 cause cleidocranial dysplasia . Cell 1997 ; 89 : 773 – 779 . Google Scholar Crossref Search ADS PubMed WorldCat 14 Komori T , Yagi H, Nomura S et al. Targeted disruption of Cbfa1 results in a complete lack of bone formation owing to maturational arrest of osteoblasts . Cell 1997 ; 89 : 755 – 764 . Google Scholar Crossref Search ADS PubMed WorldCat 15 Otto F , Thornell AP, Crompton T et al. Cbfa1, a candidate gene for cleidocranial dysplasia syndrome, is essential for osteoblast differentiation and bone development . Cell 1997 ; 89 : 765 – 771 . Google Scholar Crossref Search ADS PubMed WorldCat 16 Pratap J , Galindo M, Zaidi SK et al. Cell growth regulatory role of Runx2 during proliferative expansion of preosteoblasts . Cancer Res 2003 ; 63 : 5357 – 5362 . Google Scholar PubMed OpenURL Placeholder Text WorldCat 17 Ducy P , Zhang R, Geoffroy V et al. Osf2/Cbfa1: A transcriptional activator of osteoblast differentiation . Cell 1997 ; 89 : 747 – 754 . Google Scholar Crossref Search ADS PubMed WorldCat 18 Lee MH , Javed A, Kim HJ et al. Transient upregulation of CBFA1 in response to bone morphogenetic protein-2 and transforming growth factor beta1 in C2C12 myogenic cells coincides with suppression of the myogenic phenotype but is not sufficient for osteoblast differentiation . J Cell Biochem 1999 ; 73 : 114 – 125 . Google Scholar Crossref Search ADS PubMed WorldCat 19 Fujita T , Azuma Y, Fukuyama R et al. Runx2 induces osteoblast and chondrocyte differentiation and enhances their migration by coupling with PI3K-Akt signaling . J Cell Biol 2004 ; 166 : 85 – 95 . Google Scholar Crossref Search ADS PubMed WorldCat 20 Lee KS , Kim HJ, Li QL et al. Runx2 is a common target of transforming growth factor beta1 and bone morphogenetic protein 2, and cooperation between Runx2 and Smad5 induces osteoblast-specific gene expression in the pluripotent mesenchymal precursor cell line C2C12 . Mol Cell Biol 2000 ; 20 : 8783 – 8792 . Google Scholar Crossref Search ADS PubMed WorldCat 21 Bae JS , Gutierrez S, Narla R et al. Reconstitution of Runx2/Cbfa1-null cells identifies a requirement for BMP2 signaling through a Runx2 functional domain during osteoblast differentiation . J Cell Biochem 2007 ; 100 : 434 – 449 . Google Scholar Crossref Search ADS PubMed WorldCat 22 Reinhold MI , Naski MC. Direct interactions of Runx2 and canonical Wnt signaling induce FGF18 . J Biol Chem 2007 ; 282 : 3653 – 3663 . Google Scholar Crossref Search ADS PubMed WorldCat 23 Ling L , Dombrowski C, Foong KM et al. Synergism between Wnt3a and heparin enhances osteogenesis via a phosphoinositide 3-kinase/Akt/RUNX2 pathway . J Biol Chem 2010 ; 285 : 26233 – 26244 . Google Scholar Crossref Search ADS PubMed WorldCat 24 Stewart M , Terry A, Hu M et al. Proviral insertions induce the expression of bone-specific isoforms of PEBP2alphaA (CBFA1): Evidence for a new myc collaborating oncogene . Proc Natl Acad Sci USA 1997 ; 94 : 8646 – 8651 . Google Scholar Crossref Search ADS PubMed WorldCat 25 van der Deen M , Akech J, Wang T et al. The cancer-related Runx2 protein enhances cell growth and responses to androgen and TGFbeta in prostate cancer cells . J Cell Biochem 2010 ; 109 : 828 – 837 . Google Scholar PubMed OpenURL Placeholder Text WorldCat 26 Pande S , Browne G, Padmanabhan S et al. Oncogenic cooperation between PI3K/Akt signaling and transcription factor Runx2 promotes the invasive properties of metastatic breast cancer cells . J Cell Physiol 2013 ; 228 : 1784 – 1792 . Google Scholar Crossref Search ADS PubMed WorldCat 27 Blyth K , Vaillant F, Jenkins A et al. Runx2 in normal tissues and cancer cells: A developing story . Blood Cells Mol Dis 2010 ; 45 : 117 – 123 . Google Scholar Crossref Search ADS PubMed WorldCat 28 Won KY , Park HR, Park YK. Prognostic implication of immunohistochemical Runx2 expression in osteosarcoma . Tumori 2009 ; 95 : 311 – 316 . Google Scholar Crossref Search ADS PubMed WorldCat 29 Ozaki T , Wu D, Sugimoto H et al. Runt-related transcription factor 2 (RUNX2) inhibits p53-dependent apoptosis through the collaboration with HDAC6 in response to DNA damage . Cell Death Dis 2013 ; 4 : e610 . Google Scholar Crossref Search ADS PubMed WorldCat 30 van der Deen M , Taipaleenmaki H, Zhang Y et al. MicroRNA-34c inversely couples the biological functions of the runt-related transcription factor RUNX2 and the tumor suppressor p53 in osteosarcoma . J Biol Chem 2013 ; 288 : 21307 – 21319 . Google Scholar Crossref Search ADS PubMed WorldCat 31 Pfaffl MW . A new mathematical model for relative quantification in real-time RT-PCR . Nucleic Acids Res 2001 ; 29 : e45 . Google Scholar Crossref Search ADS PubMed WorldCat 32 Kresse SH , Rydbeck H, Skarn M et al. Integrative analysis reveals relationships of genetic and epigenetic alterations in osteosarcoma . PLoS One 2012 ; 7 : e48262 . Google Scholar Crossref Search ADS PubMed WorldCat 33 Dahl JA , Collas P. A rapid micro chromatin immunoprecipitation assay (microChIP) . Nat Protoc 2008 ; 3 : 1032 – 1045 . Google Scholar Crossref Search ADS PubMed WorldCat 34 Trapnell C , Pachter L, Salzberg SL. TopHat: Discovering splice junctions with RNA-Seq . Bioinformatics 2009 ; 25 : 1105 – 1111 . Google Scholar Crossref Search ADS PubMed WorldCat 35 Trapnell C , Williams BA, Pertea G et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation . Nat Biotechnol 2010 ; 28 : 511 – 515 . Google Scholar Crossref Search ADS PubMed WorldCat 36 Guenther MG , Levine SS, Boyer LA et al. A chromatin landmark and transcription initiation at most promoters in human cells . Cell 2007 ; 130 : 77 – 88 . Google Scholar Crossref Search ADS PubMed WorldCat 37 Heintzman ND , Hon GC, Hawkins RD et al. Histone modifications at human enhancers reflect global cell-type-specific gene expression . Nature 2009 ; 459 : 108 – 112 . Google Scholar Crossref Search ADS PubMed WorldCat 38 Bernstein BE , Kamal M, Lindblad-Toh K et al. Genomic maps and comparative analysis of histone modifications in human and mouse . Cell 2005 ; 120 : 169 – 181 . Google Scholar Crossref Search ADS PubMed WorldCat 39 Bernstein BE , Mikkelsen TS, Xie X et al. A bivalent chromatin structure marks key developmental genes in embryonic stem cells . Cell 2006 ; 125 : 315 – 326 . Google Scholar Crossref Search ADS PubMed WorldCat 40 Hassan MQ , Gordon JA, Beloti MM et al. A network connecting Runx2, SATB2, and the miR-23a∼27a∼24-2 cluster regulates the osteoblast differentiation program . Proc Natl Acad Sci USA 2010 ; 107 : 19879 – 19884 . Google Scholar Crossref Search ADS PubMed WorldCat 41 Huang S , Wang S, Bian C et al. Upregulation of miR-22 promotes osteogenic differentiation and inhibits adipogenic differentiation of human adipose tissue-derived mesenchymal stem cells by repressing HDAC6 protein expression . Stem Cells Dev 2012 ; 21 : 2531 – 2540 . Google Scholar Crossref Search ADS PubMed WorldCat 42 Ziyan W , Shuhua Y, Xiufang W et al. MicroRNA-21 is involved in osteosarcoma cell invasion and migration . Med Oncol 2011 ; 28 : 1469 – 1474 . Google Scholar Crossref Search ADS PubMed WorldCat 43 Ouyang L , Liu P, Yang S et al. A three-plasma miRNA signature serves as novel biomarkers for osteosarcoma . Med Oncol 2013 ; 30 : 340 . Google Scholar Crossref Search ADS PubMed WorldCat 44 Peng X , Guo W, Liu T et al. Identification of miRs-143 and −145 that is associated with bone metastasis of prostate cancer and involved in the regulation of EMT . PLoS One 2011 ; 6 : e20341 . Google Scholar Crossref Search ADS PubMed WorldCat 45 Liu H , Lin H, Zhang L et al. miR-145 and miR-143 regulate odontoblast differentiation through targeting Klf4 and Osx genes in a feedback loop . J Biol Chem 2013 ; 288 : 9261 – 9271 . Google Scholar Crossref Search ADS PubMed WorldCat 46 Zhao G , Cai C, Yang T et al. MicroRNA-221 induces cell survival and cisplatin resistance through PI3K/Akt pathway in human osteosarcoma . PLoS One 2013 ; 8 : e53906 . Google Scholar Crossref Search ADS PubMed WorldCat 47 Wang X , Guo B, Li Q et al. miR-214 targets ATF4 to inhibit bone formation . Nat Med 2013 ; 19 : 93 – 100 . Google Scholar Crossref Search ADS PubMed WorldCat 48 Qi H , Aguiar DJ, Williams SM et al. Identification of genes responsible for osteoblast differentiation from human mesodermal progenitor cells . Proc Natl Acad Sci USA 2003 ; 100 : 3305 – 3310 . Google Scholar Crossref Search ADS PubMed WorldCat 49 Kulterer B , Friedl G, Jandrositz A et al. Gene expression profiling of human mesenchymal stem cells derived from bone marrow during expansion and osteoblast differentiation . BMC Genomics 2007 ; 8 : 70 . Google Scholar Crossref Search ADS PubMed WorldCat 50 Lynch MP , Capparelli C, Stein JL et al. Apoptosis during bone-like tissue development in vitro . J Cell Biochem 1998 ; 68 : 31 – 49 . Google Scholar Crossref Search ADS PubMed WorldCat 51 Yu YL , Chou RH, Chen LT et al. EZH2 regulates neuronal differentiation of mesenchymal stem cells through PIP5K1C-dependent calcium signaling . J Biol Chem 2011 ; 286 : 9657 – 9667 . Google Scholar Crossref Search ADS PubMed WorldCat 52 Tan J , Lu J, Huang W et al. Genome-wide analysis of histone H3 lysine9 modifications in human mesenchymal stem cell osteogenic differentiation . PLoS One 2009 ; 4 : e6792 . Google Scholar Crossref Search ADS PubMed WorldCat 53 Takanashi M , Oikawa K, Fujita K et al. Heterochromatin protein 1gamma epigenetically regulates cell differentiation and exhibits potential as a therapeutic target for various types of cancers . Am J Pathol 2009 ; 174 : 309 – 316 . Google Scholar Crossref Search ADS PubMed WorldCat 54 Hsieh J , Nakashima K, Kuwabara T et al. Histone deacetylase inhibition-mediated neuronal differentiation of multipotent adult neural progenitor cells . Proc Natl Acad Sci USA 2004 ; 101 : 16659 – 16664 . Google Scholar Crossref Search ADS PubMed WorldCat 55 Duncan EM , Muratore-Schroeder TL, Cook RG et al. Cathepsin L proteolytically processes histone H3 during mouse embryonic stem cell differentiation . Cell 2008 ; 135 : 284 – 294 . Google Scholar Crossref Search ADS PubMed WorldCat 56 Xiao ZS , Thomas R, Hinson TK et al. Genomic structure and isoform expression of the mouse, rat and human Cbfa1/Osf2 transcription factor . Gene 1998 ; 214 : 187 – 197 . Google Scholar Crossref Search ADS PubMed WorldCat 57 Liu JC , Lengner CJ, Gaur T et al. Runx2 protein expression utilizes the Runx2 P1 promoter to establish osteoprogenitor cell number for normal bone formation . J Biol Chem 2011 ; 286 : 30057 – 30070 . Google Scholar Crossref Search ADS PubMed WorldCat 58 Pregizer S , Baniwal SK, Yan X et al. Progressive recruitment of Runx2 to genomic targets despite decreasing expression during osteoblast differentiation . J Cell Biochem 2008 ; 105 : 965 – 970 . Google Scholar Crossref Search ADS PubMed WorldCat 59 Rodriguez JP , Gonzalez M, Rios S et al. Cytoskeletal organization of human mesenchymal stem cells (MSC) changes during their osteogenic differentiation . J Cell Biochem 2004 ; 93 : 721 – 731 . Google Scholar Crossref Search ADS PubMed WorldCat 60 Little GH , Noushmehr H, Baniwal SK et al. Genome-wide Runx2 occupancy in prostate cancer cells suggests a role in regulating secretion . Nucleic Acids Res 2012 ; 40 : 3538 – 3547 . Google Scholar Crossref Search ADS PubMed WorldCat 61 van der Deen M , Akech J, Lapointe D et al. Genomic promoter occupancy of runt-related transcription factor RUNX2 in Osteosarcoma cells identifies genes involved in cell adhesion and motility . J Biol Chem 2012 ; 287 : 4503 – 4517 . Google Scholar Crossref Search ADS PubMed WorldCat 62 Lengner CJ , Steinman HA, Gagnon J et al. Osteoblast differentiation and skeletal development are regulated by Mdm2-p53 signaling . J Cell Biol 2006 ; 172 : 909 – 921 . Google Scholar Crossref Search ADS PubMed WorldCat 63 Donovan J , Abraham D, Norman J. Platelet-derived growth factor signaling in mesenchymal cells . Front Biosci (Landmark Ed) 2013 ; 18 : 106 – 119 . Google Scholar Crossref Search ADS PubMed WorldCat 64 Kubota K , Sakikawa C, Katsumata M et al. Platelet-derived growth factor BB secreted from osteoclasts acts as an osteoblastogenesis inhibitory factor . J Bone Miner Res 2002 ; 17 : 257 – 265 . Google Scholar Crossref Search ADS PubMed WorldCat 65 Pountos I , Georgouli T, Henshaw K et al. The effect of bone morphogenetic protein-2, bone morphogenetic protein-7, parathyroid hormone, and platelet-derived growth factor on the proliferation and osteogenic differentiation of mesenchymal stem cells derived from osteoporotic bone . J Orthop Trauma 2010 ; 24 : 552 – 556 . Google Scholar Crossref Search ADS PubMed WorldCat 66 Cawley S , Bekiranov S, Ng HH et al. Unbiased mapping of transcription factor binding sites along human chromosomes 21 and 22 points to widespread regulation of noncoding RNAs . Cell 2004 ; 116 : 499 – 509 . Google Scholar Crossref Search ADS PubMed WorldCat 67 Grosso AR , de Almeida SF, Braga J et al. Dynamic transitions in RNA polymerase II density profiles during transcription termination . Genome Res 2012 ; 22 : 1447 – 1456 . Google Scholar Crossref Search ADS PubMed WorldCat 68 Antonell A , Del CM, Magano LF et al. Partial 7q11.23 deletions further implicate GTF2I and GTF2IRD1 as the main genes responsible for the Williams-Beuren syndrome neurocognitive profile . J Med Genet 2010 ; 47 : 312 – 320 . Google Scholar Crossref Search ADS PubMed WorldCat 69 Lazebnik MB , Tussie-Luna MI, Hinds PW et al. Williams-Beuren syndrome-associated transcription factor TFII-I regulates osteogenic marker genes . J Biol Chem 2009 ; 284 : 36234 – 36239 . Google Scholar Crossref Search ADS PubMed WorldCat 70 Ota M , Sasaki H. Mammalian Tead proteins regulate cell proliferation and contact inhibition as transcriptional mediators of Hippo signaling . Development 2008 ; 135 : 4059 – 4069 . Google Scholar Crossref Search ADS PubMed WorldCat 71 Hong JH , Hwang ES, McManus MT et al. TAZ, a transcriptional modulator of mesenchymal stem cell differentiation . Science 2005 ; 309 : 1074 – 1078 . Google Scholar Crossref Search ADS PubMed WorldCat 72 Bhat KP , Salazar KL, Balasubramaniyan V et al. The transcriptional coactivator TAZ regulates mesenchymal differentiation in malignant glioma . Genes Dev 2011 ; 25 : 2594 – 2609 . Google Scholar Crossref Search ADS PubMed WorldCat © 2014 AlphaMed Press This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model) TI - The Regulatory Landscape of Osteogenic Differentiation JF - Stem Cells DO - 10.1002/stem.1759 DA - 2014-10-01 UR - https://www.deepdyve.com/lp/oxford-university-press/the-regulatory-landscape-of-osteogenic-differentiation-pns5RubEuH SP - 2780 EP - 2793 VL - 32 IS - 10 DP - DeepDyve ER -