|Home | About | Journals | Submit | Contact Us | Français|
Sox9 encodes an essential transcriptional regulator of chondrocyte specification and differentiation. When Sox9 nuclear activity was compared with markers of chromatin organization and transcriptional activity in primary chondrocytes, we identified two distinct categories of target association. Class I sites cluster around the transcriptional start sites of highly expressed genes with no chondrocyte-specific signature. Here, Sox9 association reflects protein-protein association with basal transcriptional components. Class II sites highlight evolutionarily conserved active enhancers directing chondrocyte-related gene activity through direct binding of Sox9 dimer complexes to DNA. Sox9 binds through sites with sub-optimal binding affinity; the number and grouping of enhancers into super-enhancer clusters likely determines the levels of target gene expression. Interestingly, comparison of Sox9 action in distinct chondrocyte lineages points to similar regulatory strategies. In addition to providing insights into Sox family action, our comprehensive identification of the chondrocyte regulatory genome will facilitate study of skeletal development and human disease.
The mammalian skeleton is synthesized by cartilage-secreting chondrocytes and bone-forming osteoblasts. Different skeletal structures arise from distinct cell lineages: neural crest cells form much of the cranial vault and face; paraxial mesoderm derivatives generate additional head structures, vertebrae and ribs, while lateral plate mesoderm derivatives generate the sternum and limb skeleton (Olsen et al., 2000). Direct intramembranous ossification and cartilage-templated endochondral ossification represent alternate modes of bone formation for specific skeletal structures (Helms and Schneider, 2003; Kronenberg, 2003). In endochondral ossification, mesenchymal cells initially differentiate into mitotic chondrocytes that deposit extracellular matrix to form cartilage molds. Mitotic chondrocytes transit to postmitotic hypertrophic chondrocytes that eventually undergo cell death, leaving a matrix that is converted to bone by invading osteoblasts.
SRY-box containing gene 9 (Sox9) is a key regulator of vertebrate endochondral skeletal development (Akiyama and Lefebvre, 2011). Sox9 is initially expressed in mesenchymal condensations that differentiate into both chondrocytes and osteoblasts; early expression is essential for further development of both cartilage and bone (Akiyama et al., 2002; Akiyama et al., 2005; Bi et al., 1999). Subsequently, Sox9 expression resolves exclusively to chondrocytes; here Sox9 activity is essential for the chondrogenic program (Akiyama et al., 2002; Bi et al., 1999). In this, Sox9 is expressed at highest levels in mitotic and early pre-hypertrophic chondrocytes; Sox9 is down-regulated as chondrocytes undergo hypertrophic expansion and in joint forming regions and joint-associated articular cartilage (Akiyama et al., 2002; Dy et al., 2012; Wright et al., 1995). Recent studies have extended Sox9 action from chondrocyte specification and early chondrocyte differentiation to the initiation of chondrocyte hypertrophy (Dy et al., 2012). In summary, Sox9 functions at multiple stages of the cartilage program from mesenchymal condensation to chondrocyte hypertrophy. The importance of Sox9 levels to normal development is exemplified by the fact that haploinsufficiency for Sox9 is lethal in both mouse and man (Bi et al., 2001; Schafer et al., 1996), and causes campomelic dysplasia (Foster et al., 1994; Wagner et al., 1994).
Several cis-regulatory enhancer elements associated with direct, Sox9-dependent regulation have been identified by random and bioinformatically-driven analysis of non-coding regions flanking Sox9 (Mead et al., 2013), Ctgf (connective tissue growth factor) (Huang et al., 2010), and seven genes encoding cartilage matrix proteins: Col2α1 (Bell et al., 1997; Lefebvre et al., 1996; Leung et al., 1998; Zhou et al., 1998), Col9α1 (Genzer and Bridgewater, 2007; Zhang et al., 2003), Col10α1 (Dy et al., 2012), Col11α2 (Bridgewater et al., 1998; Bridgewater et al., 2003; Liu et al., 2000), Col27α1 (Jenkins et al., 2005), Hapln1 (Kou and Ikegawa, 2004), Matn1 (Nagy et al., 2011; Rentsendorj et al., 2005), and Acan (Han and Lefebvre, 2008; Hu et al., 2012). However, as there has been no systematic study of Sox9 actions, Sox9’s broad regulatory functions in cartilage development are not well understood.
To this end, we performed a detailed analysis of Sox9 binding, chromatin organization and transcriptional programs within mammalian chondrocytes isolated directly from the neonatal mouse rib. In addition, Sox9 programs were compared between rib and nasal chondrocytes to understand Sox9 actions in chondrocytes arising from distinct cell lineages. Our studies demonstrate different modes of Sox9 engagement at cartilage specific gene targets compared with expressed non-cartilage specific genes.
Rib chondrocytes were manually dissected from post-natal day 1 (P1) mouse to include proliferative and prehypertrophic zones and to exclude mature hypertrophic regions (Supplemental Experimental Procedures and Figures S1A–S1C). Chondrocytes were subjected to chromatin immunoprecipitation with a variety of antibodies (Table S1) followed by high-throughput sequencing (ChIP-seq); the key features of each ChIP-seq dataset are summarized in Table S2.
For the Sox9 rib ChIP-seq, 27,656 raw peaks met the peak-calling criterion. The first clear feature of Sox9 ChIP-seq peaks is a striking enrichment around the TSS; approximately 24.6% (6,794/27,656) of all peaks lie +/− 500 base pairs (bp) around the TSS even though this region represents only 0.001% of the genome (McLean et al., 2010) (Figure 1A). In addition, a large fraction of Sox9 bound regions lie at considerable distance from the TSS; about 50% of all peaks map to an interval between +/− 50 kb and 500 kb from the TSS (Figure 1B).
We first examined the characteristics of the TSS associated Sox9 dataset within the +/− 500 bp window around the TSS, hereafter referred to as Class I peaks (Sheet 1 of Table S3). The peak quality scores for Class I peaks ranked lower than non TSS associated Sox9 regions (Figure 2A): few Class I peaks were mapped within the top 2,000 peaks while many were found within the lowest ranked 2,000 peak regions (Figures 2B and 2C). These data suggest either a generally weaker direct binding interaction, or an indirect mode of interaction close to the TSS. Consistent with the latter argument, Sox9 motifs were not enriched in TSS proximal peak regions whereas enrichment was observed in Sox9-associated regions more distant from the promoter (Table S4 and Figure 2D).
We analyzed gene expression within the same region of the rib by microarray analysis and compared this data with predicted Class I peak regions. The Sox9 ChIP-seq signal intensity around the TSSs was highly correlated with the expression level of the associated genes, and with the binding of RNA polymerase II and p300, an activating component of the transcriptional complex (Figure 2E); RNA polymerase II and p300 show a similar level of co-enrichment (Figure S2A–S2C). To understand the biological importance of this peak set, we performed a GREAT Gene Ontology analysis (GREAT GO) (McLean et al., 2010). Although some chondrogenesis-related genes are present in the nearest genes to Class I peaks (Sheet 1 of Table S3), Class I targets fall within terms related to general cellular processes, not chondrocyte specific activities (Figure 2F) though comparison of MGI mouse expression data does show a significant correlation with rib gene expression (Figure 2F).
In summary, for Class I targets, Sox9 associates around the TSS of what comprises a large set of genes controlling general cellular functions not specific to chondrocytes. The level of Sox9 binding reflects the levels of core transcriptional complex engagement, and with this, the levels of expression of the associated gene.
Next, we examined the features of the remaining Class II peaks (20,862 regions, Sheet 2 of Table S3). Class II peaks showed a highly enriched recovery of Sox9 motifs, even amongst the most distant peaks sets relative to the TSS (Figure 3A), indicative of direct Sox9 binding, and a high conservation score across vertebrate species around the peak center consistent with an expected conservation of the cis-regulatory genome (Figure 3B). A clear association of Class II regions with enhancer signatures was evident; Sox9 bound regions showed an H3K4me2high/H3K4me3low enhancer signature (He et al., 2010), a strong association of p300 and RNA polymerase II consistent with active enhancers (Heintzman et al., 2007; Kim et al., 2010; Visel et al., 2009) and peaks of H3K27 acetylation (H3K27Ac) flanking Sox9, indicative of open chromatin (Heintzman et al., 2007; Rada-Iglesias et al., 2011; Zentner et al., 2011) (Figures 3C and 3D). In striking contrast to the Class I data set, GREAT GO showed a highly significant recovery of expected terms for a Sox9-regulated skeletal program (Figure 3E).
We observed substantial Sox9 ChIP-seq signal at 7 of 9 cartilage enhancers that have been identified through rigorous transgenic analysis, though the one site inAcan was below the statistical cutoff of the stringent peak calling process (Table S5). In contrast, we observed Sox9 binding at only one out of 4 regulatory elements associated with Sox9 regulation in testis cells (Table S5). Figures 3F and 3G show screen shots of ChIP-seq data around two essential cartilage matrix protein encoding genes, Acan and Col2α1.
At a higher order level, a marked clustering of Sox9 binding was evident around key chondrocyte genes expressed at high levels in chondrocytes (Figure S2D). This is supported by further analysis using the algorithm to detect so-called super-enhancers (Whyte et al., 2013). Among 13,648 adjusted rib chondrocyte enhancers in this analysis, we identified 502 super-enhancer groupings around 422 genes. Expression correlation analysis showed that expression levels of super-enhancer-associated genes were significantly higher than those of typical enhancer-associated genes (p < 2.2e-16) (Figures 3H and 3I, and Sheet 3 and 4 in Table S3), suggesting that multiple cooperative Sox9 bound enhancers may have evolved to ensure appropriate expression of chondrocyte genes.
We connected Class II Sox9 binding events to nearest genes and intersected this data set with microarray generated gene expression profiles of comparable chondrocytes. Class II peak associated-genes were positively expressed in chondrocytes relative to mouse embryo fibroblasts (MEF). Further, summing multiple sites of engagement around a putative target, it is clear that expression levels of that target increase with increasing levels of Sox9 engagement (Figure 4A).
We utilized a distance weighted regulatory potential analysis tool, S-score (Tang et al., 2011) (Supplemental Experimental Procedures), to predict direct target genes for Sox9. The analysis extracted 735 genes with high confidence (Figure 4B and Sheet 5 in Table S3), including all known Sox9 chondrocyte target genes mentioned in the Introduction. David GO analysis on this gene list showed that skeletal system development and cartilage related terms were highly enriched (Figure 4C). Together, the data are consistent with a central role of Sox9 in the activation of chondrocyte targets. Furthermore, the sum of Sox9 engagement around a gene rather than Sox9 binding at any individual site provides the strongest prediction of gene expression levels (Figures S3A and S3B).
Sox9 has been shown to activate transcription of chondrocyte specific genes through its binding to the inverted repeat of a quite variable Sox recognition sequence and that dimerization is required for activation (Bernard et al., 2003; Bridgewater et al., 2003; Mead et al., 2013; Sock et al., 2003). We performed de novo motif analysis for Class II Sox9 peak regions using two algorithms, CisGenome (Ji et al., 2008) and MEME-ChIP (Bailey et al., 2009). The primary motifs recovered were similar to a previously predicted Sox-dimeric motif, a pair of the Sox consensus sequences (A/T)(A/T)CAA(A/T)G oriented head to head and separated by 4 nucleotides (Bridgewater et al., 2003; Sock et al., 2003) (Figures 4D and S3C). As expected for direct binding of Sox9 at Class II targets, this dimeric motif is highly enriched at the predicted center of Sox9 ChIP-seq peaks (Figure 4E). A Sox-monomer motif was also recovered though these exhibited a lower enrichment score and were only weakly centered about Sox9 peaks (Figures 4F, 4G, and S3C). Thus, the Sox monomer motif does not likely represent the preferred primary site for Sox9 engagement in chondrocytes. We also recovered significant enrichment of motifs predicting engagement of AP-1, Nfat, Fox, Runx and Hox transcription factor families, consistent with the integration of multiple regulatory inputs through Sox9 directed enhancer modules (Figure S3C).
To verify Sox9-regulated enhancers from our data, we selected a set of predicted enhancers around likely Sox9 target genes, analyzing enhancer activity in vivo in a Zebrafish Enhancer Detection (ZED) system (Bessa et al., 2009). Given strong vertebrate conservation of enhancer modules, we expected that many would likely operate outside of the mouse, in a system better suited for large-scale transgenic approaches. Fourteen of 17 tested regions showed chondrocyte specific reporter gene expression (Figure 5A and Table S6). While some enhancers were active as a single copy, others required multimerization to reveal enhancer activity. Most enhancers were not active in all chondrocytes but displayed mosaic activity in individual skeletal elements. Increasing the number of copies of an individual element increased both the likelihood of enhancer activity and the level of reporter gene expression, and decreased mosaicism within the observed expression pattern (Figure S4A).
To explore Sox9 interactions in more depth, we focused on an enhancer in the intron of the cartilage matrix-encoding gene, Col9α1. Sox9 binding was associated with an evolutionary conserved motif predicted to bind Sox9 dimers within the identified chondrocyte enhancer element (Figures 5A and 5B). Enhancer-driven chondrocyte specific gene expression (Figure 5C, WT) was lost on mutation of the Sox9 site (Figure 5C, MT). Interestingly, when the endogenous motif was substituted with the optimum consensus predicted from the analysis of the entire Sox9 Class II data set, we observed stronger, less mosaic chondrocyte specific reporter activity (Figure 5C, OP).
Direct Sox9 interaction with target DNA was examined by an Electrophoresis Mobility Shift Assay (EMSA) on the predicted Sox dimeric motif for the wild-type Col9α1 enhancer (WT, Figure 5D) and one substituting the optimized Sox9 motif from whole data analysis motif recovery (OP, Figure 5D). Sox9 bound to both sequences; however, binding was markedly stronger with the optimized Sox9 target sequence (Figure 5D). As expected, enhancer binding was dependent on the HMG-box of Sox9 (Figure 5D). Mutations in either half site abolished Sox9 binding (MT1 and MT2 in Figure 5D). Further, a mutated Sox9 oligomer failed to compete with Sox9 binding to the wild-type Sox9 motif (MT, Figure 5E). As expected, the wild-type motif competed effectively but only at a high (five-fold) molar excess (WT, Figure 5E). In contrast, the optimal Sox9 site containing oligonucleotide competed effectively at equimolar levels leading to the loss of most Sox9-binding (OP, Figure 5E).
To determine whether other endogenous Sox9 sites showed sub-optimal binding relative to the recovered optimum Sox9 motif, we examined Sox9 binding to oligonucleotides incorporating 10 additional predicted endogenous Sox9 binding motifs and assayed each by EMSA. All 10 probes were bound by Sox9 protein with different levels of binding ability in the EMSA assay (Figure S4B). All were competed as effectively, or more effectively by 3-fold molar excess of optimal oligonucleotide compared with oligonucleotides for each wild-type binding site (Figure S4C).
Together these data make a strong case for a direct mode of Sox9 engagement through a dimeric binding motif. Further, the data suggest that regulation may have favored weaker binding regions than those predicted from optimal nucleotide positioning, and the reiteration of multiple weaker, cis-regulatory modules to maximize output. Consistent with this view, a statistical analysis of the occurrence of two optimal half-sites (ACAAAG, 298 in 2,474 Sox motifs) predicts 34 optimal motifs in this dataset, while the observed number is 13 (binomial p-value 4e-5), and none of these are associated with any likely chondrocyte target gene.
We sought functional roles for Class I and Class II engagements through two approaches: (1) examining gene expression changes of targets upon Sox9 manipulation and (2) intersection of human SNP data with our data set.
Regarding the first approach, we ectopically expressedSOX9 in human dermal fibroblasts (hDFs) with an adenoviral vector; hDFs are reported to acquire chondrocytic phenotypes upon the introduction of SOX9 (Ikeda et al., 2004). We confirmed no endogenous SOX9-activity in hDFs by comparing activation of a 48-bp COL2A1 intron 1 enhancer reporter construct (Kan et al., 2009) in hDFs with two Sox9+ chondrogenic cell lines (Figure 6A). Next, we infected hDFs with Sox9 containing virus, and performed RNA sequencing (RNA-seq) 2 days post infection. A box plot analysis of expression changes revealed that SOX9 overexpression induced both Class I and Class II targets (Class I targets, p<2.2e-16; Class II targets, p<2.2e-16), though as expected Class II targets showed much larger expression changes compared to Class I targets (p=3.4e-16) upon SOX9 overexpression (Figure 6B). Furthermore, the top 200 highly expressed genes within hDFs displayed a significant elevation in mRNA levels when SOX9 overexpression was compared with an control (p=1.7e-05, Figure 6C). GO analysis revealed that skeletal system development related terms were the most significantly enriched in gene sets displaying a 2-fold or greater increase on SOX9 overexpression (p= 7.2e-06, Figure 6D and Sheet 3 in Table S7). In contrast to the above gene sets comprising 235 genes (fold change rank 1 to 235), the next set of 235 genes (fold change rank 236 to 470), a low responder gene sets to SOX9 input, showed little association with skeleton related terms though ossification related terms were enriched with a low p value (p=0.048, Figure 6D and Sheet 4 in Table S7). In summary, Sox9 is likely to up-regulate a subset of Class II chondrocyte targets in hDFs. Further, Sox9 can broadly elevate transcript levels for non-chondrocyte genes consistent with a general interaction with basal transcriptional components.
To examine the relationship of Sox9 binding with 73 million human SNPs in the NCBI SNP database we mapped 523,295 SNPs to the mouse genome (mm9), and 100 SNPs were then matched to homologous mouse Sox9 binding sites in rib chondrocytes (Sheet 6 in Table S3). We analyzed the closest genes to the SNPs overlapping with Sox9 binding sites, and 8 genes (Nfatc1, Ucma, Col9α3, Itga5, Grasp, Col2α1, Glt25d2, and Mia1) were identified from the list of Sox9 putative targets identified by the S-score analysis of Class II targets. This set of SNPs may be useful foci for future exploration in human disease.
To investigate whether Sox9 regulatory networks are similar between chondrocytes derived from neural crest and mesoderm, we performed expression profiling and Sox9 ChIP-seq on cartilage of E18.5 neural crest-derived nasal septum. Nasal chondrocytes at this time are proliferative and Sox9+, Runx2− (Figures S5A and S5B).
We observed a strong correlation in gene expression between rib and nasal chondrocytes (Figure 7A). Sox9 bound Class II regions recovered from rib chondrocytes were largely bound by Sox9 in nasal chondrocytes and vise versa (Figure 7B). The major differences reflect differences in the weakest component of the dataset (Figure 7B, intensity plots); the higher ranking Sox9 peaks in rib chondrocyte data show a higher relative association with the top peaks in nasal chondrocytes than the entire Sox9 Class II peak set (Figure S5C). GREAT GO showed skeletal related terms highly enriched in nasal chondrocytes (Figure S5D). S-score analysis predicted a similar set of Sox9 target genes to those identified in rib chondrocytes (Figure S5E and Sheet 7 in Table S3) and skeletal development-related terms are significantly enriched in the target gene set (Figure S5F). These findings suggest that when Sox9 is activated in distinct skeletal lineages, Sox9 engages a similar set of target genes through a highly similar set of cis-regulatory modules.
Interestingly, despite the high correlation amongst predicted Sox9 regulated enhancers some differences are observed. We found a small number of peaks with differential signal intensities for Sox9 ChIP-seq between these two populations. To distinguish the biological significance of regions that had significantly higher Sox9 peak intensity in one population compared to the other, we performed GREAT GO analysis on such peak regions identified by peak intensity MA plots (Figure 7C). Remarkably, when nasal chondrocytes were compared to rib chondrocytes, Sox9 peak regions that had significantly higher peak intensity (log2 (nasal/rib) > 2, −log10P > 10) associated with disease terms reflecting abnormal palate and craniofacial development (Figure 7D, Blue panel; Sheet 8 in Table S3). Conversely, regions with higher Sox9 binding intensity in rib than in nasal chondrocytes (log2(rib/nasal) > 2, −log10P > 10) were associated with genes whose mutation caused abnormal skeletal or cartilage development (Figure 7D, Red panel; Sheet 9 in Table S3). For example, among genes associated with each term in Figure 7D (Sheet 10 in Table S3), the Eya1 gene showed a higher expression level and stronger Sox9 association in nasal chondrocytes whereas Hoxc8 expression and Sox9 association in flanking regions were higher in rib chondrocytes (Figure 7E and Sheet 11 in Table S3). In summary, general Sox9 binding profiles and identified Sox9 targets are similar between rib and nasal chondrocytes though variation in the levels of Sox9 binding, and a small number of qualitative differences, may play a significant role in regulatory outcomes within distinct chondrocyte populations.
Whereas a significant body of information has accrued on epigenetic control of mammalian gene activity from systematic study of some well-studied cell types in culture, much less is known about relevant cell types in vivo (Consortium et al., 2007; Maher, 2012; Rosenbloom et al., 2010). The mammalian skeleton shows a relatively simple, spatially defined differentiation program where cell-lineage and genetic studies have established cellular and regulatory hierarchies. Chondrocytes and osteoblasts arise from multiple lineages where a small number of transcriptional regulators play a predominant role in specifying skeletal cell types (Akiyama et al., 2002; Bi et al., 1999; Dy et al., 2012; Komori et al., 1997). In particular, Sox9 has emerged as a key regulator of chondrogenic programs (Akiyama et al., 2002; Bi et al., 1999; Dy et al., 2012) and Runx2 and Sp7/Osx critical determinants of osteoblast development (Komori et al., 1997; Nakashima et al., 2002). The skeletal system and these cell-type specific regulatory factors are attractive targets for developing a deep understanding of regulatory interactions shaping the mammalian body plan.
Our studies provide many important insights into the Sox9-directed process of chondrogenesis. The analysis of Class I peaks shows broad, non-specific engagement of Sox9 at promoter regions within chondrocytes that is likely mediated through protein-protein interactions with the basal transcription apparatus. Ectopic expression studies suggest that Sox9 engagement may elevate overall expression levels. The analysis of Sox9 at non-promoter associated Class II peaks indicates direct DNA binding through dimeric Sox9 recognition within evolutionarily conserved enhancer elements. These enhancers cluster around highly expressed cartilage-related genes, multiple enhancers displaying active chromatin signatures extending for many 10’s to 100’s of kbs from the TSS of the target gene. In vitro analysis of recovered Sox9 motifs and in vivo enhancer organization argues for a mode of regulation favoring the use of multiple enhancers each with a sub-optimal Sox9 binding ability that form super-enhancer-like groupings, which likely underpins the high levels of gene expression observed for the many chondrocyte matrix encoding genes. Interestingly, Sox9 utilizes a qualitatively similar enhancer set regardless of chondrocyte origin though where there are key genes that distinguish cranial and axial chondrocytes, notably Hox genes in the latter, novel patterns of Sox9 engagement are observed in the chondrocyte lineages.
These data provide a comprehensive, genome-scale assessment of Sox9’s action within the nucleus of the normal chondrocyte. A limited analysis of Sox9 engagement around promoter regions of a selected sub-set of targets has been reported from the study of a rat chondrosarcoma line (Oh et al., 2010). These studies were recently extended by RNA-seq and Sox9 ChIP-seq and included some analysis of primary chondrocytes (Oh et al., 2014). Unfortunately, we were unable to obtain the primary data underpinning these studies so we have been unable to corroborate, compare or incorporate the findings of these authors into our work.
The two classes of Sox9 binding on the genome point to distinct function of Sox9 in the regulation of gene activity in chondrocytes. Class I binding correlates strikingly, both quantitatively and spatially, with engagement of p300, RNA polymerase II and presumably a spectrum of other components linked to the basal transcriptional apparatus. Further, the absence of enrichment of Sox9 binding motifs in this data set suggests an indirect mode of engagement with one or more of the above. Sox9 is known to associate with p300; consequently, their interaction is likely to account for some, and possibly all, of the Class I data (Furumatsu et al., 2005; Tsuda et al., 2003). This data, and analysis of gene expression following ectopic Sox9 expression in fibroblasts, suggests a model wherein Sox9 engagement with the basal transcription apparatus generally elevates gene expression. As the levels of Sox9 binding correlate with transcriptional activity rather than cell type specificity of the target gene, it is likely that Sox9’s action on Class I targets would be most significant for highly expressed genes, including house-keeping genes, independent of any particular function within the cell.
Class II binding associates with broad range of skeletal enhancers; the major cartilage matrix producing genes are prominent targets of the enhancer network. Amongst this group, we confirmed Sox9 engagement in most established enhancers previously linked to Sox9’s regulation of early chondrocyte differentiation (Akiyama et al., 2002; Akiyama and Lefebvre, 2011; Bell et al., 1997), and Sox9’s auto-regulation (Mead et al., 2013); Sox9 appears as a major target of its own regulation in our data with discrete binding sites extending for 100’s of kilobases. Two other Sox factors, Sox5 and Sox6 act downstream of Sox9 (Akiyama et al., 2002) and are also predicted to be direct targets of Sox9 regulation in the current data.
More contentious is Sox9’s function in maturation of hypertrophic chondrocytes. Initial studies based on mouse genetics and in vitro studies suggested Sox9 inhibits the hypertrophic program (Akiyama et al., 2002; Ikeda et al., 2004). Further, Sox9 is reported to repress transcription of Col10a1, a key gene in the production of a collagen-type restricted to post-mitotic and hypertrophic chondrocyte fates. Here, repression was reported to map to a Sox motif upstream of the Col10a1 gene (Leung et al., 2011). However, in contrast to these findings, recent studies disrupting Sox9 within the growth plate suggest Sox9 promotes hypertrophy, acting cooperatively with Mef2c, a master regulator of hypertrophy, to promote transcription of Col10a1 (Dy et al., 2012). We observe a weak Sox9 peak adjacent to the region identified upstream of Col10α1 but no Sox9 binding to the regions around Mef2c. Our ChIP population is predominantly comprised of mitotic chondrocytes, the prehypertrophic zone is a small fraction and hypertrophic chondrocytes are absent. Clearly, Sox9 protein levels remain elevated within Col10α1expressing hypertrophic chondrocytes so Sox9 may play a continuing role in maturing chondrocytes. Thus, given potential heterogeneity of our ChIP population, future ChIP-seq studies that can overcome the current barrier to extensive dissection of primary cell types and the requirement for substantial amounts of nuclear material to generate comparable robust data to that presented here, will help in clarifying Sox9 action in hypertrophic chondrocyte development.
Sox9 acts predominantly as a homodimer in chondrocytes (Bridgewater et al., 2003; Coustry et al., 2010; Genzer and Bridgewater, 2007; Jenkins et al., 2005). Our findings lend strong support to Sox9 dimers as the key regulatory input into chondrocyte cis-regulatory modules, which is consistent with the identification of mutations that disrupt the dimerization capability of Sox9 in patients with campomelic dysplasia (Bernard et al., 2003; Sock et al., 2003). Biochemical and statistical analysis of Sox9 motifs indicates that enhancers favor lower affinity Sox9 sites compared with the highly symmetric homodimer motif that represents the motif preference from the entire Sox9 bound enhancer dataset. The reason for this is not clear but lower affinity interactions may favor cooperative engagement at regulatory regions. As an example, Sox5 and Sox6, members of the distinct SoxD sub-family are known to act downstream of Sox9 in chondrocyte regulation. While these members can dimerize with each other, they do not dimerize with Sox9, but engage with Sox9 in the cooperative regulation of cartilage enhancers for Col2α1 and Acan (Han and Lefebvre, 2008; Lefebvre, 2010; Lefebvre et al., 1998).
In examining the enhancer landscape around target genes, and the results of enhancer analysis in transgenic studies, collectively our data suggest that the transcriptional outcome in terms of the probability of expression and the level of expression of a given target gene is likely determined through the use of multiple enhancer elements. Though our transgenic studies are limited to facial structures in the zebrafish all enhancers have broad, general cartilage specificity. Thus, at this level of analysis, enhancers are likely “level-determining modules” that sum to give appropriate levels of expression of the skeletal target gene rather than “spatial-determining modules” each acting in a unique range of skeletal elements that sum to cover the entire skeleton.
Sox9 is a master regulator of male sex determination and skeletal development. Further, Sox9 action is critical in different aspects of lung, hair, kidney, ear, hair and gut development (reviewed in Lefebvre et al., 2007). How are these specific outcomes realized? A study of Sox9 interactions in male sex determination focused on ChIP-chip analysis of promoter regions (Bhandari et al., 2012). However, on the basis of data here, this will likely have selected for Class I interactions and consistent with this conclusion motif recovery did not provide compelling evidence for a Sox9 binding motif. A comparative analysis of the limited set of cis-regulatory modules engaged by Sox9 in the testis that include Sox9-itself (Sekido and Lovell-Badge, 2008), miR202 (Wainwright et al., 2013), and Amh (Arango et al., 1999; De Santa Barbara et al., 1998), none of these enhancers were bound by Sox9 in chondrocytes.
Interestingly, reporter constructs comprising 5 copies of optimized Sox9 dimer motifs show no chondrocyte activity in the zebrafish assay. This may reflect accessibility of the reporter construct if for example specific histone modification of an enhancer is required to enable accessibility of Sox9 (Zhou et al., 2011). Indeed, the testis enhancers linked to Sox9’s testis regulatory program did not display poised or active histone marks (H3K4me2 or H3K27acetyation, respectively) in chondrocytes (data not shown). Sox9 binding may also be insufficient for chondrocyte expression (data not shown). Cooperative partnerships between Sox9 and Sox5/6 discussed earlier may be one important factor in determining a specific regulatory outcome (Han and Lefebvre, 2008; Lefebvre, 2010; Lefebvre et al., 1998). Cooperative interaction with non Sox-family members is also suggested by our study. Motifs for AP-1, NFAT, Fox, Runx, and Hox transcription factor families were all significantly enriched in the Sox9 ChIP-seq datasets. Given that several of these families have been implicated in formation and maintenance of cartilage (Behrens et al., 2003; Gross et al., 2012; Karreth et al., 2004; Kimura et al., 2010; Nifuji et al., 2001; Rodova et al., 2011; Wang et al., 2009; Yoshida et al., 2004; Yueh et al., 1998), their actions may modulate Sox9 engagement, and cooperate with Sox9 to regulate enhancer activity. The identification and validation of the core Sox9 network provides a regulatory framework for extending our understanding of the transcriptional mechanisms at play in mammalian skeletogenesis.
Rib chondrocytes isolated from newborn mice were immediately subjected to chromatin preparation and ChIP as described (Odom et al., 2004; Vokes et al., 2007). ChIP-seq libraries were constructed using ChIP-seq DNA Sample Prep Kit (IP-102-1001; Illumina Inc., San Diego, CA). Sequencing was performed on Genome Analyzer II (Illumina).
RNA was isolated form rib or nasal chondrocytes using TRIZOL reagent (15596-026; Life Technologies) and RNeasy Mini Kit (74104; QIAGEN) according to manufacturers’ instruction. Microarray analysis were performed on GeneChip mouse Gene 1.0 ST array (901169; Affymetrix, Santa Clara, CA). RNA-seq analysis was performed on the Illumina platform with Nextseq 500.
We thank Drs. Henry M. Kronenberg, Clifford J. Tabin, Ung-il Chung, Kevin Peterson, Xiaoxiao Zhang, and Anton Valouev for their helpful inputs; Laurie Chen, Joe Vaughan, Jill McMahon, Christian Daly, Jennifer Couget, Toshiyuki Ikeda, Taku Saito, Fumiko Yano, Rie Yonemoto, and Nozomi Nagumo for providing technical assistance; and Drs. Peng Huang, Alex Schier (Harvard University), Gage Crump (University of Southern California) and his lab members for help with the zebrafish work. Work in A.P.M.’s laboratory was supported by a grant from the NIH (DK056246). Work in the University of Tokyo was funded by Grant-in-Aid for Scientific Research (23689079, 26713054, 15K15732, 24240069, and 26221311) from the JSPS, Takeda Science Foundation Research Grant, the Graduate Program for Leaders in Life Innovation, and the Core-to-Core Program A. Advanced Research Networks.
Microarray, ChIP-seq, and RNA-seq data reported in this paper are available in the Gene Expression Omnibus (GEO) under accession numbers GSE69108, GSE69109, and GSE 69110, respectively.
AUTHOR CONTRIBUTIONSS.O. conceived and performed experiments on ChIP-seq and expression profiling, analyzed and interpreted data, and wrote the manuscript. X.H. conceived and performed zebrafish and biochemical experiments, analyzed and interpreted data, and wrote the manuscript. H.H. analyzed and interpreted data. A.P.M. conceived and designed the project, interpreted data, and wrote the manuscript. S.O. and X.H. contributed equally to this work.
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.