HFD-feeding results in diabetic peripheral neuropathy

C57BL/6J mice were fed a high-fat diet (HFD) over a 12 weeks period (from 4 to 16 weeks old). DXA analysis of body composition demonstrated increased body weight, fat body mass and percentage of fat mass over animals fed a normal diet (ND (Fig. 1a–d, 1.5-, 4- and 2.2-fold increase, respectively).21,22 Next, serologic studies were performed to confirm the induction of glucose metabolic dysfunction (Fig. 1e–g).22 Fasting blood glucose was significantly elevated in HFD-fed animals (Fig. 1e), and the HFD-fed cohort demonstrated glucose intolerance (Fig. 1f) and insulin resistance (Fig. 1g). The presence of a small peripheral neuropathy was first assessed by quantitative analysis of intraepidermal nerve fiber (IENF) density and behavioral paw withdrawal testing23 (Fig. 1h–j). Quantification of the immunostaining for the pan-neural marker protein gene product 9.5 (PGP 9.5) revealed a 32% reduction in IENF density in the HFD-fed cohort (Fig. 1i), concomitant with a significantly delayed paw withdrawal time during a hotplate test24 (Fig. 1j). These results confirmed the development of a type 2 diabetic polyneuropathy in HFD-fed mice.25

Fig. 1

High-fat diet (HFD) feeding induces obesity, glucose intolerance, insulin resistance, and peripheral neuropathy. Normal diet (ND) or HFD feeding in C57BL/6J mice was instituted on week 4 of life, with analysis up to week 16 of life. a Representative Dual-energy X-ray absorptiometry (DXA) images of body fat. b Body weight after 0, 4, 8 and 12 weeks dietary treatment. n = 8 mice per group. c Fat mass and d Percentage of fat mass after 12 weeks dietary treatment. n = 8. e Fasting blood glucose after 12 weeks feeding. n = 8. f Blood glucose concentrations 0–120 min after glucose injection (1 g/kg body weight ip) in fasted mice following 12 weeks dietary treatment. n = 8. g Blood glucose concentrations 0–120 min after insulin administration (0.5 unit/kg ip) in fasted mice following 12 weeks dietary treatment. n = 8. h, i Images and quantification of Protein Gene Product 9.5 (PGP 9.5) immunofluorescent staining (green) in sagittal sections of mice hind paw skin after 12 weeks dietary treatment. n = 5. Scale bar: 100 µm. j Paw thermal withdrawal latency measurements at 0, 4, 8 and 12 weeks of dietary treatment. n = 8. Graphs represent average values ± 1 SD, *P < 0.05, **P < 0.01, and ***P < 0.001. Comparisons between groups were analyzed by unpaired Student’s t test. In subfigures (b and j), comparisons between groups were analyzed versus D7

Diabetic peripheral neuropathy affects long bone periosteum

Diabetic polyneuropathy is well-known to affect cutaneous innervation,26 but to our knowledge there are no reports on skeletal innervation in the context of type 2 diabetes. As the most innervated skeletal location,27 the long bone periosteum was evaluated. Nerve density within the periosteum was quantified within three representative bones, including the femur, tibia and the 1st metatarsus. Immunostaining for the pan-neuronal marker β III-tubulin (TUBB3) was used to quantify total innervation. Calcitonin gene-related peptide (CGRP) and tyrosine hydroxylase (TH) were used to quantify sensory and sympathetic nerve fibers, respectively. Quantification of TUBB3+ nerve fibers showed a 76% decrease at the femoral midshaft under HFD-fed conditions (Fig. 2a). CGRP+ (Fig. 2b) and TH+ (Fig. 2c) nerve fibers also decreased by 62% and 55%, respectively. The distal femur periosteum was also assessed for total innervation, as well as sensory and sympathetic innervation. The results revealed similar significant decreases in TUBB3+ (Fig. 2d), CGRP+ (Fig. 2e) and TH+ (Fig. 2f) nerve fibers by 69%, 38% and 51%, respectively. The innervation at the tibial bone was similarly affected by HFD-feeding. Quantitative evaluation of periosteal innervation showed a 67% reduction of TUBB3+ nerve fibers at the midshaft and a 58% reduction at the distal tibia (Fig. S1). The innervation at the 1st metatarsus was also affected in HFD-fed mice. Quantitative evaluation of periosteal innervation showed 43% and 41% at the midshaft and distal metatarsal periosteum, respectively (Fig. S2). These findings indicated that HFD-induced diabetic polyneuropathy extended to skeletal innervation affecting both more proximal long bones and distal small tubular bones of the feet.

Fig. 2

High-fat diet (HFD) feeding induces long bone neuropathy of the periosteum. Normal diet (ND) or HFD feeding in C57BL/6J mice was instituted on week 4 of life, with analysis at week 16 of life. a The innervation of the femoral midshaft periosteum visualized through the pan-neural marker β III-tubulin (TUBB3) immunostaining, represented in green. b Sensory innervation in the same region highlighted using Calcitonin gene-related peptide (CGRP) immunostaining, shown in yellow. c Sympathetic innervation observed through Tyrosine hydroxylase (TH) immunostaining, shown in red. n = 5. (Scale bar: 100 µm). d The distal femoral periosteum exhibits similar patterns for TUBB3+ total innervation, e CGRP+ sensory innervation and f TH+ sympathetic innervation. n = 5. Scale bar: 100 µm. Graphs represent average values ± 1 SD, *P < 0.05 and **P < 0.01. Comparisons between groups were analyzed by unpaired Student’s t test

Bone alterations are associated with diabetic peripheral neuropathy

Next, quantitative analysis of cortical and trabecular bone in the femur, tibia and 1st metatarsus were analyzed with microCT (Fig. 3 and Figs. S3–S5). Quantitative analysis of the microstructural parameters of the mid-shaft cortical bone including cortical bone area (Ct.Ar), cortical bone perimeter (Ct.Pm), cross-sectional thickness (Cs.Th) and polar moment of inertia (pMOI) were significantly decreased among HFD-fed mice showing a 11.9%, 4%, 8.3%, and 16.3% reduction in femur, respectively (Fig. 3a, b). Similarly, significant reductions in these parameters were observed in both the tibia and the 1st metatarsus (Figs. S4a, b, S5a, b). For trabecular bone in the femur, fractional bone volume (BV/TV), trabecular thickness (Tb.Th), and trabecular number (Tb.N) were reduced by 42.6%, 11.4%, and 35.1%, respectively (Fig. 3c, d). These reductions were also evident in the tibia and the 1st metatarsus (Figs. S4c, d and S5c, d). Additionally, trabecular spacing (Tb.Sp) increased by 7.7% in the femur (Fig. 3d), with similar increases observed in the tibia and the 1st metatarsus (Figs. S4d and S5d) within the HFD-fed cohort.

Fig. 3

High-fat diet (HFD) feeding results in cortical and trabecular bone alterations. Normal diet (ND) or HFD feeding in C57BL/6J mice was instituted in week 4 of life, with analysis at week 16 of life. a µCT images of femoral midshaft cortical bone. b µCT quantifications of cortical area (Ct.Ar), cortical perimeter (Ct.Pm), cross-sectional thickness (Cs.Th) and polar moment of inertia (pMOI). n = 8 mice per group. c µCT images of distal femoral trabecular bone. d µCT quantifications of Bone volume per total volume (BV/TV), Trabecular thickness (Tb.Th), Trabecular number (Tb.N) and trabecular separation (Tb. Sp). n = 8 mice per group. Graphs represent average values ± 1 SD, *P < 0.05, **P < 0.01 and ***P < 0.001. Comparisons between groups were analyzed by unpaired Student’s t test

These quantitative analyses were further utilized to compare changes in cortical and trabecular bone parameters among the femur, tibia, and 1st metatarsal bone. The comparative analysis revealed that reductions in cortical bone parameters were more pronounced in the more distally located small tubular bones of the feet (Table S3). In contrast, changes in trabecular bone parameters did not exhibit substantial differences between bones (Table S4).

Multi-tissue single cell RNA-Sequencing infers somatosensory neuron-to-periosteal signaling

Communication between sensory nerve and bone cells has been demonstrated at a cellular level, forming a neuro-osteogenic network within the periosteum.28,29 We next set out to decipher potential intercellular communication between DRG sensory neurons and skeletal resident cells within the periosteum (Fig. 4a). Long bone periosteum (combined femoral and tibial periosteum) was isolated and examined by scRNA-Seq analysis in ND conditions. Six cellular clusters were identified, including mesenchymal cells (n = 434), endothelial cells (n = 127), pericytes (n = 32), neutrophils (n = 2 287), macrophages (n = 521), and T cells (n = 461) (Fig. 4b). A previously established scRNA-Seq dataset from mouse lumbar DRG neurons was reanalyzed, with 16 identified neuronal subclusters30 (Fig. 4c). These included 7 subclusters of CGRP neurons, low threshold mechanosensory neurons (LTMRs), nonpeptidergic nociceptors, proprioceptors, somatostatin+ (SST) neurons, and transient receptor potential M8 (TRPM8) cold sensitive neurons. Next, interaction analyses between DRG neurons and periosteal cells was performed using NicheNet (Fig. 4d–f),31 for which, cells from periosteum were defined as receiver cells and categorized into two distinct niches: the periosteal stromal/vascular niche (mesenchymal cells, endothelial cells, and pericytes) and the periosteal immune cell niche (neutrophils, macrophages, and T cells). With the neuronal derived ligands as senders and receiver cells defined, a prediction for ligand-receptor interactions was performed by combining gene expression with existing knowledge of signaling pathways and gene regulatory networks. The goal was to identify predicted neuronal ligands secreted from DRG neurons and their downstream signaling genes involved specifically in periosteum cellular niche regulation.31

Fig. 4

Cellular communication by interactome analysis between DRG neurons and long bone periosteum. a Schematic illustrating the interaction between dorsal root ganglia (DRG) sensory neurons and periosteal cells within their microenvironment. b UMAP projection of periosteum cell clusters via single-cell RNA sequencing (scRNA-Seq), accompanied by violin plots of known gene markers for each cluster. n = 4 042 total cells. Cell number per cluster shown in parentheses. c UMAP visualization of mouse lumbar dorsal root ganglia (DRG) neurons by scRNA-Seq, re-analyzed from ref. 102n = 1 757 total DRG neurons. Cell number per cluster shown in parentheses. d Interaction plot depicting the cell-cell communication between DRG neurons and periosteal stromal/vascular cells (including mesenchymal cells, endothelial cells, and pericytes). DRG neuron-derived ligands and receptors expressed on stromal/vascular cells are shown under normal diet conditions. n = 593 stromal/vascular cells. e Expression profile of ligands in each DRG neuron cell cluster, highlighting specificity to stromal/vascular cells. f Analysis of the regulatory potential of DRG neuron-derived ligands and their downstream target genes in stromal/vascular cells. g Expression profile of ligands in each DRG neuron cell cluster and corresponding receptors within immune cells of the periosteum (including macrophage, neutrophils and T cells). n = 3 269 immune cells. h, i Analysis of the regulatory potential of DRG neuron-derived ligands and their downstream target genes in immune cells

The ligand-receptor interaction analysis between DRG neurons and the periosteum stromal/vascular niche revealed enriched expression of numerous neural ligands, including growth factors (Vegfa, Nrg1, Gmfb, Pdgfa and Pdgfb) and neuropeptides (Gal, Calca, Bdnf, and Calcb), predicted to pair with stromal/vascular cell receptors32,33 (Fig. 4d). Among the growth factors, Vascular Endothelial Growth Factor A (VEGFA) was highly expressed by Aβ−Field−LTMR neurons, a critical regulator of osteoblast function, bone homeostasis and skeletal repair.15,34,35,36,37,38 The membrane glycoprotein Neuregulin 1 (NRG1), was highly expressed across multiple DRG neuron types, including proprioceptors, Aβ−Field−LTMR and Aβ–RA−LTMR neurons. NRG1 has been implicated during osteoblastogenesis39 and in cartilage and bone regeneration in zebrafish.40 Another growth and differentiation factor enriched in the interaction analysis was Glia Maturation Factor Beta (GMFB). Gmfb was highly expressed in nonpeptidergic nociceptors, CGRP-ζ and CGRP-η neurons. A recent study showed that Gmfb knockout mice are protected against osteoporosis by suppressing osteoclast hyperactivity in diabetic conditions.34 Platelet-derived growth factors A and B (PDGFA and PDGFB) were also identified as key soluble growth factors in the interaction analysis. Pdgfa exhibited high expression in CGRP α, ε, and Trpm8 neurons, while Pdgfb was most highly expressed in CGRP-η neurons. PDGFA has been implicated in nerve-regulated digit tip regeneration,41 and our group implicated its potential role as a nerve-derived factor in regulation of heterotopic ossification (HO).42 Galanin (GAL) was highly expressed by CGRP-γ, α and β neurons. Several studies showed that Galanin (GAL) was present in bone marrow mesenchymal cells, endothelial cells and nerve fibers within the periosteum suggesting a role in bone development and repair. Calca and Calcb genes encoding Calcitonin gene-related peptide (CGRP), a crucial neuropeptide for bone metabolism regulation, were also enriched in the interaction analysis, primarily in CGRP clusters. CGRP was shown to play a pleiotropic effect on bone cells, promoting osteoblast differentiation, inhibiting osteoclast activity,43 and enhancing reparative bone formation in response to mechanical loading and fracture.44 Brain-Derived Neurotrophic Factor (BDNF) was also enriched in the interaction analysis. Bdnf was highly expressed by CGRP-α, Sst and CGRP-ε neurons. Several studies have shown the positive effects of BDNF on bone formation by enhancing osteoblast differentiation, new bone formation and maturation.45,46 Moreover, BDNF and its receptor TrkB are present at various stages of the bone formation process and are upregulated in human osteoblasts, suggesting a role in bone development47,48 (Fig. 4e). Although the interaction analysis did not predict previously described bone-related neural ligands such as Fibroblast growth factor 1 and 9 (FGF1, FGF9),49,50,51 Follistatin Like 1 (FSTL1)17 and Sonic hedgehog signaling molecule (SHH),52 their expression was present within the present dataset of DRG neurons.

Following up the ligand-receptor interaction, various downstream signaling events were predicted to regulate target genes involved in cell cycle regulation (Cdkn1a and Dusp1), cell growth and proliferation (Ackr3, Ergr1, Fosb and Socs3), cell differentiation (Fos, Fosb and Ier2) and transcriptional regulation (Jun, Junb, Jund, Nr4a1 and Cited2). Interestingly, several downstream target genes were described to play an important role in regulating osteoblast differentiation and bone formation among them, Fos and Jun family members,53 the zinc finger protein Zfp36,54 the heme oxygenase Hmox1,55 the cyclin-dependent kinase inhibitor Cdkn1a,53 the suppressor of cytokine signaling Socs3,56 and the nuclear receptor Nr4a153 (Fig. 4f).

Next, a similar interaction analysis was performed between DRG neurons and immune periosteal niche (including neutrophils, macrophages, and T cells). The ligand-receptor interaction analysis revealed enriched expression of numerous neural ligands including, neuropeptides (Nrtn) and neurotransmitters (Sst, Pnoc, Nmu, and Gal), cell adhesion molecules (l1cam), cytokines and chemokines (Ccl2 and Cx3cl1) predicted to pair with immune cells receptors (Fig. 4g). Neurturin (Nrtn) gene was highly expressed by Aδ-LTMR neurons. It has been suggested that NRTN may play several context-dependent roles in modulating immune cell function and inflammation.57 Somatostatin (SST), a neurotransmitter known for the neuroendocrine inhibitory effects across multiple systems58 was mainly expressed by the Sst neurons. Sst is known to exert inhibitory effects on immune cell functions, cytokine production, proliferation, and inflammatory responses.59,60 Another neurotransmitter, Prepronociceptin (Pnoc) gene was specifically expressed by TrpM8 neurons. PNOC is the precursor protein for nociceptin, a neuropeptide involved in regulating inflammation, pain, and arousal.61 Neuromedin U gene (Nmu) was specifically expressed by C-LTMR neurons. Nmu is an immunoregulator described to be a potent activator of immune cells, particularly in the context of type 2 inflammation, allergic responses, and autoimmune conditions.62,63 Galanin gene (Gal) was highly expressed by CGRP-γ, CGRP-α, and CGRP-β neurons. Gal was described to modulate the expression of cytokines suggesting a role in regulating inflammatory responses.64,65 The L1 Cell Adhesion Molecule gene (l1cam) was highly expressed in CGRP-θ and Sst neurons. L1CAM was shown to contribute to the migration and extravasation of immune cells.66 Chemokines and cytokines genes such as C-C Motif Chemokine Ligand 2 (Ccl2) and C-X3-C Motif Chemokine Ligand 1 (Cx3cl1) were highly expressed in CGRP-θ and CGRP-α, respectively. CCL2 was shown to be a key mediator of the crosstalk between immune and bone cells, influencing bone homeostasis67,68 (Fig. 4h). Interestingly, some ligands, such as SST, DUSP18, CRH, NTS, GAL and GMFB, were found to activate both mesenchymal and immune cell populations, suggesting shared signaling pathways in the periosteal niche. The ligand-receptor analysis was followed by a downstream signaling inquiry highlighting the activated genes within the periosteal immune niche involved in cell cycle regulation (Cks2 and Dusp1), cell growth and proliferation (Erg1, Fosb, Jun family and Pim1), cell differentiation (Klf6 and Zfp36), and immune response (Nfkbi family and Tnf) (Fig. 4i).

The analysis model presented here describes secreted neural ligands regulating downstream genes under homeostatic conditions. Interestingly, distinct sets of ligands were predicted to activate downstream genes in periosteal mesenchymal cells versus periosteal immune cells, predicted to regulate cell cycle and proliferation, as well as cell function.

HFD-induced bone neuropathy disturbs long bone periosteal cellular signalization

Next, scRNA-Seq analysis focused on signaling changes in the transcriptome among periosteal cells exposed to either ND or HFD-feeding conditions (Fig. 5 and Fig. S6). First, changes in the cellular composition of the six periosteal cell clusters (Fig. S6a), were analyzed, revealing subtle differences in the distribution of clusters between ND and HFD conditions (Fig. S6b, c). Next, overall changes in cellular proliferation were assessed using a proliferation module score ratio (HFD/ND) calculated across periosteum cell clusters (Fig. S6d). A significant reduction in proliferative score was seen under HFD conditions among several cell types including pericytes, T cells and mesenchymal cells (Fig. S6d). Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis using all periosteal cells highlighted increased terms expression under ND conditions including ossification and angiogenesis, as well as signaling pathways such as mitogen-activated protein kinase (MAPK), Vascular endothelial growth factor (VEGF), and Transforming growth factor beta (TGFβ). Conversely, representative GO terms enriched within the HFD periosteal cells included bone resorption, hypoxia, negative regulation of angiogenesis, and cellular immune response (Fig. S6e). Further transcriptional analysis focused on the mesenchymal cell cluster, dissecting their cell subcluster frequency, phenotypes, and molecular signaling pathway modifications associated to HFD feeding conditions (Fig. 5). Mesenchymal cells were analyzed and distributed across three subclusters defined by characteristic gene markers, subcluster 1: mesenchymal progenitor cells (Pdgfrα+Ly6a+, n = 813), subcluster 2: Pre-osteoblasts (Ly6a+Runx2+Lepr+, n = 320), and subcluster 3: Osteoblasts (Runx2+Alpl+Bglap+, n = 177) (Fig. 5a). Next, characteristic gene markers were examined across the three subclusters as shown by heatmap of differentially expressed gene (DEG) profiles (Fig. 5b). Complete DEGs are provided in (Fig. S7) along with characterize periosteal cell markers such as Gli1,69Postn,70Ctsk,71Acta272 and Nes.73

Fig. 5

High-fat diet (HFD) feeding disturbs the skeletal cell phenotype and neural-skeletal cell communication within the periosteal microenvironment. 12 weeks after the initiation of the dietary treatment, four left femurs and tibias from ND and HFD mice were dissected. Periosteal cells were isolated, and scRNA-sequencing and analysis were performed. a UMAP of mesenchymal cell subclusters including mesenchymal progenitors, pre-osteoblasts, and osteoblasts, along with known gene markers for each by violin plot, and the cell number ratio among ND and HFD treated groups. Black indicates ND cells. Gray indicates HFD cells. b Heatmap for differentially expressed genes (DRGs) in each mesenchymal cell subcluster. c UMAP showing the pseudotime trajectory of the mesenchymal cell subclusters along with progenitor gene markers (Pdgfrα and Ly6a) as well as osteoblast markers (Bglap and Alpl). Black line represents the trajectory graph. Plots were generated using Monocle 3.0.1.2. d–g’ Linear graph analysis of phenotype changes across pseudotime including gene modules for stemness, proliferation, adipogenesis and osteoblastogenesis among ND (blue line) and HFD (black line) fed mesenchymal cells, and module index scoring of mesenchymal cell subclusters (mesenchymal progenitors, pre-osteoblasts and osteoblasts). Dashed lines indicate delineates early, mid and late pseudotime. h–k’ Linear graph analysis of dysregulated signaling pathways, including MAPK, TGFβ, Wnt and mTor signaling, across pseudotime among ND (blue line) and HFD (black line) fed groups and module index scoring of mesenchymal cell subclusters. Dashed gray lines in module score graphs represent early, mid, and late pseudotime. Graphs represent average values ± 1 SD. Module scoring data was analyzed using the Kolmogorov-Smirnov test. *P < 0.05, **P < 0.01, ***P < 0.001 and ****P < 0.00 01 in comparison to ND control. 1 255 total mesenchymal cells analyzed

Next, pseudotime trajectory analysis revealed a continuous cell differentiation path,74 starting from a mesenchymal progenitor cell expressing Pdgfrα+ and Ly6a+ genes. As cells progressed along the pseudotime axis, they transitioned through an intermediate pre-osteoblast state, characterized by the expression of Ly6a, Runx2, and Lepr genes. The terminal branch represented fully differentiated osteoblast cells, expressing Runx2, Alpl, and Bglap (Fig. 5c).75 Here, the trajectory analysis laid the foundation for further investigation of the phenotypical changes and regulatory signaling under ND and HFD-fed conditions. Using the R package Escape,76 we evaluated th