Diverse mesenchymal heterogeneity captured by single-cell sequencing. To understand the cellular composition of the embryonic coronal suture, we performed single-cell RNA sequencing at E15.5 and E17.5 on dissected coronal sutures, including small amounts of frontal and parietal bone, after removing the skin and brain (Fig. 1a). We filtered using Seurat 3 R-Package23 and obtained 8279 cells at E15.5 (median of 2460 genes per cell) and 8682 cells at E17.5 (median of 3200 genes per cell) (Fig. 1b). We identified 14 cell clusters through unsupervised graph clustering of the two datasets combined (Supplementary Table 1). Osteogenic and mesenchymal cell types were identified based on the expression of broad mesenchyme/fibroblast (Col1a1) and osteoblast (Sp7) markers (Fig. 1b, dotted line; Supplementary Fig. 1a-b). The identities of clusters outside the osteogenic/mesenchymal subset were resolved using previously reported markers, and included chondrocytes, myeloid cells, mast cells, lymphocytes, pericytes, osteoclasts, endothelial cells, neurons, and glia (Fig. 1b, c). All the major cell types in our analysis were present at E15.5 and E17.5 (Supplementary Fig. 2a). Chondrocytes were especially abundant at E15.5 (Supplementary Fig. 2b), highlighting the close proximity of the E15.5 coronal suture to the chondrocranium. Myeloid cells were more abundant at E17.5, consistent with reports of increased myeloid differentiation during late embryonic stages24 (Supplementary Fig. 2b).
To analyse the osteogenic and mesenchymal cell types that might comprise and support the coronal suture, we re-clustered the osteogenic/mesenchymal population and obtained 14 clusters present at both E15.5 and E17.5 (Fig. 1d-f). Analysis of enriched genes for each cluster allowed us to assign probable identities to each cell type (Fig. 1g, Supplementary Table 2), which we validated by in situ hybridization as described below. We observed one ectocranial cluster strongly over-represented at E15.5 and a meningeal cluster over-represented at E17.5. Osteogenic cells were also more abundant in the E17.5 dataset, although it is unclear whether this reflects true biological differences versus differing cell capture between the dissections (Fig. 1a, f).
Diversity of meningeal layers below the coronal suture. The meninges are involved in the development of the calvaria and underlying brain that they separate12. They comprise dura mater, arachnoid mater, and pia mater. To determine the identity of meningeal tissues included in our dissections, we utilized a recent transcriptomic study of murine E14 meninges as a guide25. Markers associated with the pia mater (Ngfr, Lama1, Rdh10) were not co-enriched in any of our clusters, consistent with the pia mater being removed with the brain during dissections (Fig. 2a). In contrast, markers enriched within the arachnoid mater (Aldh1a2, Cldn11, and Tbx18) were abundant in MG4, and dura mater markers (Gja1, Fxyd5, and Crabp2) in MG3 and MG4, and to a lesser extent MG1 and MG2 (Fig. 2a).
To resolve the identity of MG4, we performed in situ experiments for a highly specific MG4 marker, Gjb6, in combination with the arachnoid/dura mater marker, Crabp2 (Fig. 2b). Crabp2 and Gjb6 overlapped below the bone with Crabp2 single-positive cells (MG3) found above the Crabp2+/Gjb6+ domain (Fig. 2c). Immunofluorescence for Crabp2 and Gja1 at E17.5 confirmed that these arachnoid/dura mater markers are excluded from the pia mater (Supplementary Fig. 3a). Rgs5 + pericytes were interspersed with Gjb6 + cells in the MG4 arachnoid layer, consistent with the prominent vasculature extending from the border of the dura mater and through the arachnoid mater to the pia mater26 (Supplementary Fig. 3c). In the Crapb2+/Gjb6− layer (MG3), we also observed co-expression of Crabp2 with Nppc, with a zone of Nppc+/Crabp2− cells (MG2) above (Fig. 2d). In addition, Nppc was co-expressed with the known dura mater marker Fxyd525 (Supplementary Fig. 3b, d). MG3 may represent the dural border cell layer that has been described using EM in adult meninges27,28.
The Nppc+/Crabp1− population above Nppc+/Crabp2+ cells is consistent with MG2 identity in our UMAP analysis (Fig. 2b, d). This domain was also positive for Matn4 and Ctgf, markers that only overlap with Nppc in MG2 (Fig. 2f, Supplementary Fig. 3b, e). Matn4 expression does not overlap with the MG3/MG4 marker Crabp2 (Fig. 2e), and chondrocytes in other sections express high levels of Matn4 but not Ctgf (Supplementary Fig. 3f). We also detected some Matn4 cells that were low or absent for Nppc, and these were situated closer to the suture and calvarial bones than the Nppc-high cells (Fig. 2f). Our UMAP analysis indicates that these Matn4+/Nppc− cells represent MG1 (Fig. 2b). MG1 expressed higher levels of the chondrocyte markers Col2a1 and Acan than MG2, but not at the level seen in bona fide chondrocytes (Supplementary Fig. 3g). These findings are consistent with MG1/MG2 representing a periosteal dura mater population29 primed to form cartilage, with the MG1 population having a stronger chondrogenic-like program. These data reveal a diversity of mesenchyme cell types within the meninges, with those closest to the suture sharing features with chondrogenic cells and likely functioning as a specialized connective tissue bridging the calvarial bones to the meninges (Fig. 2g).
Diversity of ectocranial layers above the coronal suture. Ectocranial mesenchyme guides the migration of early osteogenic cells and helps establish suture boundaries7,21,22. Consistent with analysis of the mouse frontal suture30, we identified clusters enriched for markers of the dermis (Ly6a, Dpt; EC1-3) (Fig. 1g, Fig. 3a). Ly6a is expressed in hypodermis of neonatal mouse skin31 and we observe co-expression of Ly6a with additional EC1-3 markers, Pi16 and Dpt, within a single layer of mesenchyme above the skull bones (Fig. 3b, c; Supplementary Fig. 4a, b). Our previous work had identified expression of Jag1 within both an ectocranial layer and suture mesenchyme cells22, and scRNAseq analysis shows high expression of Jag1 in EC1, and to a lesser extent EC2 and EC3 (Fig. 1f-g, Supplementary Fig. 4g). Co-expression of Jag1 with Pi16 confirms Jag1 expression within the hypodermal layer (Supplementary Fig. 4h).
In addition to hypodermal EC1-3 layers, we also noted a prominent EC4 population defined by C1qtnf3 and Tnmd expression in our UMAP analysis (Fig. 3a). In situ experiments revealed two C1qtnf3+ ectocranial layers on either side of the Pi16+ hypodermis (Fig. 3d, Supplementary Fig. 4c). However, only the C1qtnf3+ layer between the calvarial bones and hypodermis expresses Tnmd, suggesting that this corresponds to EC4, which we term “suprasutural mesenchyme” (Fig. 3a, e; Supplementary Fig. 4d). It is likely that the outer C1qtnf3+ layer was removed along with the skin during dissection and thus not included in our single-cell analysis. Consistently, Epha4, which has also been shown to display ectocranial expression7, was most strongly expressed in this outer C1qtnf3+ layer despite not being particularly enriched in clusters EC1-4 in UMAP analysis (Supplementary Fig. 4i). Thus, Jag1 and Epha4 appear to label distinct ectocranial layers.
Interestingly, just underneath the C1qtnf3+/Tnmd+ layer (EC4), and sitting on top of the suture, we observed a C1qtnf3−/Tnmd+ population (LIG) that was enriched for a number of genes associated with tendon/ligament development (e.g. Scx, Tnmd, Mkx, Thbs2, Bgn), as well as markers of smooth muscle (e.g. Acta2, Tagln, Myl9) (Supplementary Fig. 4j). Two of the most specific markers for LIG were Tac1, a gene that encodes four neuropeptides including Substance P that is implicated in tendon mechanosensation32, and Chodl, a membrane protein with expression in tendons of humans33 and mouse limbs34 (Fig. 3a). In situ hybridization for Tac1 and Chodl revealed highly restricted expression in mesenchyme overlying the coronal suture and connecting the edges of the frontal and parietal bones (Fig. 3f, h, Supplementary Fig. 4e). Although broader, in situ experiments showed that Scx and Tnmd were also expressed within this layer (Fig. 3g; Supplementary Fig. 4f, k), and we observed only minimal overlap between C1qtnf3 and Tac1 expression (Supplementary Fig. 4l). Thus, the LIG cluster represents a transcriptionally distinct and spatially confined subset of the ectocranium overlying the coronal suture and connecting the adjacent bones.
Distinct osteoblast trajectories in the suture versus periosteum. The suture is a source of new osteoblasts for skull bone expansion. To better understand the trajectories of osteogenesis within and around the sutures, we identified six osteogenic clusters (OG1, OG2, OG3, OG4, PO1, PO2) based on expression of the osteoblast transcription factors Runx2 and Sp735,36,37, as well as Dlx5 and Dlx638,39,40, (Supplementary Fig. 5). Pseudotime analysis of these clusters using Monocle3 revealed OG1 to be the earliest lineage cells, followed by two proliferative branches (PO1, PO2) and then distinct OG2 and OG3 trajectories that converged on the osteoblast cluster OG4 (Fig. 4a, b). These pseudotime trajectories were reflected by increased expression of Runx2 and Sp7 from OG1 to OG4 and accumulation of the mature osteoblast markers Ifitm5 and Dmp1 in OG4 (Fig. 4c). Prominent cluster markers included Erg and Pthlh (OG1), Lef1 and Inhba (OG2), Mmp13 and Podnl1 (OG3), and Ifitm5, Dmp1, and Sost (OG4) (Fig. 4e). Pseudotime visualization shows Erg (OG1) appearing earliest and shutting down as Lef1 (OG2) and Mmp13 (OG3) become expressed, followed by their extinguishment and progressive appearance of Ifitm5, Dmp1, and Sost (OG4) (Fig. 4d).
In situ validation at E15.5 showed co-expression of OG1 markers Erg and Pthlh in suture mesenchyme and extending along the edges of the frontal and parietal bone tips (Fig. 4f, Supplementary Fig. 6a). Interestingly, Erg but not Pthlh was asymmetrically distributed along the bones, with stronger expression above the frontal and below the parietal bone. Gli1 and Prrx1 showed similar asymmetric expression above the frontal bone, and Gli1 and Six2 below the parietal bone (Supplementary Fig. 7a-c). The different asymmetric mesenchyme patterns were not attributable to disproportionate populations of mesenchyme above and below the bones, as ectocranial and meningeal layers were generally evenly distributed along the bones (Figs. 2, 3). These findings highlight asymmetric distribution of the earliest osteogenic cells around the bone fronts, which may play a role in ensuring the later reproducible overlap of the parietal over the frontal bone.
We next examined expression of markers for OG2 and OG3, as pseudotime analysis suggested that these represent alternative pathways to osteoblasts (Fig. 4a). Markers for both OG2 (Lef1, Inhba) and OG3 (Podnl1, Mmp13) displayed expression along the outer domains of the Sp7+ bone surface, consistent with a pre-osteoblast identity. Whereas Lef1+/Inhba+ cells were most abundant at the edges of the growing bones near the suture, Podnl1+/Mmp13+ cells were enriched on bone surfaces further away from the suture (Fig. 4g, h; Supplementary Fig. 6b, c). OG2 may therefore represent more specialized suture-resident pre-osteoblasts, and OG3 periosteal pre-osteoblasts more generally found on bone surfaces.
Clusters PO1 and PO2 are enriched for markers of proliferation (Mki67, Cenpf, Top2a; Supplementary Fig. 6e), as well as expression of genes from suture-resident clusters, for example Erg (OG1) and Lef1 and Inhba (OG2) (Fig. 1a, Fig. 4e). In situ validation revealed Cenpf expression around the tips of the growing bones and extending into the sutures, though in comparison to Erg it appeared largely absent from the central part of the suture (Supplementary Fig. 6f). These data suggest that proliferative pre-osteoblasts are concentrated at the bone tips, consistent with previous direct evidence for proliferative osteogenic cells at the leading edges of the calvarial bones20.
OG4 represents mature osteoblasts and osteocytes, as revealed by markers such as Ifitm5, Dmp1 and Sost and its terminal position in the pseudotime analysis41 (Fig. 4a, c). We observe co-expression of Sp7 with Ifitm5 throughout most of the bone, except for the leading tips and the surfaces of bones which are Sp7-positive only, consistent with bone growth from both the leading edges and the periosteal surfaces (Fig. 4i, Supplementary Fig. 6d). Similarly, we observe protein expression of the osteoblast marker Cd200 in bone and the Wnt-responsive transcription factor Tcf7, a pre-osteoblast marker, along the bone tips and periosteal surfaces43,44,45, consistent with Cd200 expression in cluster OG4 and Tcf7 expression in both OG4 and the suture-specific pre-osteoblast cluster OG2 (Fig. 1g, Supplementary Fig. 8a-c). The mature osteocyte marker Sost, which labels the oldest cells in pseudotime, was expressed in only a few osteocytes distant from the suture (Fig. 4i). The spatial organization between osteoblasts, pre-osteoblasts, and progenitors was confirmed with double in situs between Podnl1 and Pthlh or Ifitm5 (Supplementary Fig. 6g-h). These findings are consistent with bone elongation occurring through osteogenesis at the suture and bone tips, with a distinct osteogenic pathway in the periosteum contributing to bone thickening (Fig. 4j).
Signaling interactions between mesenchymal layers captured from single-cell analysis. To gain insights into potential signaling between osteogenic cells and adjacent ectocranial and meningeal mesenchyme, we surveyed our datasets for expression of ligands and receptors (Fig. 5a-b). In OG4 osteoblasts, we observed enrichment of Tgfb1, Bmp3, Bmp4, Ihh, and Pdgfa, with Ihh46 and Pdgfa47,48 known to be secreted by osteoblasts. In addition to osteoprogenitors (OG1), Pthlh ligand was expressed in several meningeal clusters. In the ectocranial clusters EC4 and LIG, which are situated closest to the suture, we detected enrichment of Tgfb3, Fgf9, Fgf18, Wnt5a, Wnt9a, Wnt11, and Igf1 ligands. In the meningeal clusters closest to the suture (MG1, MG2), we also detected Tgfb3 and Igf1 ligand expression, as well as Tgfb2 that has previously been shown to be a critical meningeal-derived factor for suture morphogenesis49. Reciprocally, we detected expression of the Tgfβ receptor Tgfbr3, the Fgf receptors Fgfr1 and Fgfr2, the Wnt receptors Fzd1, Fzd2, and Fzd6, the Ihh receptor Ptch1, the Pthlh receptor Pth1r, and the Pdgfa receptor Pdgfra in various osteogenic clusters (Fig. 5b). Interestingly, Fgfr1 was expressed more strongly in suture-resident pre-osteoblasts (OG2) and Fgfr2 in periosteal pre-osteoblasts (OG3). To capture signalling interactions in an unbiased approach, we interrogated our data using the CellPhoneDB package50, scoring the ligand-receptor pairings between all clusters (Fig. 6c). Interactions between osteogenic clusters were notably weak, highlighting potential roles for osteogenic - non-osteogenic interactions in coronal suture regulation. Consistently, we noted that the mesenchymal populations closest to the suture (EC4 above, MG1/MG2 below) had the strongest interactions with osteogenic cells. The proliferative PO1 and pre-osteoblast OG3 clusters appeared to have the strongest connections with the surrounding mesenchyme, consistent with these populations residing at the surfaces of bones. Of the individual ligand-receptor interactions that were identified between EC4/MG1/M2 and osteogenic clusters, several pathways are relevant to coronal suture development, including Fgf, Tgfβ, and Wnt signalling51 (Supplementary Fig. 9).
Synostosis related genes display diverse expression patterns. To determine if genes underlying craniosynostosis in humans51,52 coalesce within a specific cell type, we mapped coronal and midline synostosis genes onto our dataset. Genes associated with midline suture synostosis were less abundant, suggesting that divergent gene expression patterns may contribute to suture-specific fusions (Fig. 6a). We noted that 8/17 coronal synostosis genes were selectively enriched in PO1 or PO2, suggesting misregulation of proliferating osteoblasts as a common mechanism of synostosis. By scoring every cell for the average expression of each coronal or midline synostosis gene, we noted that coronal synostosis scores were more restricted to the meningeal and osteogenic clusters, while midline synostosis genes were broadly distributed across clusters (Fig. 6b). Coronal synostosis genes with particular enrichment in PO1 and PO2 include Fgfr1, Fgfr2, Tcf12, Twist1, and Zeb2. Consistently, in situ validation revealed expression of Twist1 in suture mesenchyme and cells lining the bone fronts, with Twist1 expression also overlapping with the proliferative pre-osteoblast marker Cenpf (Fig. 6c, Supplementary Fig. 10). Selective enrichment of Twist1 and Tcf12 in proliferative pre-osteoblasts is consistent with functional studies showing roles for these transcription factors in regulating the balance between bone differentiation and proliferation9,20.