Spatiotemporal Characterization of Human Early Intervertebral Disc Formation at Single-Cell Resolution

Abstract

The intervertebral disc (IVD) acts as a fibrocartilaginous joint to anchor adjacent vertebrae. Although several studies have demonstrated the cellular heterogeneity of adult mature IVDs, a single-cell transcriptomic atlas mapping early IVD formation is still lacking. Here, the authors generate a spatiotemporal and single cell-based transcriptomic atlas of human IVD formation at the embryonic stage and a comparative mouse transcript landscape. They identify two novel human notochord (NC)/nucleus pulposus (NP) clusters, SRY-box transcription factor 10 (SOX10)⁺ and cathepsin K (CTSK)⁺ , that are distributed in the early and late stages of IVD formation and they are validated by lineage tracing experiments in mice. Matrisome NC/NP clusters, T-box transcription factor T (TBXT)⁺ and CTSK⁺ , are responsible for the extracellular matrix homeostasis. The IVD atlas suggests that a subcluster of the vertebral chondrocyte subcluster might give rise to an inner annulus fibrosus of chondrogenic origin, while the fibroblastic outer annulus fibrosus preferentially expresseds transgelin and fibromodulin . Through analyzing intercellular crosstalk, the authors further find that notochordal secreted phosphoprotein 1 (SPP1) is a novel cue in the IVD microenvironment, and it is associated with IVD development and degeneration. In conclusion, the single-cell transcriptomic atlas will be leveraged to develop preventative and regenerative strategies for IVD degeneration.

Keywords: annulus fibrosus; intervertebral disc formation; notochord; nucleus pulposus; single-cell RNA sequencing.

Figure 1

Single‐cell profiling of the embryonic axial skeleton cell atlas during early IVD formation. A) Schematic image of the experimental workflow. B) Representative images of human embryonic axial skeleton sections at different developmental stages. Upper lane, representative images of H&E staining; lower lane, representative images of KRT8 IHC staining. Scale bar, 500 µm. C) UMAP visualization of human embryonic axial skeleton cells identified nine different clusters. Each dot corresponds to one single cell colored according to cell cluster. D) Heatmap showing the scaled expression of differentially expressed genes for each human cluster. E) Proportion of developmental stages and cell clusters (upper lane) and fraction of human cell clusters at the Noto, Trans, and NPL stages (lower lane). F) UMAP visualization of mouse embryonic axial skeleton cells identified 10 different clusters. Each dot corresponds to one single cell colored according to cell cluster. G) Fraction of mouse cell clusters at the Noto, Trans, and NPL stages. H) Heatmap showing pairwise Pearson correlations in the global transcriptome between human and mouse cell clusters.

Figure 2

Characterization of developing NC/NP cells during early IVD formation. A) UMAP visualization (left panel) of the seven subclusters of 2336 NC/NP cells defined during human early IVD formation and cell counts (right panel) of each human NC/NP cell subcluster at the Noto, Trans, and NPL stages. B) Heatmap revealing the scaled expression of differentially expressed genes for each human (left) and mouse (right) NC/NP cell subcluster. C) UMAP visualization (left panel) of the seven subclusters of 589 NC/NP cells defined during mouse early IVD formation and cell counts (right panel) of each mouse NC/NP cell subcluster at the Noto, Trans, and NPL stages. D) Representative IHC images of SOX10 in human axial skeleton sections at the indicated developmental stages. Scale bar, 200 µm. E) Representative IF images of Sox10 in mouse axial skeleton sections at indicated developmental stages. Asterisk indicates DRG. Scale bar, 200 µm. F) Representative images of lineage tracing in Sox10‐Cre;mT/mG postnatal mouse IVDs. Scale bar, 200 µm. Asterisk indicates DRG. G) Representative images of lineage tracing in Sox10‐CreERT2;mT/mG mouse IVDs. Schematic images of tamoxifen administration are indicated on the top of each image. Asterisk indicates DRG. Scale bar, 200 µm. H) Representation analysis of GO categories showing different functions for the SOX10⁺ and Sox10⁺ clusters. I) UMAP plot color‐coded for pathway activities scored by GSVA per cell within hNC subpopulations. Notochord morphogenesis (GO:0048570): the process in which the anatomical structures of the notochord are generated and organized. Notochord development (GO:0030903): the process whose specific outcome is the progression of the notochord over time, from its formation to the mature structure.

Figure 3

Characterization of matrisome clusters during early IVD formation. A) Representative IHC images of ACAN (upper lane) and COL1A1 (lower lane) in developing human IVDs at the indicated developmental stages. Arrows indicate the increased ECM in the human NP region. Scale bar, 100 µm. B) FPKM expression of the indicated signature genes in developing human NC/NP tissues captured by LCM during early IVD formation. N of 9, 10, 11, and 13 weeks = 4; N of 12 weeks = 3. C) Heatmap showing pairwise Pearson correlations of expressed matrisome genes in human developing NC/NP cells. D) The number of expressed genes associated with six matrisome patterns in the indicated human (upper) and mouse (lower) NC/NP cell subclusters. E) Representation analysis of GO categories showing different functions for TBXT⁺/T⁺ (left), CTSK⁺/Ctsk⁺ (middle), and KRT15⁺/Krt15⁺ clusters (right). Asterisks indicate the GO categories only enriched in human NC/NP subclusters. F) Violin plots showing the expression levels of representative genes associated with six matrisome patterns in the indicated human (left) and mouse (right) NC/NP cell subclusters. G) Representative IHC images of CTSK in human axial skeleton sections at the indicated developmental stages. The red arrow indicates the ECM in the NP region at the NPL stage. Scale bar, 200 µm. H) Representative images of lineage tracing in the Ctsk‐Cre;mT/mG mouse IVD. Scale bar, 200 µm. I) Dot plot showing the mean expression of selected mechanosensitive ion channel genes among seven human NC/NP cell subpopulations. The dot size indicates the percentage of cells in subclusters with detected expression.

Figure 4

Characterization of vertebral chondrocytes during early IVD formation. A) UMAP visualization (upper) of the four vertebral chondrocyte subclusters defined during human early IVD formation at the Noto, Trans, and NPL stages (lower). B) UMAP plots and violin plots showing the expressions of MATN4, COL2A1, ACAN, MATN1, CYTL1, COL1A1, SFRP2, HHIP, VEGFA, GDF5, PIEZO2, and ADGRG6 on UMAP. C) Representative IHC images of MATN1 (upper) and POSTN (lower) in human axial skeleton sections at the indicated developmental stages. Asterisks indicate AF. Scale bar, 200 µm. D) Heatmap revealing binary regulon activities analyzed with SCENIC in each vertebral chondrocyte subpopulation. “ON” indicates active regulons, and “OFF” indicates inactive regulons. E) Representative fluorescent images of the Ccn2‐GFP mouse IVD section at P0. OAF, outer annulus fibrosus; IAF, inner annulus fibrosus; VB, vertebral body; EP, end plate; NP, nucleus pulposus. Scale bar, 200 µm. F) Violin plot showing the expression of stem/progenitor markers, including THY1, PDPN, CD164, PDGFRB, MCAM, and CTSK. G) Monocle pseudotime trajectory axis revealing the progression of human vertebral chondrocytes. ECO, endochondral ossification; EP/IAF, end plate/inner annulus fibrosus. H) Expression of selected differentially expressed genes for each branch. The cells were ordered in pseudotime analyzed with Monocle and colored in developmental stages. The full line represents the ECO, while the dotted line represents the EP/IAF. ECO, endochondral ossification; EP/IAF, end plate/inner annulus fibrosus.

Figure 5

Characterization of OAF cells during early IVD formation. A) UMAP visualization of the OAF cluster during human early IVD formation at the Noto, Trans, and NPL stages. B) Representative violin plots showing the expressions of PAX1, PAX9, MKX, COL5A1, SOX9, PDGFRL, TAGLN, FMOD, and BGN during human early IVD formation. C) Representative IHC images of TAGLN (upper) and FMOD (lower) in human axial skeleton sections at the indicated developmental stages. Scale bar, 100 µm. D) UMAP visualization of the OAF cluster during mouse early IVD formation at the Noto, Trans, and NPL stages. E) Representative violin plots showing the expressions of Pax1, Mkx, Col5a1, Sox9, Pdgfrl, and Fmod during mouse early IVD formation. F) Representative IF images of Fmod (left) and Tagln (right) in mouse axial skeleton sections at E14.5. The dotted line indicates the presumptive OAF. The asterisks indicate the presumptive IAF. Triangles indicate NP. VB, vertebral body. Scale bar, 200 µm. G) Dot plot showing the mean expression of selected mechanosensitive ion channel genes in the human OAF along with early IVD formation. Dot size indicates the percentage of cells with detected expression. H) Representative images of lineage tracing in Piezo2‐Cre;Ai9 mouse IVDs. Scale bar, 200 µm.

Figure 6

Overview of the crosstalk networks between developing NC/NP and vertebral chondrocytes during early IVD formation. A) Overview of the cellular network between developing NC/NP cells and vertebral chondrocytes during human early IVD formation. Dots indicate cell subpopulations. The dot size indicates the relative quantity of each subpopulation. The thickness of the directed line indicates the relative quantity of significant ligand–receptor pairs between any two pairs of subpopulations. B) Dot plot showing the communication probability of the indicated ligand–receptor pairs between four vertebral chondrocyte subclusters (sending signals) and seven NC/NP subclusters (accepting signals). C) Dot plot showing the communication probability of the indicated ligand–receptor pairs between seven NC/NP subclusters (sending signals) and four vertebral chondrocyte subclusters (accepting signals). D) UMAP plots showing the expression of SPP1 and CD44 on the UMAP. E) Representative IF images of SPP1 in human axial skeleton sections at the indicated developmental stages. Orange dotted squares are magnified under each image. The white dotted line indicates the presumptive OAF. The asterisks indicate the presumptive IAF. AF, annulus fibrosus; OAF, outer annulus fibrosus; IAF, inner annulus fibrosus; VB, vertebral body; NC, notochord; NP, nucleus pulposus; OC, ossification center. Scale bar, 100 µm. F) Representative IF images of CD44 in human axial skeleton sections at the indicated developmental stages. Scale bar, 500 µm. G) Representative IF images of EMCN in human axial skeleton sections at the indicated developmental stages. The arrows indicate the immunofluorescent signals surrounding the VB. The asterisks indicate the cartilage canals. AF, annulus fibrosus; VB, vertebral body; NP, nucleus pulposus; OC, ossification center. Scale bar, 500 µm. H) Western blot results showing the expression of Vegfa, Cd9, Hif1a and Sox9 in ATDC5 chondrogenic progenitor cells after rSpp1 treatment. I) Representative IHC images of Spp1 in control and degenerated rat IVDs. Scale bar, 200 µm.

Figure 7

Spatial transcriptome analysis of developing IVD cells. A) UMAP visualization of the six connective tissue subclusters during human early IVD formation at the Noto, Trans, and NPL stages. B) Dot plot revealing the expressions of signature genes among the six human connective subclusters. C) Representation analysis of GO categories showing different functions for the hCT1 (left), hCT2 (middle), and hCT3 (right) clusters. D) UMAP visualization of the seven connective tissue subclusters during mouse early IVD formation at the Noto, Trans, and NPL stages. E) Dot plot revealing the expressions of signature genes among the seven mouse connective subclusters.