The regulated expansion of chondrocytes within growth plates and joints ensures proper skeletal development through adulthood. Mutations in the transcription factor NKX3.2 underlie spondylo-megaepiphyseal-metaphyseal dysplasia (SMMD), which is characterized by skeletal defects including scoliosis, large epiphyses, wide growth plates and supernumerary distal limb joints. Whereas nkx3.2 knockdown zebrafish and mouse Nkx3.2 mutants display embryonic lethal jaw joint fusions and skeletal reductions, respectively, they lack the skeletal overgrowth seen in SMMD patients. Here, we report adult viable nkx3.2 mutant zebrafish displaying cartilage overgrowth in place of a missing jaw joint, as well as severe dysmorphologies of the facial skeleton, skullcap and spine. In contrast, cartilage overgrowth and scoliosis are absent in rare viable nkx3.2 knockdown animals that lack jaw joints, supporting post-embryonic roles for Nkx3.2. Single-cell RNA-sequencing and in vivo validation reveal increased proliferation and upregulation of stress-induced pathways, including prostaglandin synthases, in mutant chondrocytes. By generating a zebrafish model for the skeletal overgrowth defects of SMMD, we reveal post-embryonic roles for Nkx3.2 in dampening proliferation and buffering the stress response in joint-associated chondrocytes.
In much of the developing vertebrate body, cartilage is a transient tissue that progressively remodels to bone through a process known as endochondral ossification. In other contexts, such as cartilage of the nose and ear and the articular cartilage lining the bony surfaces of healthy joints, cartilage is permanent. The growth and function of both transient and permanent cartilage relies on the stratification of chondrocytes into distinct zones. In the growth plates of endochondral bones, chondrocytes are arranged into a resting zone containing stem cells (Mizuhashi et al., 2018; Newton et al., 2019), a zone of round proliferative chondrocytes, a zone of proliferative flattened chondrocytes that merges into a pre-hypertrophic zone, a hypertrophic zone in which chondrocytes enlarge and calcify, and a transitional zone in which chondrocytes undergo apoptosis or transdifferentiate into osteoblasts (Giovannone et al., 2019; Jing et al., 2015; Kronenberg, 2003; Yang et al., 2014; Zhou et al., 2014). Joint cartilage is also zonally arranged: superficial chondrocytes of a flattened morphology line the synovial cavity and produce specialized molecules such as Prg4 (lubricin) (Askary et al., 2016; Kozhemyakina et al., 2015), with deeper chondrocytes arranged in zones reminiscent of growth plates and transitioning into the underlying subchondral bone (Lui et al., 2015). Defects in the specification and maintenance of these cartilage zones can lead to dwarfism and other skeletal dysplasias in the context of the growth plate, and arthritis in the context of the joint. The mechanisms for maintaining the correct proportions and identities of cartilage layers at either growth plates or joints remain, however, incompletely understood. In particular, we still know little about how the proliferative expansion of cartilage is zonally regulated to meet the differing demands of endochondral bone and joint growth.
Humans with spondylo-megaepiphyseal-metaphyseal dysplasia (SMMD; OMIM #613330), a rare skeletal dysplasia linked to homozygous frameshift mutations in the homeobox transcription factor NKX3.2 (also called BAPX1 or NKX3-2), display skeletal overgrowth in the wrists and digits that is accompanied by supernumerary bones (pseudoepiphyses) as well as scoliosis (Hellemans et al., 2009). In these patients, the cartilaginous growth plates of long bones are abnormally wide. Whereas work in animal models has revealed requirements for Nkx3.2 in regulating chondrocyte development within both the growth plate and joints, the skeletal overgrowth seen in SMMD has not been observed. In chick and mouse, Nkx3.2 is expressed in the proliferating and pre-hypertrophic chondrocytes of the growth plate but largely excluded from the less proliferative hypertrophic zone and articular surfaces of joints (Church et al., 2005; Provot et al., 2006; Tribioli et al., 1997). Nkx3.2−/− mice have defects in their vertebrae and cranial skeleton, which have been attributed to decreased cartilage formation (Akazawa et al., 2000; Tribioli and Lufkin, 1999) and is exacerbated by further mutation of the paralog Nkx3.1 (Herbrand et al., 2002). However, despite expression of Nkx3.2 in the long bones of the developing limbs, Nkx3.2−/− mice do not display the elongated wrists and digits of SMMD patients at birth. As these mice die shortly after birth, possibly because of asplenia, the lack of limb overgrowth phenotypes could reflect postnatal roles for Nkx3.2 in restricting cartilage growth. Nkx3.2 is thought to form a positive regulatory loop with Sox9 (Yamashita et al., 2009; Zeng et al., 2002), a master regulator of chondrogenesis (Bi et al., 1999), consistent with their co-expression in growth plate chondrocytes (Provot et al., 2006). Studies in the chick limb and mammalian cell culture have suggested that Nkx3.2 may delay hypertrophic maturation of growth plate chondrocytes through repression of Col10a1 (Provot et al., 2006) and Runx2 (Kawato et al., 2011; Yamashita et al., 2009), a transcription factor required for the differentiation of cells into mineral-producing osteoblasts and hypertrophic chondrocytes (Takeda et al., 2001).
In addition to roles in growth plate regulation, work in zebrafish has revealed essential roles for Nkx3.2 in formation of the jaw joint (Miller et al., 2003). Inhibition of Nkx3.2 function with an antisense morpholino results in a complete loss of the zebrafish jaw joint, including the bony retroarticular process that is a site of insertion for the jaw adductor (adductor mandibulae). The role of Nkx3.2 in specifying the jaw joint appears to be shared between zebrafish and amphibians but not mammals, potentially due to positive regulation of joint-promoting Gdf5/6 family members only in non-mammalian species (Lukas and Olsson, 2018). Although the malleus-incus joint (the evolutionary homolog of the zebrafish jaw joint) or the temporomandibular joint (the mammalian jaw joint) are normal in Nkx3.2−/− mice, the malleus is narrower and the associated gonial bone is lost (Tucker et al., 2004).
Here, we characterize a zebrafish nkx3.2 mutant that lacks the jaw joint (see also Miyashita et al., 2020 for a detailed analysis of adult cranial defects). Similar to SMMD patients, mutants are adult viable and develop cartilage outgrowths, albeit in the jaw joint region rather than wrists and digits, as well as scoliosis. By comparing mutants to animals with transient, early knockdown of nkx3.2 function, we find that Nkx3.2 functions at post-embryonic stages to prevent cartilage overgrowth and pattern the spine. Single-cell RNA-seq of facial chondrocytes reveals a mutant-enriched population of chondrocytes with high levels of cell cycle and stress-induced pathway genes, including the prostaglandin D2 synthase gene ptgdsb.1 and the mTORC1 pathway genes sestrin1 (sesn1) (Budanov and Karin, 2008) and eukaryotic translation initiation factor 4E binding protein 3 (eif4ebp3) (Tsukumo et al., 2016; Yogev et al., 2013). We validate that, in place of the missing jaw joint, adult nkx3.2 mutants develop abnormally proliferative chondrocytes with elevated levels of ptgdsb.1 and sesn1, yet do not upregulate genes associated with hypertrophic maturation. Rather than suppressing hypertrophic maturation, our results indicate post-embryonic roles for Nkx3.2 in restricting chondrocyte proliferation, possibly linked to prostaglandin and mTORC1 regulation, to ensure proper skeletal proportions, including in the limbs of SMMD patients.
Mutant nkx3.2 zebrafish lack a jaw joint and display skull and spine defects as adults
Given concerns over efficacy and specificity of morpholino knockdown (Eve et al., 2017; Gentsch et al., 2018; Joris et al., 2017; Kok et al., 2015; Law and Sargent, 2014), we first aimed to confirm requirements for nkx3.2 in zebrafish jaw joint development (Miller et al., 2003) using genetic mutants. To do so, we generated two mutant alleles, nkx3.2el802 and nkx3.2ua5011 using TALEN (Huang et al., 2011; Sander et al., 2011) and CRISPR (Hwang et al., 2013) mutagenesis, respectively (Fig. 1A). nkx3.2el802 contains a 14 bp deletion that removes the ‘A’ of the start codon and is predicted to prevent translation of nkx3.2 mRNA, and nkx3.2ua5011 introduces a 20 bp deletion and frameshift mutation predicted to disrupt the DNA-binding homeobox domain. These are both predicted to be severe loss-of-function alleles; lack of a zebrafish-reactive Nkx3.2 antibody did not allow us to test this directly.
Consistent with the nkx3.2 morpholino knockdown phenotype (Miller et al., 2003), animals homozygous for either nkx3.2 mutant allele fail to develop a jaw joint, as assessed by Alcian Blue labeling of cartilage at 7 days post-fertilization (dpf) (n=3/3; Fig. 1B; Miyashita et al., 2020). In wild-type animals at 11 dpf, trps1:GFP marks articular chondrocytes that are flanked by sox10:dsRed+ chondrocytes of the upper and lower jaw cartilages. Live imaging of nkx3.2 mutant embryos revealed loss of trps1:GFP expression and a continuous field of sox10:dsRed+ chondrocytes across the presumptive jaw joint region (n=5/5; Fig. 1C). In contrast to the vast majority of nkx3.2 morpholino-treated embryos that die by 7 dpf, zebrafish homozygous for either mutant allele are viable and fertile into adulthood, especially when raised separately from their wild-type siblings, and display a fixed open jaw (n=6/6; Fig. 1D; Miyashita et al., 2020).
Adult nkx3.2 mutant zebrafish display a number of other cranial skeletal defects beyond jaw joint loss, as well as vertebral abnormalities. Mutants display altered orientation of facial bones (e.g. ceratohyal), possibly due to mechanical restriction from jaw joint fusion (n=6/6; Fig. 1D). The paired pterosphenoid and orbitosphenoid bones are articulating endochondral bones of the neurocranium (Cubbage and Mabee, 1996) and these are abnormally fused in nkx3.2 mutants (n=6/6; Fig. 1E). Micro-computed tomography (µCT) imaging further reveals that 11/11 nkx3.2 mutants have abnormal ossification and a shortening of the rostral spine (Fig. 2B,F), misalignment of ribs (Fig. 2C,G) and rotation and fusion of the caudal vertebrae (Fig. 2D,H) at 60 dpf, with spinal defects apparent as early as 30 dpf (n=11/11; Fig. S1). Overall, the combined rostral defects and abnormal curvature are associated with shortening of the spine (Fig. 2A,E). Spinal deformities are consistent with the reported expression of nkx3.2 in the neural and hemal arches that contribute to vertebrae (Crotwell et al., 2007), as well as the scoliosis and shortening of the neck in SMMD patients.
Ectopic jaw cartilage overgrowth in nkx3.2 mutants
We next examined the long-term fate of nkx3.2 mutant cells that would have normally formed the jaw joint. In whole-mount preparations of the 60 dpf skull stained with Alcian Blue (cartilage) and Alizarin Red (bone), we observed greatly increased Alcian-positive tissue in the region where the jaw joint would have been (n=6/6; Fig. 1D). This cartilage overgrowth was particularly evident upon dissection of the presumptive jaw joint tissue (insets in Fig. 3A,B). Hematoxylin and Eosin (H&E) staining of frontal sections through the presumptive jaw joint region further demonstrated a dramatic expansion of cartilage tissue in mutants (n=3/3; Fig. 3C,D). We next assessed proliferation in the presumptive jaw joint region at 21 dpf by subjecting animals to a 1 h pulse of bromodeoxyuridine (BrdU), which is incorporated into newly replicating DNA (Fig. 3E,F). In 3/3 wild types, we observed very few BrdU+ cells near the jaw joint, and sparse BrdU+ chondrocytes in the ceratohyal growth plate. In 3/3 mutants, we observed increased numbers of BrdU+ chondrocytes in place of the missing jaw joint, consistent with cartilage overgrowth. To validate the BrdU incorporation analysis, we also assessed expression of Proliferating Cell Nuclear Antigen (PCNA), a marker for DNA replication, in chondrocytes (Fig. 3G-J). Around wild-type joints, we observed a few PCNA+ cells within the lower jaw Meckel's cartilage, but almost no PCNA+ cells at the joint articular surface. PCNA+ chondrocytes were more abundant in the proliferative zone of the ceratohyal growth plate. In nkx3.2 mutants, numerous PCNA+ chondrocytes are evident in the fused jaw joint region. Quantification of the percentage of chondrocytes expressing PCNA (Fig. 3K) revealed a dramatic increase in proliferative chondrocytes at the nkx3.2 mutant jaw joint region (identified by nkx3.2 expression in adjacent sections) compared with the wild-type jaw joint (n=3 per genotype, P=0.009), similar to the levels of proliferative chondrocytes within the wild-type ceratohyal growth plate. We also observed a modest increase in proliferation rates within the mutant ceratohyal growth plate (n=3 per genotype, P=0.07). These findings indicate that chondrocytes in the presumptive jaw joint region, and likely to a lesser extent in the ceratohyal growth plate, become increasingly proliferative in nkx3.2 mutants, likely contributing to the cartilage overgrowth seen in juvenile and adult mutants.
Nkx3.2 function is required in post-embryonic stages to suppress cartilage overgrowth
We next sought to determine whether the cartilage overgrowth seen in adult nkx3.2 mutants was solely the consequence of not forming the jaw joint, or alternatively reflected a continued requirement for nkx3.2 in chondrocyte regulation. To uncouple early versus continual requirements for nkx3.2, we took advantage of the transient nature of antisense morpholino knockdown – dilution during embryonic cell divisions is thought to limit the window of efficacy of most morpholinos to the first few days of embryogenesis (Eisen and Smith, 2008). As described earlier, the jaw joint loss seen in both nkx3.2 mutant alleles is nearly identical to that reported for nkx3.2 morpholino knockdown (Miller et al., 2003). We therefore injected one-cell-stage embryos with the same dose of nkx3.2 morpholino reported by Miller et al. (2003). Using sox10:dsRed to visualize chondrocytes in live animals, we confirmed loss of the jaw joint in the majority of morpholino-injected animals at 7 dpf and only raised animals lacking jaw joints to adulthood. From 225 animals lacking the embryonic jaw joint, only four inflated their swim bladders and survived until adulthood (n=2 at 60 dpf, n=2 at 90 dpf), in marked contrast to the viability of jaw-joint-less nkx3.2 mutants. The four viable morpholino-injected adults exhibited a similar open gape of the jaw as seen in nkx3.2 mutants (Fig. 4A-F). In wild types, permanent Meckel's cartilage within the lower jaw is separated from the cartilaginous growth plate of the upper jaw palatoquadrate bone by the jaw joint and associated subchondral bone (Fig. 4D,G). In mutants, we observed fusion of lower jaw cartilage to the palatoquadrate cartilage growth plate (Fig. 4H; see also Fig. 3D). In contrast, in morpholino-injected animals, bone continues to separate lower jaw cartilage and the palatoquadrate cartilage growth plate despite the lack of a jaw joint, and no cartilage overgrowth was observed in any of the four adult-viable animals (Fig. 4I). Whereas 5/6 nkx3.2el802 mutant animals displayed severe caudal tail curvature defects (Fig. 4K), all four morpholino-injected adults had normal spines and tails (Fig. 4J-L). These findings indicate post-embryonic roles for Nkx3.2 in restricting joint-associated cartilage growth and ensuring normal spine development.
Single-cell RNA-seq reveals upregulated mitotic and stress response genes in nkx3.2 mutant chondrocytes
In order to understand how Nkx3.2 might regulate subtypes of joint-associated chondrocytes, we used single-cell RNA-seq to profile chondrocytes in wild-type and nkx3.2 mutant juvenile heads. In zebrafish, fli1a:GFP and sox10:dsRed transgenes co-localize predominantly in chondrocytes (Askary et al., 2017; Giovannone et al., 2019). We therefore performed fluorescence activated cell sorting (FACS) of fli1a:GFP+/sox10:DsRed+ cells (Fig. S2) from pooled 21 dpf zebrafish heads (n=5 per genotype), followed by single-cell barcoded cDNA synthesis using the 10x Genomics platform and Illumina sequencing. After quality control, we obtained 1641 cells from wild-types and 1699 cells from nkx3.2el802 mutants (Fig. 5A). We used Seurat (Butler et al., 2018) on aggregated data and Uniform Manifold Approximation and Projection (UMAP) for dimension reduction and visualization (McInnes et al., 2018) to identify 11 distinct cell clusters (Fig. 5B). Five clusters are defined by high levels of cartilage collagen-encoding genes, col2a1a and col9a1a (Col2hi/Col9hi), with two of these representing chondrocytes undergoing S-phase DNA replication (pcna+) or mitosis (ube2c+/pcna+) (Fig. 5C,E). We also observed clusters with hypertrophic chondrocyte features (col10a1a+/spp1+/col2a1alo), articular chondrocyte features (f13a1b+/prg4b+/col2a1alo), perichondrial features (foxp4+/col2a1alo), two with mesenchymal/fibroblast connective tissue features (ifitm1+/col5a1+/col2a1alo) and one with an osteoblast signature (col1a1a+/ifitm5+) (Askary et al., 2016; Moffatt et al., 2008; Zhao et al., 2015).
All of the cell clusters present in wild types are also present in nkx3.2 mutants. However, we observed two clusters comprised almost entirely of mutant cells (Fig. 5D). The first represents mitotic chondrocytes (94% mutant cells), consistent with the abnormal chondrocyte proliferation seen in place of the missing jaw joint in mutants (Fig. 3E-K). A second chondrocyte cluster (98% mutant cells) was characterized by higher levels of sesn1 and eif4ebp3 (Fig. 5C, E), which encode inhibitors and repressed targets, respectively, of mTORC1-driven cellular growth and translation (Budanov and Karin, 2008; Tsukumo et al., 2016; Yogev et al., 2013). Both genes are induced by cellular stress in other systems (Budanov and Karin, 2008; Sukarieh et al., 2009; Taba et al., 2000). Compared with other chondrocytes, this mutant chondrocyte cluster also displayed higher levels of the prostaglandin D2 synthase gene ptgdsb.1, another stress-response pathway (Fig. 5E). However, many other regulators of chondrocyte biology, such as sox9a and hypertrophic genes runx2b and col10a1a, were unchanged in mutants.
Upregulation of ptgdsb.1 and sesn1 in nkx3.2 mutant chondrocytes
We next sought to understand the in vivo location of the mutant chondrocyte cluster identified in the single-cell analysis. In wild types at 21 dpf, we observed nkx3.2 expression in chondrocytes on either side of the jaw joint but not in the superficial zone (Fig. 6A,G,I). Expression of nkx3.2 was also seen in pre-hypertrophic chondrocytes of the bidirectional ceratohyal growth plates, in a largely non-overlapping pattern to expression of the hypertrophic marker col10a1a (Fig. 6C,K,M), similar to what has been observed in chick and mouse (Church et al., 2005; Provot et al., 2006; Tribioli et al., 1997). In nkx3.2el802 mutants, nkx3.2 mRNA transcripts were present in chondrocytes spanning the fused jaw joint (Fig. 6B,H,J), consistent with single-cell RNA-seq analysis showing nkx3.2 expression in the mutant-enriched chondrocyte cluster (Fig. 5E), yet col10a1a was not expressed (Fig. 6D). In wild types, ptgdsb.1 was expressed in the late hypertrophic zone in the upper jaw palatoquadrate and ceratohyal cartilages, in the proliferative zone of the ceratohyal growth plate, at modest levels in the lower jaw Meckel's cartilage and at lower levels in jaw joint articular cartilage, particularly on the upper jaw side of the joint (Fig. 6E,G,K,M). Although sesn1 was expressed on either side of the wild-type jaw joint, it was expressed at much lower levels at the articular surface, with only rare sesn1+/nkx3.2+ double-positive cells detected; in the ceratohyal growth plate sesn1 was expressed at high levels in all regions except the proliferative zone (Fig. 6I,O). In nkx3.2 mutants, ptgdsb.1 and sesn1 were upregulated across most chondrocytes within the fused jaw joint region, which we identified by co-expression of nkx3.2 (Fig. 6E-J). In the mutant ceratohyal growth plates, we observed largely normal zones of ptgdsb.1, col10a1a, sesn1 and nkx3.2 expression, with a potential modest upregulation of sesn1 in proliferative zones (Fig. 6K-P). Quantification of fluorescence intensity in the nkx3.2-positive jaw joint region confirmed upregulation of ptgdsb.1 and sesn1 in mutants (Fig. 6Q; ptgdsb.1: n=6, P=0.003; sesn1: n=4, P=0.03). These results show upregulation of stress response and mTORC1 pathway genes in the nkx3.2 mutant chondrocytes that fail to form a jaw joint and contribute to cartilage overgrowth (Fig. 6R).
Our results support temporally distinct roles for Nkx3.2 in embryonic specification of the zebrafish jaw joint and post-embryonic regulation of joint-associated cartilage growth. By comparing embryonic-only loss of nkx3.2 (morpholino knockdown) to constitutive loss (mutants), we uncovered post-embryonic requirements for Nkx3.2 in restricting subarticular cartilage growth and patterning the spine. In contrast to Nkx3.2−/− mice that die shortly after birth and do not display skeletal overgrowth (Akazawa et al., 2000; Tribioli and Lufkin, 1999), the cartilage overgrowth seen in the jaw region of adult viable nkx3.2 mutant zebrafish is a closer model to the supernumerary bones and lengthening of the wrist and digits in SMMD patients. As the number of cartilage segments in the digits appears to be proportional to cartilage length, with joint-derived signals such as Wnt9a inhibiting the formation of another joint within a certain distance (Hartmann and Tabin, 2001), the increased numbers of digit bones could be explained by cartilage overgrowth analogous to that seen in the zebrafish nkx3.2 mutant jaw. Spinal curvature defects are also more pronounced in zebrafish nkx3.2 than mouse Nkx3.2 mutants, again more closely aligning to the scoliosis seen in SMMD individuals (Hellemans et al., 2009; Tribioli and Lufkin, 1999).
Previous studies have focused on a potential role for Nkx3.2 in inhibiting hypertrophic maturation of chondrocytes. As in chick (Provot et al., 2006), we observed that nkx3.2 expression largely anti-correlates with col10a1a expression in hypertrophic chondrocytes. However, we observed no changes in hypertrophic cartilage maturation, or expansion of hypertrophic gene expression, in zebrafish nkx3.2 mutants, similar to what has been reported for mouse Nkx3.2 mutants (Akazawa et al., 2000; Tribioli and Lufkin, 1999). Instead, we observed a marked increase in chondrocyte proliferation within the mutant jaw joint region, and a more modest trend toward increased proliferation within the mutant growth plate of the ceratohyal endochondral bone. We did not, however, observe apparent cartilage overgrowth in the cranium and spine that are also affected in zebrafish nkx3.2 mutants, suggesting that Nkx3.2 has roles beyond restricting chondrocyte proliferation.
Our expression analysis also points to important differences in zones of chondrocytes between growth plates and joints. In the bidirectional growth plates of the ceratohyal (Fig. 6S), we observed a central proliferative zone of nkx3.2lo, sesn1lo, ptgdsb.1+ chondrocytes (Fig. 3E,I; Giovannone et al., 2019). Flanking these proliferative chondrocytes, we observed nkx3.2+, ptgdsb.1− zones, which we consider pre-hypertrophic based on the absence of col10a1a expression. Toward the middle of the ceratohyal we observed successive zones of col10a1a+ hypertrophic chondrocytes and then ptgdsb.1+; col10a1alo late hypertrophic chondrocytes. In contrast to the growth plate, we did not observe a prominent nkx3.2lo, ptgdsb.1+ proliferative zone close to the jaw joint. Relative to its expression in the pre-hypertrophic zones of the growth plates, nkx3.2 expression is also higher in subarticular chondrocytes of the juvenile jaw joint, consistent with a greater role in restricting chondrocyte proliferation. In nkx3.2 mutants, upregulation of ptgdsb.1 but not the hypertrophic marker col10a1a in the fused jaw joint region suggests that chondrocytes are adopting a partial ptgdsb.1+ proliferative identity rather than a late hypertrophic identity. This is supported by single-cell transcriptomic analysis showing that mitotic chondrocytes are greatly increased in abundance in nkx3.2 mutants, and our observations of increased proliferation and cartilage overgrowth at the fused jaw joint region. However, upregulation of sesn1 in the fused jaw joint region suggests that mutant chondrocytes are not equivalent to growth plate proliferative chondrocytes that normally lack sesn1. This is borne out by mutant chondrocytes forming a distinct cluster in our single-cell analysis. Our data therefore support a role for Nkx3.2 in limiting chondrocyte proliferation, particularly at the jaw joint where it is most abundantly expressed.
The higher levels of nkx3.2 expression at the wild-type jaw joint compared with the pre-hypertrophic zone of the growth plates may reflect slower growth of the subarticular zone relative to the expansive growth of the endochondral cartilage template; this is reflected by a fivefold lower proliferation rate at the jaw joint (Fig. 3K). These findings align with recent single-cell transcriptome analysis in mouse showing that chondrocyte zones in joints and growth plates are related but not identical (Lui et al., 2015). It will be interesting to test how Bmp, Fgf, and Shh signaling pathways regulate nkx3.2 expression at joints versus growth plates, given previous studies showing important roles for each of these pathways in regulating Nkx3.2 expression (Wilson and Tucker, 2004; Zeng et al., 2002). It will also be important to delete nkx3.2 only at post-embryonic stages to further confirm its requirements independent from initial specification of the fish jaw joint. Whereas we detected a modest increase in chondrocyte proliferation at the ceratohyal growth plate of nkx3.2 mutant zebrafish, we did not detect prominent cartilage overgrowth. This cannot be attributed to compensation by nkx3.1 (which we failed to detect in our single-cell RNA-seq analysis) or a broader genetic compensation (as we saw no evidence of the nonsense-mediated decay of nkx3.2 transcripts required for transcriptional adaptation). One possibility is that the preferential requirement for nkx3.2 at the zebrafish jaw joint reflects the lower levels of baseline proliferation there, and hence the greater magnitude of cartilage overgrowth when nkx3.2 is missing.
The increased proliferation of joint-associated chondrocytes in nkx3.2 mutants was linked to increased expression of stress-induced genes such as prostaglandin synthase and mTORC1 inhibitors. The presence of stress-induced pathways in late hypertrophic chondrocytes could reflect the hypoxic stress linked with their eventual apoptosis. However, strong expression of prostaglandin D2 synthase ptgdsb.1 in the proliferative zone of the wild-type growth plate was surprising, as it is not clear what types of stress these cells would experience. Expression of Ptgdsb proteins in the adult zebrafish skeleton has also been observed by proteomic analysis (Kessels et al., 2014). In mouse, single-cell transcriptomic analysis of the growth plate revealed expression of the prostaglandin E2 synthase Ptges3 in both the early proliferative zone and the late hypertrophic zone (Li et al., 2016), similar to ptgdsb.1 in zebrafish. The functions of prostaglandins in regulating chondrocyte biology, whether they are linked to proliferation, and why zebrafish would have D2 and mouse E2 subtypes remain outstanding questions. In the zebrafish testes, treatment with the D2 analogue BW-245C from 15 to 40 days of life resulted in upregulated expression of sox9a (Pradhan and Olsson, 2014), a master regulator of chondrogenesis in the skeletal system (Bi et al., 1999; Yan et al., 2005), though whether similar regulation occurs in chondrocytes has not been examined. In nkx3.2 mutants, it will be interesting to determine whether the upregulation of stress-induced pathways such as prostaglandin production, as well as modulators of mTORC1 signaling such as sesn1, are linked to the observed ectopic cartilage outgrowth, and possibly scoliosis and other skeletal defects. In the postnatal mouse, the proliferative zone of the growth plate undergoes an mTORC1-regulated transition to become a source of clonal stem cells that fuel cartilage growth (Mizuhashi et al., 2018; Newton et al., 2019). One possibility is that Nkx3.2 may function to normally regulate such an mTORC1-regulated transition below the joint surface, and to a certain extent in the pre-hypertrophic zone, to prevent cartilage overgrowth. The supernumerary bones and lengthening of the wrists and digits in SMMD patients may reflect an analogous role for Nkx3.2 in fine-tuning the rate of stem cell-mediated expansion of cartilage. In the future, it will be informative to determine how prostaglandins and other stress-induced pathways interact with mTORC1 signaling to regulate chondrocyte proliferation, and how Nkx3.2-mediated regulation of these pathways differentially fine-tunes stem cell-mediated chondrocyte expansion in joints versus growth plates.
MATERIALS AND METHODS
The Institutional Animal Care and Use Committees of the University of Southern California, Columbia University, and the University of Alberta approved all use of zebrafish in this study. Previously reported zebrafish lines used in this study include: Tg(sox10:dsRED)el10 (Das and Crump, 2012), trps1j127aGt (RRID:ZFIN_ZBD-GENO-100809-11) (Talbot et al., 2010), Tg(fli1a:EGFP)y1 (ZIRC Cat# ZL1085 ZDB-ALT-011017-8) (Lawson and Weinstein, 2002). nkx3.2 mutant lines were generated by TALEN (el802) and CRISPR (ua5011) mutagenesis. TALEN constructs were generated using TALEN targets L: TAATGCTCCGCGACCCGGTC and R: GACATGGCTGTGCGCAGTAA. For genotyping, fin clips were performed at 14 or 90 dpf and genotyped by PCR using GoTaq DNA polymerase (Promega). Genotyping for allele el802 was performed using nkx32-F: TAACCCTAATGCTCCGCGAC and nkx32-R: TGGCGAGACTCCTCTTTTCG primers to generate 118 bp WT and 104 bp MUT bands. nkx3.2 mutants inflate their swim bladders and are recovered at expected rates for analyses at 7-21 dpf (39/161: 24.2% mutant). When raised with wild types, mutant adults are recovered at lower rates (6/90: 6.6% mutant) perhaps owing to food competition with phenotypically normal clutch-mates. To raise sufficient numbers of mutants to adulthood for experimental purposes, genotyping was performed at 14 dpf. Mutants and wild types were then raised at similar densities in different tanks and size matched for downstream analyses at the indicated ages with no survival defects noted.
Wild-type (AB; sox10:GFP) and nkx3.2ua5011/ua5011 zebrafish were scanned using MILabs μCT at the School of Dentistry, University of Alberta, Canada. For µCT scanning, 30 and 60 dpf individuals from two distinct F2 populations for each genotype (wild type: n=10; nkx3.2ua5011: n=11 at each age) were fixed in 4% paraformaldehyde (PFA) for 24 h then dehydrated in a graded ethanol series. Scanning parameters were as follows: voxel size=10 µm; voltage=50 kV; current=0.24 mA; exposure time=75 ms. Skeletal reconstruction was performed using AMIRA (Milabs) by selecting and delineating regions of high tissue densities at a voxel size of 10 μm.
Histology and proliferation analysis
Anaesthetized fish were transferred into system water containing 4.5 mg/ml BrdU (B5002, Sigma-Aldrich) and immersed for 1 h before euthanasia. Tissue was fixed overnight in 4% PFA, transferred to 20% EDTA at room temperature for decalcification for 3 days and then processed for paraffin embedding. Immunohistochemistry was performed on 5 µm sections with antigen retrieval by steaming for 20 min in citrate buffer (pH 6.0). Primary antibodies used were rat anti-BrdU (1:100, MCA2060GA, Bio-Rad), mouse anti-PCNA (1:1000, P8825, Sigma-Aldrich) diluted in serum-free antibody diluent (Dako, Agilent) overnight at 4°C. Primary antibodies were detected by incubating slides in secondary AlexaFluor antibodies (A21094, A11001, Invitrogen) diluted at 1:500 in antibody diluent for 1 h at room temperature in the dark with Hoechst 33342 nuclear stain. H&E staining was performed on 5 µm sections. Briefly, sections were de-paraffinized, stained for 2 min in Hematoxylin solution (Harris Modified, HHS16, Sigma-Aldrich), followed by a brief acid rinse and water wash and incubation for 2 min in Blueing Reagent (Scott's Tap Water Substitute, 11160, EK Industries). Slides were then stained in Eosin Y (E6003, Sigma-Aldrich) for 30 s followed by mounting using Cytoseal 60 (8310, Richard-Allan Scientific) for imaging.
In situ hybridization and RNAscope
In situ hybridization labeled probes for col10a1a (Askary et al., 2016) and ptgdsb.1 were generated by PCR. ptgdsb.1 primers: Fwd, CTGCAAACATGACGAGTGTC; Rev, AATGCTTTGCCATTCTTTCCC. Digoxigenin (DIG)-labeled antisense probes were synthesized using SP6 or T7 polymerase (Roche). In situ hybridization was performed as previously described (Askary et al., 2016). Briefly, after deparaffinization, 5 µm sections were digested in 7.5 µg/ml proteinase K for 5 min and post-fixed in 4% PFA/0.2% glutaraldehyde for 20 min. Then, 1 µg of DIG-labeled probe per slide was incubated overnight at 62°C. Following hybridization, slides were washed 3× in 50% Formamide, 1× saline sodium citrate (SSC), 0.1% Tween-20 Wash buffer and 3× in MABT [0.1M maleic acid, 0.15M NaCl, 0.1% Tween-20 (pH 7.5)]. Slides were blocked in 2% Roche Blocking buffer followed by 1 h incubation in anti-DIG-AP (Roche). Color reaction was developed with NBT/BCIP (Roche). Nuclei were counterstained with Fast Red and slides were mounted using Cytoseal for imaging. RNAscope staining was performed using the RNAscope multiplex fluorescent assay v2 using the manufacturer's protocol for formalin-fixed paraffin-embedded 5 µm sections (323100, Advanced Cell Diagnostics). RNAscope probes for col10a1a (C1), nkx3.2 (C2), ptgdsb.1 (C3) and sesn1 (C3) were synthesized by Advanced Cell Diagnostics. For RNAscope analysis, nuclei were stained with DAPI. For all in situ hybridizations and RNAscope stains we saw similar patterns of expression in at least three individuals for each experiment.
As per Miller et al. (2003), knockdown of nkx3.2 (bapx1) expression during early embryonic development was performed by injecting antisense morpholino oligos. bapx1-MO1 (5′-GCGCACAGCCATGTCGAGCAGCACT-3′; ATG start complementary sequence underlined) was purchased from GeneTools, diluted to 3 mg/ml and injected into the yolk of one-cell-stage Tg(sox10:DsRed)el10 embryos. Morpholino-injected animals to be raised to adulthood were screened for jaw joint fusion by imaging sox10:DsRed+ craniofacial cartilages at 7 dpf.
Adult specimens were fixed overnight at 4°C in 4% PFA. Following a 1 h wash in Tris/MgCl2, cartilage was stained overnight in 0.01% Alcian Blue in 10 mM MgCl2. Samples were re-hydrated through an ethanol/100 mM Tris (pH 7.5) series into 0.5% KOH and incubated overnight at 4°C. Samples were bleached in 3% H2O2/0.5% KOH for 6-8 h and then neutralized overnight in 35% NaBO4. Tissue was cleared using 1% Trypsin in 35% NaBO4 for 3-4 h followed by staining in 0.02% Alizarin Red (pH 7.5) overnight. Samples were washed in 50% glycerol/0.5% KOH to remove residual stain and stored in 100% glycerol before imaging using a Leica S8 APO stereomicroscope.
Single-cell preparation and sequencing
For single-cell profiling of juvenile chondrocytes, zebrafish double-positive for fli1a:GFP; sox10:DsRed were selected and screened for mutant joint fusion phenotype. PCR analysis to confirm genotypes was performed at 14 dpf. At 21 dpf, wild-type and nkx3.2el802 mutant craniofacial skeletons were micro-dissected and processed into a single-cell suspension using mechanical and chemical dissociation for 30 min at 28°C in protease solution [0.25% Trypsin, 1 mM EDTA (pH 8.0), PBS] containing collagenase-D (Roche Life Sciences) (n=5 animals/genotype). GFP+/DsRed+ live cells were sorted into suspension medium (1% calf serum, 0.8 mM CaCl2, 50 U/ml penicillin, 0.05 mg/ml streptomycin) and immediately loaded onto a 10x Genomics microfluidic chip. Barcoded cDNA libraries were generated using 10x Genomics Chromium Single Cell 3′ Library and Gel Bead Kit v.2 according to the manufacturer's instructions. Library sequencing on Illumina NextSeq was performed to a depth of at least 1,000,000 reads/cell for each library. Cell Ranger software (v. 3.1.0, 10x Genomics) was used for barcode recovery, genome alignment (Ensembl GRCz11) and to generate gene-by-cell count matrices with default parameters for each library. RNA-seq files have been deposited in the NCBI Gene Expression Omnibus and are available under accession number GSE151354.
R software with Seurat Version 3 was used for downstream analysis. Cells expressing <200 and >2500 unique RNA Features and >5% mitochondrial RNA were excluded. Datasets were aggregated using the merge function. Aggregate data was log normalized, centered and scaled using default settings of NormalizeData and ScaleData functions. Linear dimensional reduction was performed with principal component analysis using the RunPCA function. Cell clusters were determined using Find Neighbors and FindClusters functions. Cluster marker genes were identified using FindMarkers function. In order to focus the analysis on the craniofacial chondrocytes affected in nkx3.2 mutants, cells of the gill were identified based on the expression of ucmaa and ppp1r1c (P.F. and J.G.C., unpublished) and fin cells based on expression of posterior hox genes (hox9-13) (Fig. S3) (Ahn and Ho, 2008). We then excluded gill and fin cells from downstream analysis using the subset function and re-clustered the craniofacial-specific cells as above. UMAP non-linear dimensional reduction was used to visualize clusters and gene expression.
Imaging and statistical analysis
Fluorescence imaging of live zebrafish and sections was performed using a Zeiss LSM800 confocal microscope and Zen software or Leica SP8 confocal and Leica LAS software. Skeletal preparations were imaged using Leica S8APO and DM2500 microscopes. Image processing was performed with the same settings used among all images from each experiment using Fiji (Schindelin et al., 2012) and Adobe Photoshop. Proliferating PCNA+ cells and total Hoechst+ nuclei at the jaw joint and ceratohyal growth plates were counted using the cell counter function of Fiji and expressed as a %PCNA+/total cells. The regions selected for proliferation analysis were identified by wild-type joint morphology and/or by performing RNAscope staining for nkx3.2 in adjacent 5 µm sections (e.g. Fig. 3G,H PCNA staining versus Fig. 6I,J nkx3.2 expression). Fluorescence intensity of sesn1 and ptgdsb.1 transcripts was measured in at least two sections per animal in Fiji using the ‘Measure’ feature of a region of interest defined by co-expression with nkx3.2 and shown as mean intensity/area. Statistical analysis was performed in GraphPad Prism software using unpaired two-tailed t-test function (proliferation analysis) or unpaired two-tailed t-test with Welch's correction (fluorescence intensity analysis) and represented in the scatter dot plot as individual data points with mean±s.d.
We thank Megan Matsutani and Jennifer DeKoeyer Crump for fish care; Jeffrey Boyd at the University of Southern California Stem Cell Flow Cytometry Core for FACS; Nellie Nelson for single-cell RNA-seq library preparation; the University of Southern California's Center for High-Performance Computing (HPC) for single-cell RNA-seq data alignments; Amjad Askary for generating the nkx3.2el802 allele; and Ted Allison for his helpful feedback and facilitating our collaboration.
Conceptualization: J.S., J.G.C.; Methodology: J.S., N.N., T.M., P.B.; Software: J.S., P.F.; Validation: J.S., N.N., A.N.K., P.F.; Formal analysis: J.S., N.N., A.N.K., T.M., P.B., D.G., J.G.C.; Investigation: J.S., N.N., A.N.K., T.M., P.B., P.F.; Resources: P.F., D.G., J.G.C.; Data curation: J.S., P.F.; Writing - original draft: J.S., J.G.C.; Writing - review & editing: J.S., P.F., J.G.C.; Visualization: J.S., P.B.; Supervision: D.G., J.G.C.; Project administration: J.G.C.; Funding acquisition: J.S., D.G., J.G.C.
This work was supported by National Institutes of Health grants R00DE027218 (to J.S.) and R35DE027550 (to J.G.C), and Natural Sciences and Engineering Research Council of Canada grant RGPIN-2014-06311 (to D.G.). Deposited in PMC for release after 12 months.
RNA-seq files have been deposited in GEO under accession number GSE151354.
Peer review history
The peer review history is available online at https://dev.biologists.org/lookup/doi/10.1242/dev.193409.reviewer-comments.pdf
The authors declare no competing or financial interests.