Loss of Foxc1 and Foxc2 function in chondroprogenitor cells disrupts endochondral ossification

Endochondral ossification initiates the growth of the majority of the mammalian skeleton and is tightly controlled through gene regulatory networks. The forkhead box transcription factors Foxc1 and Foxc2 regulate aspects of osteoblast function in the formation of the skeleton, but their roles in chondrocytes to control endochondral ossification are less clear. Here, we demonstrate that Foxc1 expression is directly regulated by the activity of SRY (sex-determining region Y)-box 9, one of the earliest transcription factors to specify the chondrocyte lineage. Moreover, we demonstrate that elevated expression of Foxc1 promotes chondrocyte differentiation in mouse embryonic stem cells and loss of Foxc1 function inhibits chondrogenesis in vitro. Using chondrocyte-targeted deletion of Foxc1 and Foxc2 in mice, we reveal a role for these factors in chondrocyte differentiation in vivo. Loss of both Foxc1 and Foxc2 caused a general skeletal dysplasia predominantly affecting the vertebral column. The long bones of the limbs were smaller, mineralization was reduced, and organization of the growth plate was disrupted; in particular, the stacked columnar organization of the proliferative chondrocyte layer was reduced in size and cell proliferation was decreased. Differential gene expression analysis indicated disrupted expression patterns of chondrogenesis and ossification genes throughout the entire process of endochondral ossification in chondrocyte-specific Foxc1/Foxc2 KO embryos. Our results suggest that Foxc1 and Foxc2 are required for normal chondrocyte differentiation and function, as loss of both genes results in disorganization of the growth plate, reduced chondrocyte proliferation, and delays in chondrocyte hypertrophy that prevents ossification of the skeleton.

Endochondral ossification initiates the growth of the majority of the mammalian skeleton and is tightly controlled through gene regulatory networks. The forkhead box transcription factors Foxc1 and Foxc2 regulate aspects of osteoblast function in the formation of the skeleton, but their roles in chondrocytes to control endochondral ossification are less clear. Here, we demonstrate that Foxc1 expression is directly regulated by the activity of SRY (sex-determining region Y)-box 9, one of the earliest transcription factors to specify the chondrocyte lineage. Moreover, we demonstrate that elevated expression of Foxc1 promotes chondrocyte differentiation in mouse embryonic stem cells and loss of Foxc1 function inhibits chondrogenesis in vitro. Using chondrocyte-targeted deletion of Foxc1 and Foxc2 in mice, we reveal a role for these factors in chondrocyte differentiation in vivo. Loss of both Foxc1 and Foxc2 caused a general skeletal dysplasia predominantly affecting the vertebral column. The long bones of the limbs were smaller, mineralization was reduced, and organization of the growth plate was disrupted; in particular, the stacked columnar organization of the proliferative chondrocyte layer was reduced in size and cell proliferation was decreased. Differential gene expression analysis indicated disrupted expression patterns of chondrogenesis and ossification genes throughout the entire process of endochondral ossification in chondrocyte-specific Foxc1/Foxc2 KO embryos. Our results suggest that Foxc1 and Foxc2 are required for normal chondrocyte differentiation and function, as loss of both genes results in disorganization of the growth plate, reduced chondrocyte proliferation, and delays in chondrocyte hypertrophy that prevents ossification of the skeleton.
The majority of the mammalian skeleton forms through a mechanism known as endochondral ossification (1,2). In this developmental event, mesenchymal progenitor cells condense at the sites of newly forming bone and differentiate into chondrocytes. These chondrocytes undergo continued differentiation and organize themselves into cell layers that form the growth plate, which ultimately drives endochondral bone growth. Round resting zone chondrocytes form at the distal ends of long bones. These cells then differentiate inward to become flattened, highly proliferative columnar chondrocytes. Much of the extension of the bone length is achieved by the proliferative activities of these cells. Columnar chondrocytes exit the cell cycle to differentiate into prehypertrophic chondrocytes and then enlarge to form hypertrophic chondrocytes that form a mineralized matrix and sets the foundation for future bone growth. The fate of the hypertrophic chondrocytes is split into a number of outcomes. A portion of these cells will undergo apoptosis and are removed from the bone; alternatively, hypertrophic chondrocytes will transdifferentiate into bone-forming osteoblast cells that contribute to the ossified bone structure (3,4). Additional osteoblast cells originate from a cell layer, the periosteum, that lines the newly forming bone, and invade into the newly formed marrow spaced along with blood vessels.
The differentiation of chondrocytes to control growth plate functions is a tightly regulated process controlled by multiple signaling networks. In particular, cross-talk between Indian hedgehog (IHH) and parathyroid hormone-related peptide signals coordinate the proliferation of the columnar chondrocytes and their exit from the cell cycle to differentiate into prehypertrophic and hypertrophic chondrocytes (5,6). In addition, fibroblast growth factor (FGF) signaling networks also regulate growth plate chondrocyte function needed for proper bone growth (7,8). Disruption to these pathways can affect the formation of the skeleton and result in bone growth disorders in humans (9).
Campomelic dysplasia (OMIM # 114290) is a lethal skeleton malformation characterized by bowed limb bones and a reduced size of the rib cage (10). Mutations in the transcription factor gene SRY (sex-determining region Y)-box 9 (SOX9) cause campomelic dysplasia, and extensive functional analysis of SOX9 has defined it as a master regulator of the chondrocyte lineage (10)(11)(12). SOX9 function regulates multiple stages of chondrocyte differentiation and development including the initial acquisition of the chondrocyte fate, the proliferation of growth plate chondrocytes, and the transition to hypertrophic chondrocytes (12)(13)(14). In addition to its profound role in directing cells down the chondrocyte lineage, Sox9 acts to prevent differentiation toward other lineages (13,15,16). Although important in the formation of the chondrocyte lineage, SOX9 is dispensable for the initiation of the chondrogenic lineage and the induction of gene expression patterns associated with this fate (17). This finding suggests that additional transcription factors function along with SOX9 to control chondrocyte formation during endochondral ossification.
The forkhead box (FOX) transcription factors are candidates for early regulators of chondrocyte differentiation. Both Foxc1 and Foxc2 genes are expressed in the condensing mesenchyme of the presumptive endochondral skeleton (18)(19)(20)(21). Furthermore, Foxc1 and Foxc2 are required for proper endochondral ossification as mice deficient for Foxc1 or Foxc2 display disruptions to the formation of the endochondral skeleton (19,20,22). Homozygous null Foxc1 mouse mutants die shortly after birth and display small rib cages that lack an ossified sternum (19). The neural arches of the vertebral column are not fully mineralized in these mutants. The limbs are shorter in Foxc1-deficient mice and Foxc1 can regulate IHH signaling to control endochondral growth in the limb (21). Foxc2 homozygous null mutant mice also die before birth with patterning and ossification defects apparent in the axial skeleton (skull, rib cage, and vertebral column) (20). Although both Foxc1 and Foxc2 are expressed at high levels in chondrocytes of the developing limbs (21), loss of function mutation of either gene results in milder phenotypes than that observed in the axial skeleton. As FOXC1 and FOXC2 proteins have near-identical DNA-binding domains (23)(24)(25), it is possible that Foxc1 or Foxc2 may compensate for the loss of the other. Compound Foxc1 −/− ;Foxc2 −/− mice arrest in development before the onset of skeletal formation, and therefore, the combined functions of these genes in endochondral ossification is not known (26).
We wished to determine how Foxc1 and Foxc2 function in chondrocytes to regulate endochondral ossification. We first used in vitro assays to demonstrate that Foxc1 expression could be directly regulated by SOX9 activity and that gain or loss of Foxc1 function could positively or negatively regulate in vitro chondrocyte differentiation, respectively. We also demonstrate that the loss of both Foxc1 and Foxc2 function in early chondrocyte cells in the developing mouse disrupts normal endochondral ossification processes. These findings indicate that Foxc1 and Foxc2 gene function is required in the chondrocyte cells to correctly form the endochondral skeleton.

SOX9 directly regulates Foxc1 expression
Foxc1 and Foxc2 are expressed in condensing prechondrogenic mesenchyme cells at a time when SOX9 is active (18,19). Given that their mRNA expression is reduced in Sox9deficient chondrogenic tissues (17), we sought to determine whether Foxc1 and Foxc2 were directly regulated by SOX9. First, we used an inducible mouse embryonic stem (mES) cell line that contains a tetracycline-off inducible Sox9 gene (27). Upon the removal of doxycycline (DOX) for 48 h, we observed an elevation of SOX9 protein and mRNA levels (Fig. 1A) as well as an increase in collagen, type II, alpha 1 (Col2a) mRNA, a known SOX9 target of transcriptional regulation (28). Expression of Foxc1 was also elevated in response to SOX9 induction, but levels of Foxc2 mRNA remain unchanged (Fig. 1A). We next examined previously published SOX9 chromatin immunoprecipitation (ChiP)-seq data (29) and identified four SOX9-binding peaks near the mouse Foxc1 gene. Three peaks were located upstream of Foxc1, which we termed distal A (mm10 chr13:31,764,541-31,764,717), B (mm10 chr13:31,765,465-31,765,623), and C (mm10 chr13:31,779,560-31,779,803), respectively, and one peak located downstream of Foxc1, which we termed distal D (mm10 chr13:31,820,626-31,820,791). No peaks were found in proximity to the Foxc2 gene. We verified SOX9 binding to these sites using ChIP in ATDC5 cells and found SOX9 was associated with all four regulatory elements as well as the known SOX9-binding site in the intron 1 enhancer of the Col2a gene (Fig. 1B). Next, we cloned each regulatory region into a luciferase reporter that contains a basal promoter and tested for activation by SOX9 in ATDC5 cells. We found that only distal C was activated in response to SOX9. Although distal B was not activated by SOX9, we did detect elevated activity in ATDC5 cells compared with the empty reporter vector, suggesting this element may confer Sox9-independent chondrocyte regulatory activity for Foxc1 expression. Together, these findings indicate that Foxc1 is a direct target of SOX9 transcriptional regulatory activity.

Foxc1 regulates chondrocyte differentiation in vitro
Next, we examined whether Foxc1 functions in regulating chondrocyte differentiation. First, we used a Foxc1-inducible mES cell line (30) to assess whether Foxc1 overexpression influences chondrocyte differentiation. We used a chondrocyte differentiation protocol outlined in Figure 2A and described in (31). Dox was removed from mES cells 2 days before chondrocyte differentiation, and we confirmed that FOXC1 protein levels were elevated by Dox removal (Fig. 2B). Expression of SRY (sex-determining region Y)-box 6 (Sox6), Col2a, runt related transcription factor 2 (Runx2), and collagen, type X, alpha 1 (ColX) mRNA was elevated in response to FOXC1 protein induction at 21 days of differentiation compared with that of uninduced controls. These results indicate that enforced Foxc1 expression could enhance the differentiation capacity of mES cells.
To examine whether loss of Foxc1 function affected chondrocyte differentiation, we used Crispr-Cas9 to mutate the Foxc1 gene in ATDC5 cells. We generated a cell line (crisprmutated Foxc1 [crFOXC1]) that introduced a premature stop codon that would truncate FOXC1 protein after helix 2 in the forkhead domain. Given that this truncation occurred before helix 3, the DNA recognition helix, this mutation would result in a nonfunctioning protein. Reduced FOXC1 protein levels were observed in this cell line (Fig. 3A). We then induced chondrocyte differentiation by supplementing the culture media with insulin, transferrin, and selenium to initiate Foxc1 and Foxc2 regulate endochondral ossification chondrocyte differentiation. After 21 days of differentiation, we detected reduced Alcian blue-stained chondrogenic nodules in the crFOXC1 cells (Fig. 3B). Furthermore, loss of Foxc1 function reduced levels of genes expressed in early chondrocyte (Col2a, Sox9), proliferating chondrocytes (Fgfr3), prehypertrophic chondrocytes (Ihh), and hypertrophic chondrocytes (ColX and matrix metallopeptidase 13 [Mmp13]). Together, these findings indicate an impairment of chondrocyte differentiation when Foxc1 function is lost.
Foxc1 and Foxc2 are expressed in the perichondrium and the resting zone of the growth plate Next, we examined the localization of Foxc1 and Foxc2 expression in the developing skeleton. We focused on both Foxc1 and Foxc2 as each gene has a documented function for the proper formation of the skeleton (19)(20)(21). More importantly, a direct comparison of gene expression patterns of Foxc1 and Foxc2 mRNA in the developing skeleton has yet to be performed. We used dual-labeling in situ hybridization to simultaneously localize Foxc1 and Foxc2 mRNAs in the growth plate of the mouse hind limb. We observed both overlapping and distinct expression patterns for Foxc1 and Foxc2 mRNAs. Both Foxc1 and Foxc2 are expressed in the condensing mesenchyme and anlage in the developing hind limb bud at 12 days post coitum (dpc), although Foxc1 expression is more widely distributed than that of Foxc2 at this time point (Fig. 4A). At 13.5 dpc Foxc1 and Foxc2 expression becomes restricted to the peripheral region of the proximal skeletal elements (Fig. 4B). In the less-mature distal skeletal elements,  (29). C, Foxc1 distal regulatory elements were cloned into pGL4-luciferase reporters and transfected along with Sox9 in ATDC5 cells. Only Foxc1-distal C was activated by Sox9. Data presented are all data points generated from three biological replicates with each containing three technical replicates. The bottom and top boundaries of each box represent the 25th and 75th percentiles, respectively, whereas the upper and lower error bars represent the 90th and 10th percentiles, respectively. The solid bar inside the box represents the median value, whereas the dashed line represents the mean. Dots indicate each data point. Statistical analysis was performed using one-way ANOVA with Holm-Sidak pairwise multiple comparisons. *pvalue < 0.05. Col2a, collagen, type II, alpha 1; qRT-PCR, quantitative reverse transcriptase PCR; SOX9, SRY (sex-determining region Y)-box 9.
Foxc1 expression is widely expressed throughout the newly forming bone, whereas Foxc2 retains a more restricted expression pattern at this stage. In the tibia at 14.5 dpc, Foxc1 mRNA was strongly detected in the perichondrium, resting zone chondrocytes, and late hypertrophic chondrocytes/primary ossification center and lower expression detected in the proliferating and prehypertrophic chondrocytes (Fig. 4C). We also observed strong Foxc1 expression in the cells lying between the tibia and femur. Foxc2 mRNA expression was restricted to the perichondrium and newly emerging primary ossification center at this time point (Fig. 4C). A distinct expression pattern for Foxc2 mRNA was detected in the outer cell layer of perichondrium. Expression of Foxc1 and Foxc2 mRNAs in the proximal tibia at 16.5 dpc continued to be enriched in the surrounding perichondrium (Fig. 4, D-I). Foxc1and Foxc2-expressing cells were also detected in the resting zone chondrocytes (Fig. 4G) and in putative borderline chondrocytes cells found between the proliferating zone chondrocytes and the perichondrium (Fig. 4H). Very little Foxc1 or Foxc2 mRNA was detected in proliferating chondrocyte or hypertrophic chondrocyte cells (Fig. 4I). Together, these results indicate that Foxc1 and Foxc2 have restricted expression patterns in the developing bones: expression of both Foxc1 and Foxc2 is abundant in early stages of  (31). B, expression of FOXC1 was confirmed in mES cells containing a Dox-inducible flag-tagged Foxc1 expression cassette. C, chondrocyte differentiation in mES cells was monitored by measuring expression of Sox9, Sox6, Col2a, runt related transcription factor 2, and ColX expression. Data presented are from four biological replicates. The bottom and top boundaries of each box represent the 25th and 75th percentiles, respectively, whereas the upper and lower error bars represent the 90th and 10th percentiles, respectively. The solid bar inside the box represents the median value, whereas the dashed line represents the mean. Dots indicate each data point. Statistical analysis was performed using one-way ANOVA with Holm-Sidak pairwise multiple comparisons. *p-value < 0.05; *p-value < 0.05. Col2a, collagen, type II, alpha 1; ColX, collagen, type X, alpha 1; DOX, doxycycline; Sox6, SRY (sex-determining region Y)-box 6; SOX9, SRY (sex-determining region Y)-box 9.

Deletion of Foxc1 and Foxc2 in chondrocyte progenitors causes skeletal dysplasia
To understand the roles of Foxc1 and Foxc2 in the formation of the endochondral skeleton, we generated conditional, compound Foxc1 and Foxc2 mutant mice using the Col2-cre driver strain (32). We chose to delete both Foxc1 and Foxc2 genes to eliminate any potential genetic compensation that may occur when one Foxc paralog was deleted. First, we assessed the cre activity of this strain by crossing Col2-cre mice with ROSA26 mTmG mice that expresses enhanced GFP (EGFP) when cre recombinase is present. At 14.5 dpc, EGFP was detected throughout the developing skeleton. In the hind limb, cre activity was detected in the growth plate chondrocytes, perichondrium, and primary ossification center (Fig. 5, A and B). This expression pattern indicates that the Col-cre strain is active in cells of the developing skeleton and cell layers that express Foxc1 and Foxc2.

Foxc1 and Foxc2 are required for the formation of the growth plate
The limbs of Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos displayed reduced mineralization and reduction in length (Fig. 6C). We   therefore decided to investigate the formation of the growth plate of these mutants to identify the contribution of Foxc1 and Foxc2 to the function of this structure. Sections through the 16.5 dpc proximal tibia revealed growth plate anomalies in the Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (Fig. 7A). In particular, the characteristic stacked organization of the proliferative zone chondrocytes was not as prevalent in the Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos. The length of the resting zone chondrocytes layer was expanded, whereas the proliferating zone layer was reduced in the Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (Fig. 7B). We next assessed whether this reduction in size of the proliferating zone chondrocytes corresponded to a reduction in cell proliferation. As observed in Figure 7C, we did detect a reduction in the number of Ki67-positive cells in the growth plate of Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos. We also examined whether cell death was affected in Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos and found no changes in the number of apoptotic cells in the growth plates between control and Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (data not shown).
Endochondral ossification gene expression is reduced in Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos Given that Foxc1 and Foxc2 are transcription factors, it is expected that these genes exert their effects through regulation of gene expression. To further understand how Foxc1 and Foxc2 contribute to the formation of the skeleton, we analyzed gene expression of endochondral genes in Col2-cre;Foxc1 Δ/ Δ ;Foxc2 Δ/Δ embryos. We isolated RNA from the rib cage of 16.5 dpc embryos and monitored gene expression by quantitative reverse transcriptase PCR (qRT-PCR). We chose the rib cage as a tissue source as nonskeletal tissues could be efficiently dissected away from the skeletal elements. Expression of Foxc1 and Foxc2 was reduced in rib cage RNA from Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (Fig. 8A). We assayed expression of genes that act throughout all stages of endochondral ossification and found that all genes assessed had decreased expression levels when Foxc1 and Foxc2 were absent. For example, genes expressed early in the formation of chondrocytes (Sox9, Sox6, and Col2a) were reduced, and genes expressed in later chondrocyte differentiation stages such as Ihh, Fgfr3, and ColX, as well as genes expressed during osteoblast formation and mineralization (Sp7, Runx2, Col1a, and Spp1) were also reduced in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ rib cage RNA (Fig. 8A).
To gain a complete picture of gene expression changes in response to deletion of Foxc1 and Foxc2 in chondrocyte progenitors, we performed RNA-Seq from three additional samples of 16.5 rib cage RNA isolated from control and Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos. In total, we found 83 genes downregulated and 232 genes upregulated in rib cage RNA isolated from mutant embryos compared with controls (log2 foldchange +/-1; p value <0.05; false discovery rate <0.01; Tables S1 and S2). The 25 genes with the greatest reduction in expression in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (Fig. 8B) included many of the genes involved in endochondral ossification. Gene ontology analysis revealed that the majority of biological functions affected in our downregulated gene set were involved in endochondral ossification, including cartilage development, osteoblast differentiation, and ossification (Fig. 8C). Functional classification of the genes present in the upregulated dataset were enriched with those involved in lipid metabolic processes, fatty acid metabolism, and epidermis development (Table S3). Together, data from these gene expression analyses suggest that loss of Foxc1 and Foxc2 function in Col2-cre-expressing cells affected many stages of chondrocyte and osteoblast differentiation during the endochondral ossification processes.
Next, we examined the expression of chondrocyte differentiation genes in formation of the growth plate of Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos in more detail. We found that the localization of SOX9 and SOX6 protein was unaffected in the tibia growth plate. These proteins were found throughout the resting zone and proliferative zone in the proximal tibia at 16.5 dpc (Fig. 9, A and B). Expression of Fgfr3 was detected in the proliferative zone of Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (Fig. 9C) The extent of this expression domain appeared reduced in the Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos, likely because of the reduced size of the proliferating zone (Fig. 7). RUNX2 protein localized to the prehypertrophic chondrocytes of Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos (Fig. 9D), and Ihh mRNA was also expressed in this region, but expression intensity was reduced (Fig 9E). Localization of COLX protein in hypertrophic chondrocytes was also altered in the tibia growth plate of Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos at 16.5 dpc. Although the protein was present, its localization zone was expanded into the primary ossification in Col2-cre;Foxc1 Δ/ Δ ;Foxc2 Δ/Δ embryos (Fig. 9, F and H). MMP13 protein was also localized correctly to the hypertrophic chondrocytes at the osteochondral interface. Together, these expression patterns suggest that in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos, chondrocytes are able to correctly form and progress through their differentiation processes; however, this progression is altered as evidenced by expanded COLX protein.
Finally, we assessed the effect of loss of Foxc1 and Foxc2 function in the Col2-cre-expressing cells on the ossification process. Von Kossa staining in the 16.5 dpc tibia indicated that ossification occurred in the periosteum and osteochondral interface; however, large areas of the bone marrow space were unmineralized (Fig. 10, A, A 0 , B, and B 0 ). We next examined whether osteoblast formation occurred correctly in the Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos by monitoring osterix (OSX) and osteopontin (OPN) protein localization. In 16.5 dpc tibias, OSX-and OPN-positive cells were found in the osteochondral interface and in the primary ossification center (Fig. 10, C and  D), indicating that osteoblasts formed in the Col2-cre;Foxc1 Δ/ Δ ;Foxc2 Δ/Δ embryos. Osteoblasts containing COL1a protein were detected at the osteochondral interface, the primary ossification center, the periosteum, and groove of Ranvier in both control and Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos, which further indicated that osteoblast could form in the mutant bones (Fig. 10, E and F). Vascular endothelial growth factor A (VEGFA) protein localization was present but reduced in the hypertrophic chondrocytes in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ  Foxc1 and Foxc2 regulate endochondral ossification embryos (Fig. 10, G and H). Together, these findings indicate that although endochondral bone formation and mineralization is abnormal in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos, the molecular processes that regulate this bone formation do occur in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos albeit in a disrupted manner.

Discussion
Foxc1 and Foxc2 are required for normal skeletal development. We demonstrate that Foxc1 and Foxc2 are important regulators of chondrocyte formation and function required during endochondral ossification. Foxc1 expression, but not Foxc2, is directly regulated by SOX9 activity. Enforced expression of Foxc1 promoted the chondrocytic differentiation of mouse ES cells. Together, these results indicate that Foxc1 acts early in the formation of chondrogenic cells. This idea is further illustrated given that expression of Foxc1 in the developing growth plate was enriched in the resting zone chondrocytes, with less mRNA detected in more differentiated cells. Loss of both Foxc1 and Foxc2 function in early Col2-expressing chondrocytes resulted in a disruption of endochondral ossification, leading to severe skeletal hypoplasia.
Little is known about the mechanisms that regulate Foxc1 and Foxc2 gene expression in the developing endochondral skeleton. We provide evidence that SOX9 can directly regulate expression of Foxc1 in chondrocyte cells. Overexpression of Sox9 in mES cells was sufficient to induce Foxc1 mRNA. SOX9 binding to four distal regulatory regions in the Foxc1 gene had been suggested from ChIP-seq studies (29). We confirmed SOX9 association with regulatory elements in the Foxc1 gene and ascribed functional activity to these regions. Although SOX9 was able to bind to all four elements, only one element (distal C) elicited a transcription response to SOX9 in ATDC5 chondrocyte cells. A second region (distal B) did confer increased reporter activity in ATDC5, suggesting it contains regulatory information to activate Foxc1 expression in chondrocyte cells; however, its activity is likely independent of SOX9 as its activity was unaffected by SOX9 levels. Recently, expression of FOXC1 mRNA in human breast cancer cells was controlled by SOX9 (34). Regulation of Foxc1 by SOX9 is  Foxc1 and Foxc2 regulate endochondral ossification further supported by the reduction of Foxc1 expression levels in the limbs of Sox9-deficient embryos (17). This study also reported a decrease in Foxc2 mRNA in Sox9 mutant mice, although we did not detect any induction of Foxc2 expression using Sox9-inducible mES cells. It should be noted that these cells were not differentiated toward the chondrocyte lineage, and therefore, it is possible that SOX9 might regulate Foxc2 expression in a tissue-specific context. Foxc1 activation by SOX9 suggests a role for Foxc1 function in the early stages of chondrocyte differentiation. We observed that Foxc1 expression is spatially enriched in resting zone chondrocytes compared with proliferating and hypertrophic zone chondrocytes that represent later differentiation stages. In addition, overexpression of Foxc1 in mES cells led to an increased chondrocyte differentiation capacity, suggesting this early increase in Foxc1 levels promoted cells to the chondrocyte fate.
We also demonstrate the comparative expression patterns of Foxc1 and Foxc2 in the developing mouse limb. We determined that Foxc1 and Foxc2 have some overlapping as well as distinct expression patterns in the growing limb. Expression of Foxc1 and Foxc2 was readily detected in early condensing mesenchyme in the limb bud, and expression became restricted to the perichondrium surrounding the growth plate and in the resting zone chondrocytes later in as endochondral ossification proceeded. Together, these findings indicate a role for Foxc1 and Foxc2 in the initial formation of chondrocytes in the nascent skeletal element. These findings are consistent with our observation in Figure 1B, whereby activation of Foxc1 expression enhances chondrogenic gene expression in mES cells. Although much of the overall expression pattern in these regions for Foxc1 and Foxc2 was similar, there were areas where only Foxc1 or Foxc2 expression could be detected, suggesting these factors may not have completely overlapping functions. Very few Foxc1or Foxc2-expressing cells could be observed in the proliferating zone, prehypertrophic and hypertrophic chondrocytes at 16.5 dpc. This expression pattern suggests that the chondrodysplasia that we observed in the Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ mutants arose from either an indirect effect of Foxc1 and Foxc2 acting in the perichondrium to regulate chondrogenesis or as a consequence of loss of Foxc1 and Foxc2 function at an early stage of chondrocyte differentiation that disrupts cell function throughout later stages of development. We did detect Foxc1 expression in a subset of cells lying between the perichondrium and the proliferating chondrocytes, suggesting that Foxc1 may be expressed in borderline chondrocytes that supply skeletal stem cells later in life (35). Spatially, we observed increased Foxc1 and Foxc2 hybridization signals in distal skeletal elements compared with proximal ones. This expression pattern may reflect a role for Foxc1 and Foxc2 in the proximal-distal patterning of the limb. In addition, as endochondral ossification in the limb proceeds in a proximal to distal manner, this expression pattern may reflect a role for Foxc1 and Foxc2 in the early stages in the formation of skeletal elements and expression decreases as the element matures.
Global KO of either Foxc1 or Foxc2 results in a number of skeletal phenotypes in mice (19,20,22). These mutations predominantly affect the axial skeleton, although the bones of the appendicular skeleton do display a modest reduction in size. Loss of both Foxc1 and Foxc2 in Col2-cre-expressing cells results in a more severe disruption to the endochondral skeleton than loss of either Foxc1 or Foxc2 alone. We observe a general disorganization of the growth plate that affects endochondral ossification processes as indicated by the broad reduction in gene expression of regulators of these events. We reported a reduction in gene expression affecting all stages of chondrocyte differentiation and function as well as genes involved in osteoblast formation and mineralization. The growth plate of the tibia was disorganized in the Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ mutants, and the columnar arrangement of the proliferative chondrocytes did not form, resulting in decreased cell proliferation that likely accounts for the reduced bone size. We did observe a reduction in Ihh mRNA levels, which likely accounts for the impaired proliferation observed in Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ mutants. Chondroprogenitor cells lacking Sox9 display a reduction in the length of the columnar chondrocyte zone and a reduced cell proliferation (15). Moreover, Sox6 function is required to maintain columnar organization of proliferative zone chondrocytes (12). We observe similar phenotypes in Col2cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ embryos although SOX9 and SOX6 protein levels were unaffected, suggesting that Foxc1 and Foxc2 may function in common aspects of Sox9/Sox6-dependent chondrogenic processes.
We did observe some differences in the Col2-cre;Foxc1 Δ/ Δ ;Foxc2 Δ/Δ mutants compared with global Foxc1 or Foxc2 gene mutations. Germline deletion of Foxc2 results in fused ribs (20). The number and positioning of the ribs was not affected in our mutants, suggesting a role for Foxc2 in rib patterning and specification that occurs before rib chondrocyte formation. The sternum does not fully mineralize in Foxc1 −/− mice and completely lacks ossification of the xiphoid process (19). In the Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ mutants, the rib cage is reduced in size, but mineralization along the sternum appears unaffected. This difference in phenotype may arise from a lessefficient Cre activity in cells that form the sternum, or a portion of these cells that make up this structure arise from cells that do not express Col2a. COLX and MMP13 protein levels were also markedly reduced in germline Foxc1 −/− mutants (21), whereas COLX expression was expanded in the Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ mutants and MMP13 levels were not affected. These observed phenotype differences may reflect a function for Foxc1 and Foxc2 in progenitor cells before the onset of Col2a expression or in cells that do not express Col2a and thus did not delete Foxc1 and Foxc2.
The combined loss of both Foxc1 and Foxc2 in Col2aexpressing cells resulted in delayed endochondral ossification events in both the appendicular and axial skeleton. These disruptions were more pronounced in the axial skeleton. There are a number of reasons to account for such observations. First, Col2-cre is detected in the sclerotome of the somites before these cells migrating to and condensing at sites of future vertebral bone formation (36). In contrast, Col2-cre activity is not detected in the limbs until mesenchyme cells have migrated and condensed where the future bones will develop. This suggests that Foxc1 and Foxc2 are deleted at an earlier stage of vertebral bone formation that may prevent formation of adequate numbers of sclerotome cells needed to form the axial skeleton. Second, this disruption may reflect a difference in the chondrogenesis process occurring in the axial skeleton compared with the appendicular skeleton. In the sclerotome, cells must first receive Sonic Hedgehog signal from the notochord to initiate chondrogenesis, whereas such a signaling is not required in the limb skeleton (37,38). It is possible that Foxc1 and Foxc2 may function in processing this Sonic Hedgehog signal needed for chondrogenesis in the sclerotome. Third, additional transcription factors present in the limb skeleton may compensate for loss of Foxc1 and Foxc2, and such compensation does not occur in the axial skeleton. Experiments to address such mechanisms are being pursued in our laboratory.
In conclusion, Foxc1 and Foxc2 functions in chondrocytes are required for correct endochondral ossification events to occur. Loss of Foxc1 and Foxc2 function in Col2-creexpressing cells results in skeletal dysplasia and disrupts skeletal mineralization. The expression patterns for these factors in the growth plate suggest they act at early stages of chondrogenesis and their loss of function impacts later chondrogenic differentiation stages. Chondrocytes do form in the absence of Foxc1 and Foxc2, but they are unable to correctly differentiate, resulting in a disorganized growth plate, reduced chondrocyte proliferation, and delays in chondrocyte hypertrophy. Such disruptions have the overall effect of preventing correct ossification of the endochondral skeleton.
ATDC5 cells were purchased from European Collection of Authenticated Cell Cultures and cultured in DMEM:F12 containing 5% FBS. Mutagenesis of the Foxc1 gene was achieved by using Alt-R CRISPR-Cas9 system (Integrated DNA Technologies [IDT]) using a predesigned gRNA (5 0 -CAA-CATCATGACGTCGCTGC-3 0 ) that targeted the second helix of the FOXC1 FOX DNA-binding domain. Ribonucleoprotein complexes were transfected into ATDC5 cells using RNAi-MAX. Two days after transfection, cells were then plated as a single cell per well in a 96-well plate by diluting cells to approximately 80 cells/ml, and 100 μl cell solution was added to each well. Single cells were then expanded to larger culture volumes and were screened for a reduction in Foxc1 protein levels by immunoblotting with anti Foxc1 antibodies (Ori-Gene), followed by sequencing of the Foxc1 ORF. Two cell lines were selected (crFOXC1-1 and crFOXC1-8) and used for subsequent functional studies. Both lines displayed similar properties in chondrocyte differentiation experiments, but results from a single line (crFOXC1-1) are reported here.

RNA isolation and qRT-PCR
RNA was isolated from cell cultures and tissues using the RNeasy mini kit (Qiagen) following the manufacturer's protocol. All lysates were homogenized using QIAshredder (Qiagen). Tissue lysates were first disrupted with a microcentrifuge pestle before passing through a QIAshredder. RNA (500 ng) was then reverse-transcribed to cDNAs using the QuantiTect reverse transcription kit (Qiagen). qRT-PCR analysis was performed as described previously (39,40). All qRT-PCR experiments were performed with at least three biological replicates, and each contained three technical replicates. Primers used for analysis were obtained as predesigned PrimeTime qPCR Primer Assays (IDT).

Chondrocyte-differentiation procedures
For chondrocyte-differentiation experiments, mES cells were grown in DOX-free media for 2 days before induction of a chondrocyte differentiation protocol as described in ref ( (31)). Briefly, mES cells were grown as hanging drops (2500 cells/drop) to form embryoid bodies (EBs). After 2 days, EBs were pooled (24 EBs per 60-mm bacterial grade Petri dish) and grown in suspension culture containing 0.1 μM transretinoic acid. Media were replaced after 3 days with retinoic acid-free media and cultured for an additional 3 days. EBs were then transferred to gelatin-coated, tissue culture-grade, 60-mm plates containing low serum media (1% FBS) and transforming growth factor beta 3 (10 ng/ml) until differentiation was complete (21 days). RNA was isolated from 24 hanging drops/EB, and experiments were performed three times. Chondrocyte differentiation of ATDC5 cells was initiated by supplementing cell culture media with 1x insulin, transferrin, and selenium supplement (Cellgro). Differentiation media were replenished every 2 to 3 days.

Skeletal preps and analysis
Skeletons from embryos collected at 18.5 dpc were stained with Alizarin Red and Alcian Blue as described in (42). Images were captured with an Olympus E520 digital camera.

In situ hybridization
Embryos were collected at the desired stage and fixed in 4% paraformaldehyde overnight at 4 C. Tissues were washed in PBS overnight before embedding in paraffin. Sections (7 um) were collected onto Superfrost slides (Fisher Scientific). In situ hybridization was performed using the RNAscope multiplex in situ hybridization kit (Advanced Cell Diagnostics). The following RNAscope (Advanced Cell Diagnostic) probes were used: Foxc1 (412851-C2, lot 19155B), Foxc2 (406011, lot 18289A), Fgfr3 (440771, lot 18289A), and negative control probe (310043; lot 18197A). In situ hybridization experiments were performed using four different littermate pairs (mutant and control).
Paraffin sections were collected as described for in situ hybridization procedures. Slides were baked at 70 C for 1 h followed by rehydration through xylene and graded ethanol series. For SOX9, SOX6, OSX, OPN, KI67, RUNX2, and GFP antibodies, antigen retrieval was conducted by boiling slides in citrate buffer (10 mM trisodium citrate, pH 6.0; 0.05% Tween 20) for 20 min. For COL1, COL2a, COLX, MMP13, and VEGFA antibodies, sections were incubated in hyaluronidase. After antigen retrieval, slides were blocked in 5% donkey serum in PBS with 0.05% Triton X-100 for 1 h followed by incubation with the primary antibody overnight at 4 C. Slides were washed in PBS with 0.05% Triton X-100 before incubation in secondary antibodies for 1 h at room temperature, followed by staining with DAPI and mounted with coverslips using Pro-Long Gold. Immunofluorescence experiments were performed using six different littermate pairs (mutant and control).

RNA-Seq
RNA was isolated from the ribs of three control and three Col2-cre;Foxc1 Δ/Δ ;Foxc2 Δ/Δ 16.5 dpc embryos as described above. RNA sample quality control was performed using the Agilent 2100 Bioanalyzer. Samples with RNA integrity number >8 were used for library preparation following the standard protocol for the NEBNext Ultra II Stranded mRNA (New England Biolabs). Library construction and sequences was carried out at the Biomedical Research Centre Sequencing Facility (University of British Columbia, Canada). Sequencing was performed on the Illumina NextSeq 500 with Paired End 42 bp × 42 bp reads. Demultiplexed read sequences were then aligned to the reference sequence using STAR (https://www.ncbi.nlm.nih.gov/pubmed/23104886) aligners. Two different pipelines, HISAT2-featureCounts-DESeq2 and STAR-RSEM-DESeq2, were used for the downstream analysis through an inhouse script. For both the pipelines, mouse reference genome GRCm38 was downloaded from the NCBI (https://www.ncbi.nlm.nih.gov/assembly/ GCF_000001635.20/) while gene annotation was downloaded from GENCODE (https://www.gencodegenes.org/mouse/). Briefly, read alignment was performed using both HISAT2 (https://pubmed.ncbi.nlm.nih.gov/25751142/) and STAR (https://www.ncbi.nlm.nih.gov/pubmed/23104886). Quantification using featureCounts (https://pubmed.ncbi.nlm.nih. gov/24227677/) and RSEM (https://pubmed.ncbi.nlm.nih. gov/21816040/) while differential expression analysis was performed using DESeq2 (https://pubmed.ncbi.nlm.nih. gov/25516281/) for both. The final output was a list of genes with significant p-value (<0.05) after correction for multiple testing with a false discovery rate of less than 0.1. A combined list was obtained by averaging the genes with a filtered log2 fold change ±1. Only common genes with consistent differential expression log2 fold change values between both the pipelines were included in the list, whereas the inconsistent ones were excluded. Correlation for each sample between the pipelines was calculated using linear regression. RNA-Seq data are available on NCBI Gene Expression Omnibus under accession GSE165951.

Data availability
RNA-Seq results are available on NCBI Gene Expression Omnibus under accession GSE165951. All original data pertaining to this study will be made available upon request.
Supporting information-This article contains supporting information.