Abstract
Neonatal injury alters synaptic transmission in the spinal superficial dorsal horn (SDH), resulting in aberrant amplification of ascending nociceptive transmission. Astrocytes orchestrate synapse development and function across the CNS and play a critical role in the emergence and maintenance of persistent pain. However, little is currently known about the postnatal development of spinal astrocytes, nor about how the maturation of SDH astrocytes is impacted by early life injury. Here, we used a hindpaw incision model of postsurgical pain in postnatal day (P) 3 mice of both sexes to elucidate the effects of neonatal injury on the maturation of SDH astrocytes. Three-dimensional morphological analysis of individual astrocytes revealed that incision elicits age-dependent changes to astrocyte structure. At P4, spinal astrocytes in incised mice show increased size and complexity compared with naive controls. This is reversed at P10 and P24, as astrocytes from incised mice are smaller and less ramified compared with their naive counterparts. Transcriptomic analysis of spinal astrocytes revealed acute changes to gene expression after neonatal injury, as 76 differentially expressed genes (DEGs) were identified at P4 (such as Thbs1, Efemp1, Acta1, Acta2, Tpm2, and Fgf14), which included genes related to cell motility and cytoskeletal organization, but very few DEGs were detected at P10 and P24. Lastly, we identified that microglial engulfment of astrocyte material occurs in the developing dorsal horn and this process is altered by neonatal incision in a sex-dependent manner. These data illustrate, for the first time, that neonatal injury alters the postnatal development of spinal astrocytes.
Significance Statement
Neonatal tissue damage persistently remodels synaptic circuits in the spinal superficial dorsal horn (SDH), which has been implicated in the ability of early life injury to “prime” developing nociceptive pathways. While astrocytes clearly regulate synapse formation, pruning and function across the CNS, nothing is known about the degree to which neonatal injury modulates the properties of astrocytes within the developing SDH. The present study demonstrates that neonatal hindpaw incision evokes age-dependent transcriptional and morphological plasticity in spinal astrocytes, highlighted by a prolonged reduction in the size and complexity of astrocytes following early life injury. These findings yield new insight into the cellular mechanisms by which neonatal tissue damage can exert long-term effects on spinal nociceptive processing.
Introduction
Exposure to invasive procedures within a critical period of early life causes heightened pain responses during early childhood (Fitzgerald et al., 1989; Grunau et al., 1994; Taddio et al., 1997; Klein et al., 2009) and exacerbates hyperalgesia in response to noxious stimuli later in life (Hermann et al., 2006; Walker et al., 2009, 2018; Hohmeister et al., 2010; Vederhus et al., 2012; Walker, 2019). This may reflect the persistent sensitization of nociceptive pathways following neonatal injury, referred to as neonatal priming (Ren et al., 2004; LaPrairie and Murphy, 2009; Zhang et al., 2010; Beggs et al., 2012; Walker et al., 2016).
The superficial dorsal horn (SDH) of the spinal cord, which is a key modulatory site for ascending pain transmission (Melzack and Wall, 1965), undergoes significant functional and structural reorganization during the early postnatal period (Fitzgerald and Jennings, 1999; Baccei and Fitzgerald, 2004; Cordero-Erausquin et al., 2005; Bremner and Fitzgerald, 2008; Ingram et al., 2008). The activity of primary afferents shapes the maturation of the SDH (Jennings and Fitzgerald, 1998; Koch et al., 2012; Koch and Fitzgerald, 2013). Exposure to painful procedures during early life perturbs the balance between synaptic excitation and inhibition within the dorsal horn network (Li et al., 2013, 2015) and facilitates the potentiation of sensory synapses onto ascending projection neurons (Li and Baccei, 2016).
These changes to SDH synaptic circuits may be driven by astrocytes, which are key regulators of synaptic maturation and function (Chung et al., 2015) and react strongly to a range of pathological states (Ji et al., 2019; Diaz-Castro et al., 2021; Burda et al., 2022). Astrocyte-secreted factors regulate synaptogenesis (Baldwin and Eroglu, 2017; Irala et al., 2024) and synaptic strength by influencing presynaptic neurotransmitter release (Albrecht et al., 2012; Crawford et al., 2012) and postsynaptic receptor recruitment (Beattie et al., 2002; Stellwagen et al., 2005; Allen et al., 2012; Blanco-Suarez et al., 2018). Astrocytes also regulate neuronal excitability by buffering ions (Kofuji and Newman, 2004; Shih et al., 2013) and neurotransmitters (Danbolt, 2001; Verkhratsky and Nedergaard, 2018). The highly ramified and complex morphology of astrocytes is central to their myriad functions and allows them to make extensive multicellular contacts within the neuropil, contributing to their ability to monitor and interact with their environment (Baldwin et al., 2024). Astrocyte leaflets associate closely with synapses, seemingly filling in the spaces between the presynaptic and postsynaptic densities (Špaček, 1985; Jones and Greenough, 1996).
Astrocytes in the developing brain exhibit significant age-dependent changes in their morphological, transcriptional, and functional properties (Cahoy et al., 2008; Felix et al., 2021). Unfortunately, the regional heterogeneity in the properties of astrocytes (Ben Haim and Rowitch, 2017; Chai et al., 2017; Morel et al., 2017) poses significant challenges to extrapolating the knowledge gained from the study of developing brain astrocytes to the immature SDH. As a result, little is known about the postnatal maturation of astrocytes residing within spinal nociceptive circuits. Furthermore, the degree to which neonatal tissue damage modulates the development of SDH astrocytes has never been investigated. Collectively, these gaps in knowledge prevent a better understanding of the cellular and molecular mechanisms governing nociceptive processing in the developing spinal cord under normal and pathological conditions.
Therefore, the goal of this study was to elucidate the morphological and transcriptional plasticity of developing dorsal horn astrocytes in the absence or presence of neonatal tissue injury. Our imaging studies demonstrate that spinal astrocytes from incised mice are acutely larger and more complex than their naive counterparts, but this effect is reversed 1 week and 3 weeks after injury. RNA sequencing of dorsal horn astrocytes reveals an upregulation of genes related to cytoskeletal organization and cell migration shortly after incision, which resolves within 1 week. Overall, this study harnesses the power of modern imaging and bioinformatics to identify nuanced changes to the structure and transcriptomic signature of astrocytes after early life injury, which are predicted to influence nociceptive processing in the developing spinal cord.
Materials and Methods
Animals
All experiments were performed in accordance with animal welfare guidelines outlined by the Institutional Animal Care and Use Committee at the University of Cincinnati. To visualize astrocytes using fluorescence-based microscopy, we crossed Aldh1l1-Cre/ERT2 BAC transgenic mice which express a tamoxifen-inducible Cre recombinase under the control of the Aldh1l1 promoter (Aldh1l1-Cre/ERT2, The Jackson Laboratory; #031008; Srinivasan et al., 2016) with Rosa-CAG-LSL-tdTomato-WPRE mice (Ai9, The Jackson Laboratory, #007909) which express tdTomato (tdTOM) fluorescent protein from the Rosa26 locus under the CAG promoter in the presence of Cre recombinase. The resulting offspring (Aldh1l1-Cre/ERT2+/−; Rosa26-LSL-TdTOM+/−) will be referred to as Aldh1l1-tdTOM mice. For experiments analyzing the morphology of individual astrocytes which required sparse labeling, a single intragastric injection (Pitulescu et al., 2010) of a dilute dose of tamoxifen (0.25 mg/kg tamoxifen dissolved in corn oil vehicle containing 2.5% ethanol) was delivered at postnatal day (P) 1. For the investigation of astrocyte–microglial interactions, a single intragastric injection of full-dose tamoxifen (50 mg/kg) was delivered at P1 to induce Cre recombination in Aldh1l1-expressing cells.
For RNA sequencing experiments, Aldh1l1-Cre/ERT2 mice were crossed with mice expressing Cre-dependent CAG-Sun1/sfGFP (Jackson #021039) to produce offspring in which astrocyte nuclei were labeled with Sun1-GFP following tamoxifen induction (Aldh1l1-Cre/ERT2+/−; Rosa26-LSL-Sun1-GFP+/−); these will be referred to as Aldh1l1-GFP mice. Full-dose tamoxifen was administered via intragastric injection to Aldh1l1-GFP pups on P1 and P2.
Hindpaw incision
On P3, half of the pups in a litter of either Aldh1l1-tdTOM or Aldh1l1-GFP mice were randomly selected to undergo a hindpaw incision while the remaining littermates were exposed to isoflurane for equivalent duration as anesthesia-only naive controls. An incision was made into the plantar surface of the left hindpaw skin and underlying muscle as described previously (Brennan et al., 1996; Li et al., 2013). The skin was then closed with 7-0 suture (Ethicon). Mice were allowed to mature until the time of tissue harvest at P4, P10, or P24, corresponding to postsurgical time points of 1, 7, and 21 d. For investigation of microglia and astrocyte interactions, mice were harvested at P8, or 5 d post-incision.
Tissue preparation for three-dimensional confocal microscopy
Aldh1l1-tdTOM mice were euthanized with sodium pentobarbital (Fatal-Plus, Vortech Pharmaceuticals) and then transcardially perfused with 4% paraformaldehyde (PFA). The spinal cord was removed via ventral laminectomy and then fixed in PFA overnight at 4°C. Tissue was then transferred to 30% sucrose solution overnight at 4°C or until the tissue sank to the bottom. The right side of the spinal cord (i.e., contralateral to the incision) was marked with a needle puncture, and then thick transverse sections (100 µm) of the spinal lumbar enlargement were cut on a Leica cryostat and transferred to 0.01 M phosphate-buffered saline (PBS).
For experiments involving the visualization of microglial lysosomes, slices were incubated in blocking buffer (0.3% Triton X-100 and 10% donkey serum in PBS; Sigma Millipore) for 1 h, followed by primary antibodies against the microglial marker Iba1 (1:1,000 dilution, 500 µg/ml stock concentration, Novus Biologicals #NB100-1028) and the lysosome-associated membrane protein CD68 (1:500 dilution, 100 µg/ml stock concentration, Abcam #ab125212) overnight at 4°C. Sections were then washed in PBS and incubated in a species-specific secondary antibody (Alexa Fluor-488, Thermo Fisher Scientific #A-11055; and Alexa Fluor-647, Thermo Fisher Scientific #A-21244) at a 1:500 dilution overnight at 4°C to allow for adequate penetration into thick sections.
Sections were mounted on glass slides (Thermo Fisher Scientific #1255015), and then a thin layer of vacuum grease was applied around the sections to form a water-resistant barrier. One to two drops of EZ Clear [80% Nycodenz (Accurate Chemical & Scientific #100334-594), 7 M urea, 0.05% sodium azide prepared in 0.02 M sodium phosphate buffer; Hsu et al., 2022] were then applied to the sections for passive tissue clearing. A glass coverslip was then placed on top of the slide preparation. Great care was taken to remove any air bubbles during slide mounting to prevent urea crystallization. A cotton bud was then used to press the glass coverslip to the vacuum grease, forming an airtight seal around the tissue and EZ Clear. It is critical to compress the vacuum grease as much as possible to prevent working-distance issues during imaging. Imaging was conducted when sections were fully clear, ∼5–10 min after the tissue contacted EZ Clear. It should be noted that the original EZ Clear protocol (Hsu et al., 2022), which includes the lipid extraction step with immersion in 50% tetrahydrofuran (THF), was tested alongside the passive clearing method described above and did not yield significantly improved results. Given that the present study was conducted in the gray matter of the spinal cord, and between P4 and P24, the amount of lipid within the region of interest was not enough to warrant the lipid extraction step.
3D image acquisition
A Yokogawa SoRa W1 spinning disk confocal microscope in W1 mode equipped with a Hamamatsu Fusion BT SCMOS camera using a 50 mm pinhole with a 60× oil immersion objective (numerical aperture, 1.42; refractive index, 1.515; frame size, cropped to 1,024 × 1,024 pixels, 0.11 mm/pixel; exposure time, 50 ms) was used to acquire all images used for astrocyte morphological analysis. Z stacks capturing entire astrocytes were imaged using 0.2 µm z steps for accurate three-dimensional (3D) reconstruction. All images within a dataset were generated with the same microscope, with identical laser power and gain settings using NIS-Elements v 5.41 software. To avoid selection bias, all fluorescently labeled protoplasmic astrocytes within the medial third of the left dorsal horn (ipsilateral to the incision) were imaged (Swett and Woolf, 1985).
Images for the analysis of microglia and astrocyte interactions were acquired using the same imaging systems and parameters as described above. However, to better resolve these cell–cell interactions, 3D deconvolution was performed using the Richardson–Lucy deconvolution algorithm in NIS-Elements.
3D image analysis
Bitplane Imaris software (Oxford Instruments) was used for quantitative 3D morphological analysis of spinal astrocytes. Using Imaris Surfaces to generate a 3D model of the astrocyte, astrocyte surface area and cell volume were measured as previously described (Testen et al., 2020). The same creation parameters and background subtraction-based threshold values (surface grain size, 0.215 µm; diameter of largest sphere, 0.805 µm; threshold value empirically determined across a sample of 20 images, 7.4408) were used for all images for accurate between-group comparisons. The Imaris Filaments Tracer was then used to obtain metrics of astrocyte complexity, namely, the number of branch points, terminal points, and 3D Sholl intersections (Turk and SheikhBahaei, 2022). As with Imaris surface generation, the same creation parameters for filaments analysis were applied to each image (algorithm, Autopath with no loops; soma starting diameter, 5.37 µm; segment seed point diameter range, 0.322–3.22 µm; segment seed point threshold, 10.977). The Imaris Convex Hull XTension was used to create a 3D surface model of the convex hull of the astrocyte filaments map. The convex hull was then used to determine the astrocyte domain volume. The convex hull volume was then divided by the cell volume to determine the domain-to-cell volume ratio for each astrocyte.
Machine learning segmentation in Imaris was used to identify and quantify the volume of tdTomato-labeled astrocyte matter within CD68- and Iba1-labeled microglial lysosomes. Training of segmentation algorithm was conducted on 11 randomly selected images within the dataset. Once the creation parameters consistently and accurately generated surfaces of tdTomato inclusions within CD68 and Iba1 double-positive regions, the algorithm was applied to all images. All AI-rendered surfaces within each image were also manually validated to eliminate any false positives. To ensure unbiased analysis, the researcher was blinded to experimental condition during image acquisition and analysis. The number of astrocyte inclusions and the total volume of these inclusions were then quantified for each image.
Nuclei harvesting and fluorescence-activated nuclei sorting
Aldh1l1-GFP mice which underwent hindpaw incision (or isoflurane exposure as a control) at P3 were euthanized at P4, P10, or P24 by sodium pentobarbital overdose. Spinal cords were rapidly dissected in ice-cold RNAse-free 0.01 M PBS, and the left dorsal quadrant of lumbar segments L3–L4 was snap-frozen on dry ice. Tissue was stored at −80°C until the time of homogenization and nuclei sorting.
Isolation of astrocyte nuclei was carried out as previously described (Chamessian et al., 2018; Serafin et al., 2019, 2021). Each spinal cord sample was individually homogenized in 1 ml homogenization buffer (HB) containing the following (in mM): 250 sucrose; 25 KCl; 20 Tris-HCl, pH 8; 5 MgCl2; 1 µM DTT; 1 tablet/10 ml cOmplete Mini EDTA-free Protease Inhibitor Cocktail (Roche); and 40 U/ml RNAsin (Promega). Samples were centrifuged and filtered, and pelleted nuclei were resuspended in a final volume of 300 µl HB and stained with 0.1 µg/ml propidium iodide (PI; Thermo Fisher Scientific).
Fluorescence-activated nuclei sorting (FANS) was performed on a FACS Aria II (BD Biosciences) fitted with a 70 µm nozzle. Sorting gates were established on GFP and PI fluorescence, and GFP+/PI+ nuclei from each sample were sorted into a tube containing 500 µl of lysis buffer RL (Norgen Biotek). For qPCR experiments, nuclei positive for PI but not GFP (presumptive non-astrocyte nuclei) were also sorted into a different tube of lysis buffer, to facilitate comparison between GFP+ and GFP− nuclei from the same spinal cord sample.
Quantitative real-time reverse transcription PCR
Spinal cords from four naive Aldh1l1-GFP mice were harvested and processed as described above and nuclei sorted into two populations: (1) GFP+/PI+ and (2) GFP−/PI+. RNA was extracted from the sorted nuclei using Norgen Total RNA Purification Plus Micro Kit with on-column DNAse treatment (Norgen #48500) and reverse transcribed using SuperScript IV First-Strand Synthesis System (Thermo Fisher Scientific #18091050). qPCR was performed on a QuantStudio 3 (Applied Biosystems) using Taqman gene expression assays for the following probes: Slc1a3 (encoding GLAST; Mm00600697_m1), Aqp4 (Aquaporin-4; Mm00802131_m1), and Hprt (Mm03024075_m1).
RNA sequencing and bioinformatic analysis
RNA was extracted for sequencing from sorted PI+/GFP+ nuclei using the Norgen Single Cell RNA Purification Kit (Norgen #51800), and cDNA libraries were prepared using NEBNext Single Cell/Low Input RNA Library Prep Kit (New England Biosciences #E6420L). Low-input nondirectional RNA sequencing with 2 × 61 bp paired-end reads was carried out by the University of Cincinnati Genomics, Epigenomics and Sequencing Core using an Illumina NextSeq 2000, with a target depth of 60 M reads. At least 22 M reads were achieved for each sample, with a mean read depth of 30 M reads across all 24 samples. The resulting FASTQ sequencing files were aligned and quantified using BaseSpace Sequence Hub app RNA-Seq Alignment v2.0.2 (Illumina), which utilizes STAR (Dobin et al., 2013) to align to the GRCm38 Mus musculus genome assembly (mm10, UCSC Genome Browser) and salmon (Patro et al., 2017) to compute gene-level transcript quantification. Normalization, principal component analysis (PCA), and differential gene expression analysis were then carried out using DESeq2 (Love et al., 2014) within the R statistical programming environment.
Differential gene expression analysis
Effects of P3 incision were analyzed separately at each time point. Differential gene expression analysis between incision and naive samples was performed using the Wald test (model, ∼sex + incision) to identify genes having adjusted p value <0.05. Visualizations of injury-dependent DEGs at each time point were produced using EnhancedVolcano (Blighe et al., 2024), ggplot2 (Wickham, 2016), and pheatmap (Kolde, 2018).
Age-dependent differential gene expression was analyzed using naive samples only across all three time points (P4, P10, and P24). The likelihood ratio test (LRT; full model, ∼sex + age; reduced model, ∼sex) was used to identify age-dependent DEGs having adjusted p value <0.001. Clustering of age-dependent DEGs was performed on variance-stabilizing transformed expression data using DEGreport package (Pantano et al., 2024). Visualizations of age-dependent DEGs were produced using DEGreport and ggplot2.
For each differential gene expression analysis described above, data was prefiltered to exclude any gene for which at least 10 reads were not detected in a minimum of four samples (i.e., the size of the smallest experimental group). Sex was included as a cofactor for all models, but sex differences were not investigated directly. The false discovery rate (Benjamini–Hochberg; FDR) was used to correct for multiple hypothesis testing, producing adjusted p values (p-adj). After hypothesis testing was performed, the resulting dataset was subjected to additional quality control, excluding genes which did not have mean normalized expression ≥20 in at least one experimental group and for which the coefficient of variation (CoV) was >65%.
Gene Ontology term enrichment analysis
Gene Ontology (GO) term enrichment analysis was performed on identified DEGs using gProfiler (Reimand et al., 2007). Data sources were restricted to only Gene Ontology (GO: Molecular Function, GO: Cellular Compartment, GO: Biological Process) terms. Visualizations of the resulting enriched GO terms with p-adj < 0.01 were produced using ggplot2.
Neonatal astrocyte receptor–ligand interactome
To identify the populations of dorsal root ganglion (DRG) neurons which may communicate with neonatal spinal astrocytes, we leveraged two existing databases which form the foundational framework of (1) the ligand–receptor interactome platform developed by Wangzhou et al. (2021) and (2) CellChat (Jin et al., 2021), as well as previously published single-cell RNA-seq data obtained from P5 mouse DRG (Sharma et al., 2020). First, our gene expression data from naive P4 astrocytes were subsetted to include any gene annotated as a receptor in either interactome database, for which mean normalized expression was at least 20 across experimental samples. These astrocyte receptor genes were then classified as either “Catalytic Receptors,” “G-protein-coupled Receptors (GPCRs),” “Ligand-gated Ion Channels,” or “Nuclear Hormone Receptors” using the complete target and family list provided by IUPHAR (Harding et al., 2022). Any receptor gene appearing in either interactome database but not assigned to any of the above families by IUPHAR was designated as “Other Receptor.” Next, the interactome databases were used to compile a list of ligands for each astrocyte receptor gene.
To identify which subtype(s) of DRG neurons might express these ligands, we obtained single-cell RNA-seq data with cell-type annotations at P5, which was the closest age available in the developmental atlas published in Sharma et al. (2020) (available for download at https://kleintools.hms.harvard.edu/tools/springViewer_1_6_dev.html?datasets/Sharma2019/all). The count matrix was imported into Seurat version 5.1.0 (Hao et al., 2024), log-transformed, and rescaled, and each cell's identity (cluster) was assigned using the annotation previously established (Sharma et al., 2020). Marker genes expressed by each cluster of DRG neurons were obtained by using the FindAllMarkers function; a gene was considered a “Marker Gene” for a given cluster if (1) the minimum percentage of cells expressing the gene in the given cluster was at least 25% and (2) the minimum difference in average Log2FC for cells in the cluster compared with all other cells was at least 0.25. These thresholds are deliberately quite relaxed, since the aim was to identify all clusters which substantively express each ligand gene. Finally, these data were merged to create a spreadsheet with each astrocyte receptor, all of its putative ligands, and all of the DRG clusters for which each putative ligand is a Marker Gene. The top 10 most highly expressed astrocyte receptors in each of the five categories noted above are provided as individual supplemental tables (Tables S1–S5), while the complete merged spreadsheet is provided in Dataset S1.
In situ hybridization and immunohistochemistry
Aldh1l1-GFP pups were euthanized at either P4 (in situ hybridization, ISH) or P24 (immunohistochemistry, IHC) with sodium pentobarbital. Spinal cord lumbar enlargements were dissected out and postfixed for 1.5 h (P4) or 2 h (P24) in 4% PFA dissolved in RNAse-free 0.01 M PBS. Spinal cords were transferred to 30% sucrose dissolved in RNAse-free PBS and stored at 4°C until the time of cryosectioning. Transverse frozen sections (14 µm) were mounted on SuperFrost slides (Thermo Fisher Scientific). For experiments comparing tissue from incised versus naive pups, slides were prepared such that each slide contained sections of all biological replicates of both experimental groups.
For Sox9 and Sox10 IHC, primary antibodies were diluted 1:500 (Sox9; R&D Systems #AF3075) or 1:1,000 (Sox10; Cell Signaling Technology #69661) and applied overnight at 4°C. Secondary antibodies conjugated to Alexa Fluor 594 (Thermo Fisher Scientific #A11012 or #A11058) were diluted 1:1,000 and applied for 1 h at room temperature. DAPI (0.1 µg/ml) was used as a counterstain before coverslipping with Vectashield mounting medium (VectorLabs #H-1000). Endogenous Sun1-GFP was not amplified with antibody staining for these experiments. Images of the left spinal dorsal horn were captured on a Keyence inverted fluorescent microscope BZ-X800 (Keyence). Z-stack images were taken under 10× magnification with a z-separation of 0.7 µm and then projected as full-focus images for analysis.
ISH for Fgf14 (ACDBio #462801) was performed using RNAscope Multiplex Fluorescent v2 kit (ACDBio) following manufacturer's directions. TSA Plus Cy3 system (Akoya Biosciences #NEL744001KT) was used to develop fluorescent signal. After ISH, IHC was performed to retrieve the endogenous Sun1-GFP signal. Primary chicken anti-GFP (Thermo Fisher Scientific #A10262) was used at 1:1,000, and secondary Alexa Fluor 488 goat anti-chicken IgY (Thermo Fisher Scientific #A11039) was used at 1:1,000. Finally, DAPI counterstain was applied, and slides were coverslipped with Vectashield mounting medium. Images of the left (i.e., ipsilateral to hindpaw incision) spinal dorsal horn were captured as Z-stack images taken under 20× magnification with a z-separation of 0.5 µm. These were projected as full-focus images for analysis. Exposure parameters were kept constant across all samples to facilitate consistent downstream processing.
Image processing for ISH and IHC experiments
Quantification of fluorescent images was performed in QuPath (Bankhead et al., 2017; ACDBio Technical Note #MK 51-154/Rev A/Date 12/21/2020), using a combination of automated detection and manual classification steps. First, nuclei were identified and segmented on DAPI fluorescence (background radius, 5 µm; sigma, 2; minimum area, 10 µm2; maximum area, 200 µm2; threshold, 10; cell expansion, 1 µm). DAPI-stained nuclei exhibiting Sun1-GFP fluorescence were manually selected and classified as “Sun1-GFP” while non-GFP nuclei were classified as “negative.” For IHC experiments investigating colocalization of transcription factors Sox9 or Sox10 with Sun1-GFP nuclei, the maximum red fluorescent intensity of the identified cell nucleus was used to classify cells as “positive” or “negative” for a given target.
For ISH experiments, Fgf14 RNAscope puncta were identified within “Sun1-GFP” nuclei using subcellular object detection on Cy3 fluorescence (detection threshold, 70; expected spot size, 1 µm2; minimum spot size, 0.25 µm2; maximum spot size, 4 µm2; include clusters, yes). Results were exported for each image, and an R script was used to extract and process relevant data as follows: the number of puncta (single spots ≤4 µm2) and clusters (larger spots >4 µm2) contained within Sun1-GFP+ nuclei were counted for each image. Single spots were counted as 1 punctum while the number of puncta contained within clusters was computed by dividing total cluster area by the mean single-spot area for that image. This total number of puncta was then divided by the number of Sun1-GFP+ nuclei contained within the image, yielding a mean value of dots per Sun1-GFP+ nucleus for that image. Prevalence of Fgf14+ astrocyte nuclei was computed as the number of Sun1-GFP+ nuclei containing at least one Cy3 punctum divided by the total number of Sun1-GFP+ nuclei in that section.
Data and statistical analysis
The D'Agostino and Pearson’s test or the Shapiro–Wilk test for normality (GraphPad Prism 10) were used to determine whether nonparametric tests or data transformation were necessary. Metrics of astrocyte size and complexity were analyzed using two-way analysis of variance (ANOVA) with injury (naive vs incision) and age (P4 vs P10 vs P24) as factors. The Šidák's multiple-comparisons test was used for post hoc analyses (p-adj refers to p values generated by post hoc analyses). 3D Sholl measurements were analyzed with repeated measures (RM) two-way ANOVA with either age and distance from center as factors or injury and distance from center as factors, with the Šidák's multiple-comparisons test used for post hoc analyses. Sample sizes for these analyses are as follows: P4 naive, n = 94 cells from 4 mice; P4 incision, n = 163 cells from 5 mice; P10 naive, n = 125 cells from 5 mice; P10 incision, n = 183 cells from 4 mice; P24 naive, n = 122 cells from 8 mice; and P24 incision, n = 172 cells from 8 mice. Both male and female mice were included in all experimental groups, but data were pooled across sexes for the measurements of astrocyte size and complexity given that the sample sizes were insufficient to rigorously analyze the data using sex as an independent variable.
For qPCR experiments, the relative quantification was conducted with the 2−ΔΔCT method with Hprt as a reference gene. Each data point represents the mean of three technical replicates per biological replicate. For IHC experiments, four biological replicates and three to four tissue sections per replicate were analyzed. For ISH experiments, four biological replicates for each experimental group (i.e., incised vs naive) and four images per biological replicate were analyzed. Incision versus naive groups were compared with a nested unpaired t test. Each data point represents one tissue section. The number of astrocyte inclusions and volume of engulfed astrocyte matter within microglial lysosomes were analyzed using two-way ANOVA with injury and sex as factors. Each data point represents the arithmetic mean of 3–4 regions of interest taken from 2–3 tissue sections per animal with n = 4–5 animals per group (i.e., naive female, incision female, naive male, incision male). Data are expressed as arithmetic mean ± SEM. All studies were conducted by investigators blinded to the treatment group during data collection and data analyses.
Data availability
All data are available upon written request submitted to the corresponding author. Raw sequencing data and metadata were deposited at GEO under accession number GSE289540. The results of all differential gene expression analyses and GO term enrichment analyses are provided as supplemental content. Additionally, we have created an interactive website for exploration of our dataset: https://cellsera.shinyapps.io/AS_Dev/. R code for data analysis and visualization is available in an Open Science Framework (OSF) repository at DOI 10.17605/OSF.IO/XPK92.
Results
High-resolution confocal imaging of individual spinal dorsal horn astrocytes for three-dimensional morphological analysis
Astrocyte function is highly dependent on their complex structure (Torres-Ceja and Olsen, 2022); however, nothing is currently known about how the morphology of dorsal horn astrocytes changes during postnatal development or about how this is impacted by injury. We therefore sought to assay the morphological properties of developing spinal astrocytes using high-resolution confocal microscopy and three-dimensional (3D) image analysis techniques. Given that much of synaptogenesis and astrocyte differentiation occurs within the first 3 weeks of postnatal development (Chung et al., 2015), we focused our analysis on this critical time period. Using a low dose of tamoxifen to induce sparse labeling of astrocytes in Aldh1l1-tdTOM mice, we imaged individual astrocytes in mice which underwent a hindpaw incision at P3, and naive littermate controls, at P4, P10, or P24 (Fig. 1a). For accurate 3D morphological analysis, images capturing the full depth of fluorescently labeled astrocytes were taken within the medial third of the left spinal dorsal horn (i.e., ipsilateral to the incision) in cleared tissue (Fig. 1b, Movie 1). As anticipated based on prior reports of astrocyte morphology in the developing hippocampus (Bushong et al., 2004), we observed clear increases in spinal astrocyte complexity between P4 and P24 (Fig. 1c).
High-resolution three-dimensional imaging of individual spinal dorsal horn astrocytes using confocal microscopy. a, Schematic depicting the experimental design. Aldh1l1-tdTOM mice were given one intragastric injection of low-dose (0.25 mg/kg) tamoxifen at P1. At P3, mice were given a hindpaw incision or anesthesia only as naive controls. Spinal cords were harvested 1 d (P4), 1 week (P10), or 3 weeks (P24) after hindpaw incision. The lower lumbar enlargement was then dissected and transverse sections of 100 µm thickness were generated. These sections were cleared using the EZ Clear protocol and then imaged using the Yokogawa SoRa spinning disk confocal. b, Left panel, Representative image of transverse section from Aldh1l1-tdTOM mouse with sparse labeling. Right panels, Representative images of a single tdTomato-labeled astrocyte. Scale bars, 10 µm. Top left, 3D rendering; top right, x–y plane; bottom left: x–z plane; bottom right: y–z plane. c, Representative images of tdTomato-labeled dorsal horn astrocytes from P4 (left), P10 (middle), and P24 (right).
Neonatal hindpaw incision alters developing astrocyte morphology
The use of Imaris software to quantify the morphological properties of dorsal horn astrocytes (Fig. 2a–d, Movie 2) demonstrated that neonatal incision evokes age-dependent changes to cell surface area (n = 94–183 cells from 4–8 mice per group; injury × age interaction: p < 0.0001; two-way ANOVA; Fig. 2e), cell volume (p < 0.0001; Fig. 2f), and domain volume (p < 0.0001; Fig. 2g). At P4, spinal astrocytes in the P3 incision group have significantly increased surface area (p-adj = 0.002; Šidák's multiple-comparisons test; Fig. 2e). However, astrocytes in mice with hindpaw incision have decreased surface area at P10 and P24 (p-adj < 0.0001). Similarly, P3 incision increased astrocyte cell volume at P4 (p-adj = 0.0024; Fig. 2f) but reduced the volume at P10 and P24 (p-adj < 0.0001) compared with naive controls. Neonatal incision also increased astrocyte domain volume at P4 (p-adj = 0.0063; Fig. 2g) but decreased domain volume at P24 (p-adj < 0.0001), with no significant changes observed at P10 (p-adj = 0.1321). Our analysis of the ratio of astrocyte domain volume-to-cell volume, which measures how far the astrocyte processes extend from the soma relative to the overall density of astrocyte processes, revealed age-dependent effects of P3 incision (injury × age interaction: p = 0.0033, two-way ANOVA; Fig. 2h). P3 incision increased the astrocyte domain volume-to-cell volume ratio at P10 (p-adj = 0.0032; Fig. 2h) but did not change the ratio at P4 (p-adj = 0.09485) or P24 (p-adj = 0.5328).
Neonatal hindpaw incision alters astrocyte morphology in the developing dorsal horn. a, Representative image of tdTomato-labeled astrocyte. b, Representative image of Imaris generated surface map of astrocyte used to quantify cell volume and surface area. c, Representative image of Imaris generated filaments map of astrocyte used to create convex hull. d, Representative image of convex hull used to calculate astrocyte domain volume. e, P3 hindpaw incision increases astrocyte surface area at P4 but reduces surface area at P10 and P24 [n = 94–183 cells from 4–8 mice per group (i.e., P4 naive, P4 incision, P10 naive, P10 incision, P24 naive, P24 incision); injury × age interaction: F(2,853) = 24.19, p < 0.0001; two-way ANOVA; **p-adj = 0.002, ****p-adj < 0.0001, Šidák's multiple-comparisons test]. f, P3 incision increases astrocyte volume at P4 but decreases volume at P10 and P24 (injury × age interaction: F(2,853) = 20.39, p < 0.0001, two-way ANOVA; **p-adj = 0.0024, ****p-adj < 0.0001). g, P3 incision increases astrocyte domain volume at P4 but decreases domain volume at P24 (injury × age interaction: F(2,853) = 14.77, p < 0.0001, two-way ANOVA; **p-adj = 0.0063, ****p-adj < 0.0001). h, Incision increases the domain-to-cell volume ratio at P10 (injury × age interaction: F(2,853) = 5.752, p < 0.0033, two-way ANOVA; **p-adj = 0.0032, Šidák's multiple-comparisons test).
To determine how incision affects astrocyte complexity, we used Imaris filaments analysis to map individual astrocyte processes (Fig. 3a,b; Movie 2). Using 3D Sholl analysis (Fig. 3c), we determined the number of intersections between astrocyte processes and spheres of defined radius from the soma in both naive and P3-incised animals at P4, P10, and P24. As anticipated, we found that 3D Sholl intersections increase with age in naive mice [age × distance interaction: p < 0.0001; repeated-measures (RM) two-way ANOVA; Fig. 3d]. Similarly, 3D Sholl intersections increased between P4 and P24 in mice with neonatal hindpaw incision (age × distance interaction: p < 0.0001; Fig. 3e). To more directly examine the effects of P3 incision on astrocyte complexity, we compared the number of Sholl intersections at each developmental time point. We found that spinal astrocytes in the incised group were more highly ramified at P4 when compared with naive controls (distance × injury interaction: p < 0.0001; RM two-way ANOVA; Fig. 3f) at distances of 20–45 µm away from the soma (p-adj < 0.05, Šidák's multiple-comparisons test). Interestingly, the opposite was found at P10, in which astrocytes in the P3-incised group had reduced Sholl intersections compared with naive controls at distances of 10–40 µm from the soma (distance × injury interaction: p < 0.0001; p-adj < 0.01; Fig. 3g). We found similar differences at P24, where astrocytes from mice with neonatal incision had decreased Sholl intersections at distances of 20–45 µm away from the soma (distance × injury interaction: p < 0.0001; p-adj < 0.05; Fig. 3h). Our analysis of astrocyte process branch points, terminal points, and total filament area all demonstrated the same pattern of increased complexity at P4, and decreased complexity at P10 and P24, following neonatal incision (age × injury interaction: p < 0.0001, two-way ANOVA; p-adj < 0.05; Šidák's multiple-comparisons test; Fig. 3i–k). Taken together, these results suggest that the morphological plasticity of developing astrocytes is impacted by neonatal incision in an age-dependent manner.
Neonatal incision yields age-dependent changes to astrocyte complexity in the dorsal horn. a, Representative image of tdTomato-labeled astrocyte at P24. b, Representative image of Imaris generated filaments map used to quantify metrics of astrocyte complexity. c, Cartoon depicting Sholl analysis method. d, Astrocyte 3D Sholl intersections increase with age in naive mice [age × distance from center interaction: F(22,3718) = 51.58, p < 0.0001; repeated-measures (RM) two-way ANOVA]. e, 3D Sholl intersections increase with age in mice with P3 hindpaw incision [age × distance interaction: F(22,5654) = 44.47, p < 0.0001; RM two-way ANOVA]. f, Astrocyte Sholl intersections at P4 increase after P3 incision (distance × injury interaction: F(11,2805) = 9.033, p < 0.0001; RM two-way ANOVA; *p-adj < 0.05, **p-adj < 0.01, ***p-adj < 0.001, Šidák's multiple-comparisons test). g, P10 astrocyte Sholl intersections decrease after neonatal incision (distance × injury interaction: F(11,3366) = 25.10, p < 0.0001; RM two-way ANOVA; **p-adj < 0.01, ***p-adj < 0.001, ****p-adj < 0.0001, Šidák's multiple-comparisons test). h, P24 astrocyte Sholl intersections decrease after P3 incision (distance × injury interaction: F(11,3256) = 11.89, p < 0.0001; RM two-way ANOVA; *p-adj < 0.05, **p-adj < 0.01, ****p-adj < 0.001, Šidák's multiple-comparisons test). i, P3 incision increases astrocyte branch points at P4 but decreases branch points at P10 and P24 (injury × age interaction: F(2,858) = 28.07, p < 0.0001; two-way ANOVA; *p-adj = 0.0169; ****p-adj < 0.0001, Šidák's multiple-comparisons test). j, P3 incision increases astrocyte terminal points at P4 but decreases terminal points at P10 and P24 (injury × age interaction: F(2,858) = 28.11, p < 0.0001, two-way ANOVA; *p-adj = 0.0203, ****p-adj < 0.0001). k, P3 incision increases total astrocyte filament area at P4 but decreases filament area at P10 and P24 (injury × age interaction: F(2,857) = 26.21, p < 0.0001; ****p-adj < 0.0001).
Isolation of spinal dorsal horn astrocyte nuclei for transcriptional analysis
To investigate the effects of neonatal incision on the transcriptional signature of astrocytes in the developing dorsal horn, we performed population-level RNA sequencing on astrocyte nuclei selectively isolated from the L3–L4 dorsal horn at P4, P10, and P24, corresponding to post-incision time points 1, 7, and 21 d. Aldh1l1-Cre/ERT2 mice (Srinivasan et al., 2016) were crossed to Cre-dependent Gt(ROSA)26-CAG-Sun1-sfGFP mice (Mo et al., 2015) to generate Aldh1l1-GFP mice in which astrocyte nuclei were labeled with a fluorescent Sun1-GFP fusion protein following intragastric administration of 50 mg/kg tamoxifen at P1 and P2 (Fig. 4a). IHC for the astrocyte marker Sox9 (Sun et al., 2017) was performed on tissue sections harvested at P3 (Fig. 4b) which showed that tamoxifen-induced Cre-dependent expression was widespread, with 87.64% ± 1.37% of Sox9+ nuclei coexpressing Sun1-GFP, and selective for astrocytes, with 92.91 ± 0.73% of Sun1-GFP+ nuclei coexpressing Sox9 (Fig. 4c). The Isolation of Nuclei Tagged in specific Cell Types (INTACT) method followed by FANS enabled recovery of putative Sun1-GFP+ astrocyte nuclei from each animal (Fig. 4d). In preliminary experiments using tissue harvested at P4 from naive pups only, GFP− nuclei were also recovered from the same samples during cell sorting, and the expression of astrocyte marker genes Slc1a3 (encoding GLAST) and Aqp4 (encoding Aquaporin-4) was highly enriched in GFP+ nuclei compared with GFP− nuclei recovered from the same animal (p < 0.0001, n = 4, unpaired t test; Fig. 4e).
Isolation of RNA from spinal astrocyte nuclei. a, Intragastric administration of 50 mg/kg tamoxifen at P1 and P2 induces widespread Sun1-GFP nuclear labeling (green) of Aldh1l1-Cre-expressing cells in the spinal cord at P4. Scale bar, 100 µm. b, Most Sun1-GFP-labeled nuclei (green) coexpress the astrocyte marker Sox9 (magenta). Scale bar, 100 µm. c, Quantification of Sox9 and Sun1-GFP coexpression in the left dorsal quadrant shows that 87.64% ± 1.37% of Sox9+ cells coexpress Sun1-GFP and 92.91 ± 0.73% of Sun1-GFP+ nuclei coexpress Sox9. d, Gating strategy to selectively isolate Sun1-GFP+ nuclei using FANS. e, Expression of astrocyte marker genes Slc1a3 (encoding GLAST) and Aqp4 (Aquaporin-4) is highly enriched in GFP+ versus GFP− sorted nuclei from the same animal. Expression was analyzed by qPCR and normalized to reference gene Hprt using the 2−ΔΔCT method. Slc1a3: n = 4, t = 12.11, df = 6, ****p < 0.0001, unpaired t test; Aqp4: n = 4, t = 13.90, df = 6, ****p < 0.0001, unpaired t test. f, Principal component analysis plot of gene expression data from 24 samples of sequenced astrocyte nuclei. Red, P4; green, P10; blue, P24. Circles, P3 naive (isoflurane exposure control); triangles, P3 incision.
The INTACT/FANS approach was used to selectively isolate astrocyte nuclei from 24 pups: four naive and P3-incised pups were harvested at each of three developmental time points (P4, P10, P24), and then low-input RNA sequencing was performed. PCA of variance-stabilizing transformed gene expression data demonstrated notable segregation by age (Fig. 4f). Segregation by sex was also evident within each age group (Fig. S1), and thus sex was included as a cofactor in gene expression modeling for subsequent differential gene expression analysis.
Naive spinal astrocytes exhibit significant age-dependent transcriptional changes in early life
Since little is known about the postnatal development of spinal astrocytes, we first sought to elucidate the age-dependent changes to astrocyte gene expression. Differential gene expression analysis across P4, P10, and P24 in naive animals revealed several thousand age-dependent DEGs (Dataset S1). Therefore, a more stringent alpha of 0.999 (p-adj < 0.001) was employed to focus on the most robust findings. A total of 2,027 age-dependent DEGs (Fig. 5a) were clustered according to their developmental pattern of expression across time using R package DEGreport, which employs divisive hierarchical clustering of variance-stabilizing transformed gene expression data based on the Euclidean gene–gene distance (Pantano et al., 2024). This analysis identified four groups, and then GO term enrichment analysis was performed on the list of DEGs for each group (Dataset S1). Enriched GO terms which were unique to each group were then plotted to elucidate the processes and functions unique to that group. Group 1 genes exhibit generally increasing expression across time (Fig. 5b, left, red; 968 genes), including canonical astrocyte marker genes Slc1a2, encoding GLT-1 (p-adj = 3.07 × 10−10), and Kcnj10, encoding Kir4.1 (p-adj = 2.77 × 10−15). Enriched GO terms unique to this group are related to cellular metabolism (Fig. 5b, right), consistent with the increasing energy demand required for growth. Group 2 genes, conversely, exhibit generally decreasing expression with age (Fig. 5c, left, green; 614 genes), including the astrocyte-selective glutamate transporter GLAST (encoded by Slc1a3; p-adj = 1.04 × 10−32) and several Hox family transcription factors which are instrumental in homeotic differentiation and patterning along the anterior–posterior axis during development (reviewed in Hubert and Wellik, 2023). This is reflected in the enriched GO terms unique to this group, which are related to embryonic development, macromolecule metabolism, DNA binding and transcription, and enzymatic activity (Fig. 5c, right). The generally decreasing expression of genes belonging to this cluster could reflect a tapering off of the highly dynamic transcriptional environment of early life, as the astrocytes settle into their mature state. Groups 3 and 4 contain fewer DEGs than Groups 1 and 2 and are characterized by a “U-shaped” pattern of gene expression, wherein a local minimum exists at P10. In Group 3, the larger reduction in expression occurs between P4 and P10, with minimal recovery at P24 (Fig. 5d, left, blue; 191 genes). The few enriched GO terms unique to this group are related to mitosis (Fig. 5d, right), consistent with a prior study demonstrating that mitotic proliferation of astrocytes primarily occurs between embryonic day 15 and P3 (Tien et al., 2012). Finally, Group 4 genes exhibit a slight decrease between P4 and P10, followed by a larger increase at P24 (Fig. 5e, left, purple; 254 genes). Enriched GO terms unique to this group are focused on extracellular export, either via secretion or vesicle release, and trans-synaptic signaling (Fig. 5e, right). This could reflect increased cellular communication with neurons and other cell types beginning in early adolescence.
Naive spinal astrocytes exhibit significant age-dependent transcriptional changes in early life. a, 2,027 age-dependent DEGs were identified (likelihood ratio test, p-adj < 0.001) and clustered into four groups. b–e, Age-dependent DEGs clustered by expression pattern over P4 to P24. Left, Z-score of normalized expression of all genes assigned to each group. Right, Gene Ontology (GO) term enrichment analysis results showing enriched GO terms unique to the group shown on left. All enriched GO terms have p-adj < 0.01, and the top 10 terms in each GO category (ranked by p-adj) are plotted for each group. GO:BP, Gene Ontology:Biological Process, red; GO:CC, GO:Cellular Component, blue; GO:MF, GO:Molecular Function, green. Circle size indicates the number of age-dependent DEGs in the given group which intersect with the list of genes assigned to the given GO term.
Potential signaling mechanisms underlying interactions between DRG neurons and astrocytes within the neonatal dorsal horn
To investigate candidate signaling pathways by which sensory neurons might communicate with astrocytes within the immature dorsal horn, we identified genes that were (1) expressed by P4 astrocytes under naive conditions and (2) annotated as a receptor in either the ligand–receptor interactome platform established by Wangzhou et al. (2021)or the CellChat platform (Jin et al., 2021). These astrocytic genes were then classified as “Catalytic Receptors,” “GPCRs,” “Ligand-gated Ion Channels,” “Nuclear Hormone Receptors,” or “Other Receptors” (see Materials and Methods). Next, the above two ligand–receptor interactome databases were used to establish a list of ligands for each identified astrocyte receptor. An available single-cell RNA-seq atlas of the P5 DRG (Sharma et al., 2020) was then employed to elucidate the transcriptional subpopulation(s) of neonatal sensory neurons which express the relevant ligands.
With regard to the ligand-gated ion channels, multiple subunits of NMDA receptors (encoded by Grin1, Grin2a, Grin2b, and Grin2d), as well as the GluA1 subunit of AMPA receptors (Gria1), are expressed by neonatal spinal astrocytes and predicted to bind the leucine-rich repeat and fibronectin type III domain containing 1 protein (LRFN1) encoded by Lrfn1 (Table S1). Notably, Lrfn1 is expressed by multiple subtypes of sensory neurons in the P5 DRG, including (but not limited to) three transcriptional subpopulations of CGRP+ neurons (Table S1). Neonatal CGRP-expressing neurons could also signal to developing astrocytes through nuclear hormone receptors such as those encoded by Nr2f1, Ppara, Nr1h2, or Nr1h3 via their expression of the corresponding ligands FAT atypical cadherin 1 (encoded by Fat1), nuclear receptor coactivator 6 (Ncoa6), or endothelial differentiation-related factor 1 (Edf1; Table S2). When considering the classification of catalytic receptors, fibroblast growth factor (FGF) receptors 3 and 1 (encoded by Fgfr3 and Fgfr1, respectively) were highly expressed by neonatal dorsal horn astrocytes and may be activated by multiple subtypes of nociceptors (including both CGRP+ neurons and nonpeptidergic nociceptors), as well as pruriceptors (Sst+ neurons), via signaling initiated by a host of ligands within the FGF family (Table S3). Similarly, the integrin beta 1 subunit (CD29; encoded by Itgb1) and integrin alpha V (i.e., CD51; Itgav) are highly expressed by dorsal horn astrocytes at P4 and are predicted to respond to myriad ligands (such as various isoforms of metallopeptidases, collagens, and laminins) that are widely distributed across the transcriptional subpopulations of neonatal DRG neurons (Table S3).
The GPCRs highly expressed at the mRNA level by spinal astrocytes at P4 include endothelin receptor type B (encoded by Ednrb), which is predicted to bind the G-protein subunit Gαi2 (encoded by Gnai2), as well as well-studied receptors such as mGluR5 (Grm5), calcitonin receptor-like receptor (Calcrl), and GPR37L1 (Gpr37l1; Table S4). Notably, GPR37L1 is selectively localized to astrocytes within the dorsal horn where it constitutively suppresses astrogliosis and limits pain sensitivity (Xu et al., 2025b). This receptor can be activated by maresin 1 (Bang et al., 2024) and Neuroprotectin-1 (Xu et al., 2025a), members of the specialized pro-resolving mediator (SPM) family of bioactive lipids derived from the omega-3 fatty acid docosahexaenoic acid (DHA). Interestingly, the ligand–receptor interactome databases used here suggest that GPR37L1 in immature spinal astrocytes could also be activated by neonatal sensory neurons via signaling involving metabolites of the glycoprotein prosaposin (encoded by the Psap gene; Table S4) that are important for sphingolipid hydrolysis. Finally, an examination of the potential ligand–receptor interactions within the “Other Receptors” category (Table S5) revealed that neonatal DRG neurons are poised to signal to developing spinal astrocytes expressing neurexin I (encoded by Nrxn1) via a wide range of ligands including multiple isoforms of neuroligins (Nlgn1-3), rabphilin 3a (Rph3a), and synaptotagmin-like 3 (Sytl3). While expression of Nlgn1-3 and Rph3a have been documented in both nociceptive and nonnociceptive DRG neurons at P5, Sytl3 seems exclusively expressed by CGRP+ neurons at this age (Table S5; Sharma et al., 2020).
Neonatal hindpaw incision induces differential gene expression in spinal dorsal horn astrocytes acutely after injury
We then compared the transcriptional profile of astrocytes in the naive versus P3 incision groups using differential gene expression analysis at each time point. Seventy-six DEGs having adjusted p value <0.05 were identified 1 d after incision, at P4 (Fig. 6a, Dataset S1). Notable upregulated genes include Thbs1, encoding thrombospondin-1, and actin cytoskeleton genes Acta1, Acta2, Actg2, and Tpm2 (Fig. 6b). These genes are especially interesting considering the increased astrocyte complexity, branch points, and endpoints revealed by the morphological studies (Figs. 2, 3). Downregulated genes include Dync2h1, encoding a member of the dynein-2 complex responsible for intracellular transport; Daam1, which interacts with the actin cytoskeleton; and Fgf14, encoding an intracellular (i.e., nonsecreted) FGF. FGF signaling in general contributes to neural stem cell differentiation into astrocytes (Savchenko et al., 2019) and oligodendrocytes (Furusho et al., 2011) in the developing CNS but appears to suppress astrocyte activation in the mature CNS (Kang et al., 2014). While Fgf14 is understudied in astrocytes, genetic deletion of this isoform results in impaired neurogenesis and aberrant synaptic integration and plasticity in the adult dentate gyrus (Alshammari et al., 2016). GO term enrichment analysis using the 76 DEGs identified at P4 revealed enrichment for terms related to cell motility and migration, morphogenesis, and structural terms related to the extracellular matrix and cytoskeleton (Fig. 6c, Dataset S1).
Neonatal hindpaw incision induces differential gene expression in dorsal horn astrocytes at P4. a, 76 injury-dependent DEGs were identified 1 d following neonatal hindpaw incision (Wald test, p-adj < 0.05). Each column in heat map represents one biological replicate; color indicates z-score of normalized expression across each row. b, Volcano plot showing injury-induced fold-change in gene expression. c, GO term enrichment analysis was performed using 76 injury-dependent DEGs as the input. Enriched GO terms have p-adj < 0.01. GO:BP, Gene Ontology:Biological Process, red; GO:CC, GO:Cellular Component, blue; GO:MF, GO:Molecular Function, green. Top 10 enriched GO terms (ranked by p-adj) are plotted for each GO category. Size of circle indicates the number of injury-dependent DEGs which intersect with the list of genes assigned to the given GO term.
To validate select results from the RNA sequencing, mRNA transcripts encoding Fgf14 were quantified using ISH on tissue sections of the spinal L3-L4 dorsal horn harvested from naive or P3-incised Aldh1l1-GFP mice at P4 (Fig. 7a). Fluorescent puncta indicating Fgf14 mRNA transcripts were readily detected in most astrocyte nuclei. Although the prevalence of astrocyte nuclei expressing Fgf14 was unaffected by neonatal incision (n = 4 naive, n = 4 incision, p = 0.114, nested unpaired t test; Fig. 7b), the mean number of fluorescent puncta indicating Fgf14 transcripts was significantly lower in astrocyte nuclei from neonatally incised mice compared with naive mice (p = 0.040, nested unpaired t test; Fig. 7c–e).
Neonatal hindpaw incision significantly reduces the number of Fgf14 mRNA transcripts within spinal astrocyte nuclei at P4. a, ISH reveals widespread expression of Fgf14 mRNA in the spinal dorsal horn, including within astrocyte nuclei labeled by Sun1-GFP. Scale bar, 50 µm. b, Incision does not alter the percentage of astrocytes which express Fgf14 mRNA (t = 1.846, df = 6, p = 0.1144, nested unpaired t test, n = 4 biological replicates in each group; each data point represents 1 tissue section). c, d, QuPath was used to detect DAPI-stained nuclei (blue) which were then manually assessed for Sun1-GFP fluorescence (green; Sun1-GFP+ nuclei are described by yellow annotations). Red puncta indicating Fgf14 mRNA transcripts were then identified as subcellular objects (red) within Sun1-GFP nuclei. Scale bars, 10 µm. e, The number of Fgf14 mRNA transcripts as indicated by fluorescent in situ puncta is significantly reduced in Sun1-GFP astrocyte nuclei 1 d after P3 hindpaw incision (t = 2.612, df = 6, *p = 0.040, nested unpaired t test, n = 4 biological replicates in each group; each data point represents one tissue section).
Interestingly, minimal changes in gene expression were detected 7 or 21 d after P3 incision, with only two incision-induced DEGs identified at P10 (Fig. 8a) and 8 DEGs detected at P24 (Fig. 8b). The DEGs upregulated at P10 are Kndc1 and Kcnk2. The former encodes the Ras guanine nucleotide exchange factor v-KIND, which suppresses dendrite growth when overexpressed in vitro (Huang et al., 2007). The upregulation of this gene is interesting considering the reduction in astrocyte cell volume and complexity evident at P10. The latter encodes the potassium two-pore domain channel TREK-1, which mediates fast glutamate release in astrocytes in response to GPCR activation/G-protein binding (Woo et al., 2012). Of the eight DEGs identified at P24, most exhibited less than a twofold change in expression between P3-incised and naive mice. Among these genes is Eml1, which encodes the microtubule-associated protein (MAP) EMAP-like 1. MAPs contribute to microtubule stabilization and organization during growth (Conde and Cáceres, 2009).
Minimal persistent changes in gene expression within dorsal horn astrocytes after neonatal incision. a, Volcano plot showing two injury-dependent DEGs at 1 week after P3 hindpaw incision. b, Volcano plot showing eight injury-dependent DEGs at 3 weeks following P3 hindpaw incision. DEGs were identified by Wald test, p-adj < 0.05.
We also noted that Mobp and Mbp, which both exhibit a minimal but statistically significant downregulation at P24 following P3 incision, are canonical marker genes of oligodendrocytes. Although it is known that a minority of cells labeled by the Aldh1l1-Cre/ERT2 mouse line in adulthood are oligodendrocytes (Srinivasan et al., 2016), the prevalence of oligodendrocytes in dorsal horn nuclei labeled by this mouse line during early life was an important question to address in order to interpret our RNA-seq results. Immunostaining for the oligodendrocyte marker Sox10 (Kuhlbrodt et al., 1998) was performed on tissue sections from naive pups harvested at P24, which indicated that 13.07% ± 0.85% of Sun1-GFP labeled nuclei coexpress Sox10 (n = 4; Fig. S2). We therefore cannot completely exclude a small contribution of oligodendrocyte-expressed genes to the observed transcriptional effects of neonatal incision, although we would note that our previous Sox9 immunostaining (Fig. 4c) clearly indicates that astrocytes make up the majority of Sun1-GFP+ cells from which RNA was extracted.
Taken together, these results suggest that neonatal incision acutely modulates the transcriptional profile of dorsal horn astrocytes but does not evoke long-term changes to astrocyte gene expression within the spinal nociceptive network.
Neonatal incision has sexually divergent effects on microglial engulfment of astrocyte matter in the developing dorsal horn
One day after injury, we found an upregulation in the transcription of cytoskeletal genes (Fig. 6) along with increased astrocyte size and complexity (Figs. 2, 3). Furthermore, while we found that neonatal incision reduced astrocyte volume and branching at P10 and P24, it evoked minimal changes to astrocyte gene expression at these time points. This raised the possibility that the incision-mediated decrease in astrocyte complexity observed at these later time points resulted from the pruning of astrocytes by microglia rather than cell-autonomous processes occurring within the astrocyte itself. Notably, P3 hindpaw incision has been shown to increase microglial engulfment of both A-fiber synapses and GABAergic synapses in the developing dorsal horn (Xu et al., 2023). While the microglial pruning of astrocyte processes has not been demonstrated in the spinal cord, microglia have been shown to engulf astrocyte soma and processes in the retina, amygdala, and nucleus accumbens (Puñal et al., 2019; VanRyzin et al., 2019; Testen et al., 2024). We therefore sought to determine whether microglial engulfment of astrocytes occurs in the developing dorsal horn, and if so, the degree to which these interactions might be impacted by neonatal incision.
Aldh1l1-tdTOM mice were injected with full-dose tamoxifen at P1 to fluorescently label all spinal astrocytes, then either subjected to hindpaw incision or anesthesia only at P3. We next harvested the lumbar enlargement of the spinal cord at P8 to capture a time during the transition from an incision-evoked increase in astrocyte size and complexity (at P4) to a decrease in these metrics (beginning at P10). To visualize microglia and astrocyte interactions, we stained thick (100 µm) transverse spinal cord sections with Iba1 to label microglia and CD68 to label phagocytic lysosomes. We found clear evidence of tdTomato-labeled astrocyte matter within Iba1- and CD68-labeled microglial lysosomes (Fig. 9a,b; Movie 3). Using machine learning segmentation in Imaris, we quantified the number and volume of astrocyte inclusions within microglial lysosomes (Fig. 9c). Given the known sexual dimorphism in the regulation of spinal nociceptive processing by microglia (Sorge et al., 2015; Inyang et al., 2019; Moriarty et al., 2019; Agalave et al., 2021; Fan et al., 2025), this analysis was conducted using sex as an independent variable. The data demonstrate that P3 incision has sexually divergent effects on microglial engulfment of astrocytes (n = 7–8 mice per group; injury × sex interaction: p = 0.0007; two-way ANOVA; Fig. 9d). While neonatal incision decreases the volume of engulfed astrocyte matter in males (p-adj = 0.0494; Šidák's multiple-comparisons test), incision increases engulfed astrocyte volume in females at this time point (p-adj = 0.0109). Similarly, incision decreases the total number of astrocyte inclusions within microglial lysosomes in males (injury × sex interaction: p = 0.0003; p-adj = 0.0287; Fig. 9e) but increases the number of inclusions in females (p-adj = 0.0061). We found no differences in the total microglial volume between males and females or between naive and incision groups (Fig. 9f) suggesting that microglia fail to proliferate in the neonatal dorsal horn after early life surgical injury, although a more extensive analysis of microglial morphology, density, and gene expression would be needed to elucidate the degree to which microgliosis occurs in the developing dorsal horn after hindpaw incision. Lysosomal volume was also unchanged by incision and was similar between males and females (Fig. 9g).
Neonatal hindpaw incision increases microglial engulfment of astrocyte matter selectively in the female dorsal horn during early life. a, Representative 3D rendering of Imaris generated surface maps of microglia labeled with Iba1 (green), lysosomes labeled with CD68 (cyan), and engulfed astrocyte matter tdTomato (red). Scale bar, 2 µm. b, Raw images of the same cell in panel a shown in three different planes: x–y (top left), y–z (top right), x–z (bottom left). Microglia are labeled with Iba1 (green), lysosomes labeled with CD68 (magenta), and astrocytes labeled with tdTomato (red). Scale bar, 2 µm. c, x–y slice view of the same cell demonstrating Imaris identification of astrocyte matter engulfed by microglia. Iba 1, green; CD68, magenta; tdTomato, red; Imaris surface outline, cyan. d, P3 incision decreases the volume of engulfed astrocyte matter within microglial lysosomes in males but increases engulfed astrocyte volume in females at P8 [injury × sex interaction: F(1,26) = 14.64, p = 0.0007; two-way ANOVA; male naive vs incision: *p-adj = 0.0494; female naive versus incision: *p-adj = 0.0109, Šidák's multiple-comparisons test; n = 7–8 mice per group (i.e., male naive, male incision, female naive, female incision); each data point represents one mouse, 3–4 images from 2–3 sections averaged per mouse]. e, Neonatal incision decreases the total number of astrocyte inclusions within microglial lysosomes in males but increases the number of inclusions in females (injury × sex interaction: F(1,26) = 17.32, p = 0.0003; two-way ANOVA; *p-adj = 0.0287, **p-adj = 0.0061, Šidák's multiple-comparisons test). f, Total microglial (Iba1) volume is not impacted by P3 incision (injury: F(1,26) = 0.04081, p = 0.8415; sex: F(1,26) = 0.04828, p = 0.8278; injury × sex interaction: F(1,26) = 2.787, p = 0.1070, two-way ANOVA). g, Microglial lysosome (CD68) volume in the P8 dorsal horn is unchanged by P3 incision (injury: F(1,26) = 0.2164, p = 0.6456; sex: F(1,26) = 0.002678, p = 0.9591; injury × sex interaction: F(1,26) = 0.01824, p = 0.8936, two-way ANOVA).
High-resolution imaging of individual dorsal horn astrocytes. Representative video of a tdTomato-labeled spinal dorsal horn astrocyte visualized in 3D space. Image generated with Yokogawa SoRa spinning disk confocal microscope, z step = 0.2 µm. [View online]
Three-dimensional quantitative morphological analysis of spinal astrocytes. Video demonstrating quantification of astrocyte morphology using Bitplane Imaris software (7–14 s, Imaris surfaces rendering of astrocyte; 15–22 s, Imaris filaments rendering of astrocyte; 22–23 s, convex hull generated using Imaris filaments model). [View online]
Microglial engulfment of astrocytes. Video showing 3D rendering of representative image in Figure 9a. Green translucent object corresponds to Imaris generated surface model of Iba1-labeled microglia. Red solid object corresponds to Imaris AI-based segmentation generated surface model of astrocyte matter within microglial lysosomes. [View online]
These data demonstrate that microglial engulfment of astrocyte material does indeed occur in the developing spinal dorsal horn and that neonatal incision elicits differential effects on the engulfment of astrocyte material in male and female mice.
Discussion
The neonatal dorsal horn is characterized by marked hyperexcitability to sensory stimuli and rapidly increasing excitatory tone (Fitzgerald, 1985; Baccei et al., 2003). The demonstration that SDH astrocytes undergo profound age-dependent changes during the first three postnatal weeks raises the possibility that the maturation of spinal astrocytes shapes the physiology of the developing SDH circuit. We found that spinal astrocytes become significantly larger and more complex between P4 and P24 (Figs. 2, 3). This is consistent with prior reports characterizing the maturation of astrocytes in the brain showing that distinctive morphological features of astrocytes, such as the highly ramified and fine branchlets and the tiling of astrocytes into distinct territories, do not fully develop until 30 d after birth (Bushong et al., 2004; Felix et al., 2021). The proximity of astrocyte leaflets to the synapse is critical to their ability to buffer ions and neurotransmitters (Ghézali et al., 2016; Toman et al., 2023). Additionally, the perisynaptic astrocyte compartment forms a sort of “cradle” around the synaptic cleft, creating an important structural and functional barrier that confers spatial specificity to synaptic transmission (Ramon y Cajal and Azoulay, 1955; Verkhratsky and Nedergaard, 2014). It may be that the nascent architecture of astrocytes during early life yields increased glutamate spillover, thus contributing to the increased excitability of the neonatal SDH.
We also found extensive age-dependent changes in the expression of genes coding for transporters known to mediate astrocytic buffering of ions and neurotransmitters (Fig. 5), which may also contribute to increased excitability in the immature dorsal horn. Transcription of Slc1a2 (encoding GLT-1) and Kcnj10 (Kir4.1) significantly increases with age in spinal astrocytes. This suggests a potential developmental delay in the ability of astrocytes to efficiently uptake glutamate and K+ (Rothstein et al., 1996; Kofuji and Newman, 2004; Neusch et al., 2006; Lauriat and McInnes, 2007). In contrast, the transcription of Slc1a3 (GLAST) was found to decrease between P4 and P24. This is consistent with prior characterization of GLAST in the spinal cord showing that Slc1a3 transcripts increase between E11 and E18 and then decrease between E18 and adulthood (Shibata et al., 1997). While GLT-1 is present in most regions, GLAST is restricted to the cerebellum and superficial laminae of the dorsal horn (Rothstein et al., 1994; Torp et al., 1994; Lehre et al., 1995), supporting astrocyte heterogeneity across the CNS (Zhang and Barres, 2010; Kronschläger et al., 2021; Torres-Ceja and Olsen, 2022).
The heterogeneity of spinal astrocytes has been directly demonstrated by recent single-cell and single-nucleus transcriptomic analyses from mice, humans, and macaques which have highlighted the existence of multiple subpopulations (Rosenberg et al., 2018; Russ et al., 2021; Arokiaraj et al., 2024; Chen et al., 2025) which may play distinct roles in nociceptive processing. By pooling the isolated astrocytic nuclei prior to bulk RNA-seq analysis, the present study may have missed transcriptional changes related to age, or induced by neonatal incision, that are restricted to a specific subset of dorsal horn astrocytes. Furthermore, the sequencing of nuclear RNA may have resulted in a failure to detect injury-evoked changes in the transcription of other key astrocytic genes, including those related to metabolism and neuroinflammation, that might emerge following the sequencing of mature mRNA. Neonatal incision also likely modifies the function of developing astrocytes via a host of post-transcriptional mechanisms that would be missed with the chosen experimental approach. Finally, the conclusions drawn from this study are restricted to the first three weeks of postnatal development. Although prior work in the brain shows that much of astrocyte development occurs during the first month after birth (Felix et al., 2021), the extent to which neonatal incision impacts dorsal horn astrocytes in adulthood is not yet known.
Neonatal injury has both short-term and long-term effects on synaptic transmission in the spinal dorsal horn (Li and Baccei, 2011, 2016; Li et al., 2013, 2015). The gene encoding thrombospondin-1 (Thbs1) was strongly upregulated at P4 following neonatal incision (Fig. 6) which was interesting given that thrombospondin-1 (TSP1) drives the formation of glutamatergic synapses via binding to receptors that include the α2δ1 calcium channel subunit (Christopherson et al., 2005; Eroglu et al., 2009; Xu et al., 2010) and P3 incision increases the frequency of miniature excitatory postsynaptic currents in SDH interneurons 1–2 d postinjury (Li and Baccei, 2011). Although the upregulation of Thbs1 resolved by P10, transient increases in Thbs1 expression may contribute to the persistent changes to synaptic transmission after neonatal incision by increasing the overall number of excitatory synapses in the SDH, thus disrupting the balance of synaptic excitation and inhibition. Notably, TSP1 and TSP2 are downregulated during postnatal development across the CNS, suggesting that these factors are not required for the maintenance of glutamatergic synapses during adulthood (Wang et al., 2012). It would be interesting to determine if the administration of gabapentinoids, a class of analgesics which binds to the α2δ1 subunit (Gee et al., 1996; Field et al., 2006), during early life can prevent the priming of developing pain pathways by neonatal tissue damage. Although future studies testing the causal role of astrocytic Thbs1 will be necessary, our results suggest a potential role for astrocytes in the persistent rewiring of synaptic circuits in the spinal dorsal horn after neonatal tissue injury.
Another potential explanation for altered synaptic function in SDH circuits after neonatal injury might be the resulting changes to developing astrocyte morphology (Figs. 2, 3). Given that the highly ramified and elaborate structure of astrocytes is crucial to their many functions (Baldwin et al., 2024), injury-evoked changes to astrocyte morphology may have significant consequences for synaptic function within the SDH. Decreased astrocyte size and complexity following early life injury may suggest a compromised ability for the astrocyte to communicate with the synapse via secreted factors (Beattie et al., 2002; Stellwagen et al., 2005; Crawford et al., 2012), regulate neuronal excitability by buffering ions and neurotransmitters in the synaptic cleft (Danbolt, 2001; Kofuji and Newman, 2004), and meet the metabolic demands of the neuron (Pellerin et al., 2007; Attwell et al., 2010; Suzuki et al., 2011). Although the microscopy techniques employed in this study allowed us to resolve structures down to 200 nm, astrocyte branchlets and leaflets are 10–100 nm in diameter, which cannot be resolved using light-based microscopy. Therefore, future studies using serial block face scanning electron microscopy (Aten et al., 2022) will be necessary to identify changes to astrocyte nanostructure.
Although the effects of astrocyte-secreted factors on synaptic function are well characterized, the influence of neuronal activity on astrocyte function is poorly understood by comparison (Hasel et al., 2021). Since the developing dorsal horn is highly regulated by the activity of primary afferent inputs (Koch and Fitzgerald, 2013), and blocking peripheral nerve activity with bupivacaine prior to hindpaw incision is sufficient to prevent neonatal priming (Moriarty et al., 2018), a better understanding of how dorsal horn astrocytes respond to nociceptor activity is of clear interest. Analysis of potential interactions between the ligands expressed by neonatal sensory neurons and the receptors expressed by developing spinal astrocytes (Tables S1–S5) can yield novel hypotheses regarding the identity of the astrocytic receptors that sense changes in neuronal activity during early life. These receptors, and the downstream intracellular cascades triggered by their activation, represent important targets for future studies.
Lastly, we found that neonatal incision changes the microglial engulfment of dorsal horn astrocytes in a sex-dependent manner (Fig. 9). While P3 injury increased microglial engulfment of astrocytes in females, the opposite effect occurred in males. Although neonatal priming occurs in both males and females, the neuroimmune mechanisms by which priming occurs are different between sexes (Moriarty et al., 2019; Dourson et al., 2023). Interestingly, while microglial phagocytosis of GABAergic synapses in the SDH after early life injury occurred independently of sex (Xu et al., 2023), other evidence suggests that the microglial pruning of synapses can depend on sex (Lenz and McCarthy, 2015), highlighting the need to explore potential sexual dimorphism of microglial phagocytosis in the CNS. Sex differences in microglial phagocytosis of astrocyte progenitors have been shown (VanRyzin et al., 2019); thus, it is plausible that the microglial pruning of astrocytic processes is similarly sex specific. Whether the observed changes to astrocyte complexity following neonatal incision are a direct consequence of altered microglial phagocytosis of astrocytes is not yet clear. Astrocytic incision-evoked DEGs coding for canonical “eat me” and “don't eat me” immune signals were not detected at any developmental time points. Future studies to determine the degree to which clopidogrel, whose bioactive metabolite is a membrane-permeant P2Y12 receptor antagonist that disrupts microglial phagocytosis in the brain (Bollinger et al., 2023) and reduces activation of spinal microglia (Yu et al., 2019), prevents the structural remodeling of SDH astrocytes after P3 incision could address this issue.
In sum, our findings represent a critical first step toward elucidating the role of astrocytes in the developmental plasticity of the dorsal horn. By demonstrating that the morphology and transcriptional profile of developing spinal astrocytes are altered by early life injury, we identify astrocytes as potential mediators of injury-evoked changes to neuronal function in the dorsal horn.
Footnotes
This work was supported by the National Institutes of Health (R01NS080889 to M.L.B., F31NS135754 to J.J.Y.). We thank the Cincinnati Children’s Bio-Imaging and Analysis Facility (RRID: SCR_022628; Marina George, PhD; Sarah McLeod; and Amanda Rainey) for their assistance and expertise in image acquisition and analysis. We also thank the Research Flow Cytometry Facility in the Division of Rheumatology at Cincinnati Children's Hospital Medical Center for their assistance.
↵*J.J.Y. and E.K.S. contributed equally to this work.
The authors declare no competing financial interests.
This paper contains supplemental material available at: https://doi.org/10.1523/JNEUROSCI.1197-25.2025
- Correspondence should be addressed to Mark L. Baccei at mark.baccei{at}uc.edu.
