Single-nucleus and spatial signatures of the brainstem in normal brain and mild traumatic brain injury in male mice

Zhuang, Yuan; Liao, Xixian; Niu, Fei; Li, Ming; Yan, Yu; He, Chuanhang; Wu, Xiang; Tian, Runfa; Gao, Guoyi

doi:10.1038/s41467-025-59856-2

Download PDF

Article
Open access
Published: 31 May 2025

Single-nucleus and spatial signatures of the brainstem in normal brain and mild traumatic brain injury in male mice

Nature Communications volume 16, Article number: 5082 (2025) Cite this article

2551 Accesses
Metrics details

Subjects

Abstract

The mammalian brainstem is particularly vulnerable to mild traumatic brain injury (mTBI), which is associated with prolonged autonomic dysfunction and coma. The spatial cellular connections within the brainstem or the mechanisms underlying its response to injury have been underestimated. Here, we performed single-nucleus RNA sequencing with spatial transcriptome sequencing in both normal and mTBI brainstems in male mice, revealing thirty-five neuron and non-neuron clusters. Typically, we identified subtypes of neurons that co-release multiple neurotransmitters, especially in the sagittal midline of the brainstem. Spatially adjacent neurons sharing similar gene expression patterns. The brainstem’s response to mTBI has two features: (1) Oligodendrocytes around the fourth ventricle exhibit widespread disconnection at 1-h post-injury, and (2) Injury-related noradrenergic neurons, particularly in their interaction with neurons located in theIRt and the Sol. These findings provides a reference for further integrative investigations of cellular and circuit functions of brainstem.

Single-nucleus profiling of adult mice sub-ventricular zone after blast-related traumatic brain injury

Article Open access 05 January 2023

Brain-wide reconstruction of inhibitory circuits after traumatic brain injury

Article Open access 14 June 2022

Single-nucleus transcriptomic mapping of blast-induced traumatic brain injury in mice hippocampus

Article Open access 20 September 2023

Introduction

The mammalian brainstem represents a vital hub connecting the telencephalon and spinal cord, playing a crucial role as an interface involved in action planning¹. It also regulates the autonomic nervous system and integrates crucial physiological processes, such as respiratory regulation^2,3, cardiovascular control^4,5 and receiving ascending sensory inputs. These sensory inputs are integrated into the reticular activating system, which participates in arousal and wakefulness^6,7, attention and focus, and consciousness^8,9. Traditional anatomical studies have broadly characterized the brainstem as comprising both grey and white matter. The grey matter is involved in regulating reflexes, arousal, autonomic functions, and maintaining homoeostasis, while the white matter contains ascending sensory pathways, as well as descending pyramidal and extrapyramidal tracts responsible for motor control¹⁰. Further research suggests that the heterogeneous cell populations and the dynamic interactions between brainstem nuclei require deeper exploration to fully understand their complexity.

The primary concern regarding brainstem dysfunction arises from traumatic brain injury. Compared with severe traumatic brain injury (sTBI), which often leads to nearly irreversible neurological dysfunction and high mortality¹¹, mild traumatic brain injury (mTBI), often caused by rotational acceleration and shearing forces, is more common and can result in persistent symptoms that significantly impact quality of life^12,13. The brainstem, which is distinct from other cerebral regions, is more vulnerable to mTBI because its primary axis intersects with the rostro-caudal axis of the brain. Previous studies have highlighted that the characteristic symptoms observed after mTBI predominantly stem from disruptions to the autonomic nervous system¹⁴. Furthermore, extensive imaging research has identified brainstem injury as a central factor mediating these impairments¹⁵. However, owing to its intricate inter-structural connections and meticulous arrangement of nuclei, studying the brainstem, particularly following mTBI, remains exceedingly challenging.

Owing to advancements in single-cell (single-nucleus) and spatial transcriptome technologies, atlas studies of the mammalian cortex have been reported, revealing cortical functional heterogeneity and regional functional connections among nuclei^16,17. However, to our knowledge, no specialized single-cell or spatial transcriptome studies have targeted the entire brainstem. Therefore, we applied single-nucleus and spatial transcriptome methods to identify unique marker genes specifying each neuron and glial cell subtype and to investigate the distribution of these cell types. Furthermore, we explored the interactions between neurons and gene expression models to reveal neuronal connection preferences and gene expression coordination among specific cell-type identities. Next, we introduced an mTBI model and detected changes in glia and noradrenergic neurons in the brainstem post-injury, addressing questions about how the brainstem responds to traumatic injury.

Results

snRNA-seq profiling of cell types in the brainstem

To characterize cellular heterogeneity in the brainstem, we performed snRNA-seq using cells dissociated from mouse brainstems in normal, 1 h post-injury, and 24 h post-injury states (Fig. 1a). Ultimately, we generated 29,035 nuclei after quality filtering and detected 2216 genes per nucleus on average. We combined the cells from the control and injury groups to increase the number of cells and detected low-abundance cell populations for clustering analysis. On the basis of the unsupervised method, we identified 35 clusters (Fig. 1b), and with the expression of the pan-neuronal marker Snap25¹⁸, the 35 clusters were further divided into 18 neuronal and 17 non-neuronal clusters. Among the non-neuronal clusters, 10 specific cell types, such as astrocytes, microglia, oligodendrocytes, OPC, ependymal cells, endothelial cells, fibroblasts, pericytes, and the choroid plexus were recognized by their unique markers: Aldh1a1, Ctss, Mbp, Vcan, Ccdc170, Cldn5, Col8a1, Pdgfrb, and Ttr^19,20,21 (Fig. 1c, d).

**Fig. 1: Single-nucleus RNA-sequence transcriptional profile of the whole brainstem.**

We reclassified the neuronal clusters according to neurotransmitter-synthesizing enzymes or transporter genes. Neurons can be further divided into seven main types, such as Glutamatergic neurons (Slc17a6, Slc17a7, Slc17a8), GABAergic neurons (Gad1, Gad2), Glycinergic neurons (Slc6a5, Slc7a8), Cholinergic neurons (Chat, Ache), Serotonergic neurons (Fev, Slc6a4), Dopaminergic neurons (Th, Slc18a2), and Noradrenergic neurons (Slc16a2, Slc6a2). As expected, the three most abundant neuron types were Glutamatergic neurons (7919, comprising 49%), GABAergic neurons (5081, comprising 32%), and Glycinergic neurons (2571, comprising 16%) (Fig. 1e).

To uncover additional neuronal signatures, we further subdivided the main neuron types. Through unsupervised analysis and manual examination, we verified neurons co-expressing two or more unique markers characteristic of other neuron types. However, no special clusters were found typically in the injury group (Fig. S1). Among the Glutamatergic clusters, 12 sub-clusters were further distinguished (Fig. 2a). Most sub-clusters expressed Slc17a6, except for sub-cluster 5, which predominantly expressed Slc17a7. However, Slc17a8 was not abundantly expressed in any of the sub-clusters (Fig. 2b). Furthermore, among the sub-clusters, we identified four groups containing markers characteristic of Glycinergic neurons (Glu-Gly), four groups expressing only glutamatergic markers (glutamatergic neurons), one group expressing the typical GABAergic marker (Slc6a1, Glu-GABA), and one group co-expressing noradrenaline and the typical GABAergic marker (Reln) (Glu-Net-GABA (Reln)). Additionally, we observed that one group, which co-released glycine, also expressed the typical GABAergic marker (Pvalb) (Glu-Gly-GABA (Pvalb)). Furthermore, only one group co-released GABA, noradrenaline, and acetylcholine (Glu-GABA-Net-choline) (Fig. 2a). The three largest sub-populations within the Glutamatergic clusters were the glutamatergic, Glu-GABA, and Glu-Gly sub-clusters (Fig. 2c). Unlike the classification in the cortex, in the GABAergic clusters, we did not detect neurons that significantly expressed Sst and Vip^22,23 (Fig. 2e). By contrast, neurons in the 10 groups expressed traditional GABAergic markers, such as Gad1 or Gad2 (GABAergic neurons), with three of them co-expressing glycinergic markers (GABA-Gly), one group co-expressing glycinergic and glutamatergic markers (GABA-Gly-Glu), and one group co-expressing glycinergic and cholinergic markers (GABA-Gly-choline). Additionally, we identified three groups of GABAergic neurons expressing Reln, with two of which co-releasing glutamate (GABA(Reln)-Glu) and one of which co-released glycine (GABA (Reln)-Gly). Furthermore, we recognized only one group of GABAergic neurons expressing Pvalb, which co-released glycine (GABA (Pvalb)-Gly). Among the sub-clusters, the ‘GABAergic neurons’ cluster contained the most GABAergic neurons (Fig. 2d, f). The Glycinergic neurons contained nine sub-clusters, among which four groups of neurons expressed only glycine (glycinergic neurons) and contained the greatest number of glycinergic neurons, and three groups of neurons co-released GABA (Gly-GABA). Additionally, two groups co-expressed the glutamatergic marker (Slc17a6) (Gly-Glu), and one group co-expressed GABA and glutamate (Gly-Glu-GABA) (Fig. 2g, h). Interestingly, although we detected only 95 noradrenergic neurons, 73% of these cells co-expressed glutamate (Net-Glu), whereas the remaining 27% highly expressed only typical noradrenergic markers (noradrenergic neurons) (Fig. 2i).

**Fig. 2: Subdivision of neurons co-expressing markers characteristic of other neuron types.**

Spatial distribution and gene expression pattern heterogeneity revealed the functional diversity of neurons in the normal brainstem

The brainstem, especially the reticular activating system, executes complex functions related to vital signs. However, owing to its vague spatial distribution and plentiful neuronal connections, illustrating neuronal heterogeneity at the spatial level remains challenging. On the basis of previous studies and the Allen Mouse Brain Atlas, we confirmed that the core area of the reticular activating system is located near the midline; we then conducted spatial transcriptomics analysis using sagittal brainstem tissue located 830 μm lateral to the midline (Fig. 3a). We located the neurons spatially using robust cell type decomposition (RCTD). Compared with Glutamatergic neurons, GABAergic neurons were much more widely distributed in the mesencephalon, pons, and medulla oblongata, especially GABAergic neurons expressing Reln and co-releasing glutamate (Fig. 3b). Additionally, we identified several typical regions with relatively clear boundaries and clusters of specific neuron types. For example, there was a clustering of Dopaminergic neurons in the ventral tegmental area (VTA) of the mesencephalon (Figs. 3c–e, and S2a, b), whereas the neurons expressing only noradrenaline were clustered at the locus coeruleus (Figs. 3f–h and S2c, d). Cholinergic neurons were clustered in the facial nerve nucleus and olivary nucleus in the ventral medulla oblongata (Figs. 3i–k and S2e, f), whereas Glu-Net-GABA (Reln), Gly-Glu, Gly-GABA, and Gly-Glu-GABA neurons are co-localized along the spinal trigeminal nucleus. (Fig. S4g–i) The subtypes of Glutamatergic neurons, GABAergic neurons, and Glycinergic neurons were scattered throughout the middle of the sagittal brainstem (Fig. S4a–f), and GABA-Gly-choline, GABA-Gly-Glu, Glu-GABA-Net-choline, and Gly-GABA (Reln) neurons were distributed mainly in reticular structures, such as the pontine reticular nucleus (PnO) and intermediate reticular nucleus (IRt) (Fig. 3l).

**Fig. 3: Spatial distribution of the neuron clusters.**

By performing high-dimensional weighted gene co-expression network analysis (hdWGCNA), we identified eight neuronal functional gene expression patterns (Fig. 4a–c). Each similar pattern represented the synergistic effects of biological functions. We found that the patterns of expression of the hub genes were distinct across different neuron groups, even within the subtypes of a particular neuron type. Therefore, we further scored the hub genes at the spatial level and found that spatially neighbouring neuronal clusters presented similar gene patterns. For example, Noradrenergic, Serotonergic, and Dopaminergic neurons were located around the fourth ventricle, dorsal pons, and medulla oblongata. These neurons up-regulated hub genes in Model 3 (M3) and slightly in Model 5 (M5), which was consistent with the adjacent neuronal groups scattered in reticular structures, such as GABA-Gly-choline, Glu-GABA-Net-choline, and Net-Glu neurons, and were related to GABAergic synapse and retrograde endocannabinoid signalling pathways (Fig. 4d). Moreover, we found that neurons co-located in the spinal trigeminal nucleus, such as Gly-Glu and Gly-Glu-GABA, were down-regulated in Model 2 (M2) (Fig. 4e). The GABA (Reln)-Gly group, which was primarily distributed in the cerebellum, down-regulated genes in Models 2 (M2), 5 (M5), 6 (M6), and 8 (M8), and slightly up-regulated genes in Model 1 (M1). Nevertheless, the widespread GABA (Reln)-Glu neurons highly up-regulated genes in Model 1 (M1) and Model 7 (M7). These genes were enriched in cell communication and neurotransmitter–receptor activity processes, neuroactive ligand–receptor interactions, and cell adhesion molecule pathways. Conversely, these neurons down-regulated genes in Models 2 (M2), 3 (M3), 4 (M4), and 6 (M6) (Fig. 4f).

**Fig. 4: Neuron clusters exhibited spatial co-expression patterns.**

Oligodendrocytes exhibited widespread disconnection from neurons 1 h post-traumatic brain injury

Acute brain injury results in widespread neurological dysfunction which influences the brainstem. Using the drop-weight impact model, we found that mice had acutely decreased respiratory and heart rates; and moreover, went into a coma. By calculating the cell ratios among the three groups, we detected slight changes in the number of astrocytes and microglia, indicating a limited inflammatory response. However, the number of oligodendrocytes was significantly lower in the injury–1 h group than in the control group (20.1% vs. 23.8%) and the injury–24 h group (20.1% vs. 24.3%) while the number of OPCs remained unchanged. (Figs. 5a and S3). To explore cell heterogeneity in the different injury groups, we re-clustered the 6763 oligodendrocytes without eliminating batch effects and identified nine distinct cell populations. Ultimately, on the basis of previous studies, nine clusters were distinguished according to the hallmarks of typical oligodendrocyte states (Fig. 5b). We identified one differentiation-committed oligodendrocyte precursor (COP) cluster co-expressing Pdgfra and Cspg4, two newly formed oligodendrocyte (NFOL) clusters expressing genes related to early stages of differentiation, such as Tcf7l2 and Sox10; one myelin-forming oligodendrocyte (MFOL) cluster highly expressing genes responsible for myelin formation, such as Mog, Plp1, and Opalin; and five mature oligodendrocyte (MOL) clusters expressing genes involved in late oligodendrocyte differentiation, such as Il33, Apod, and Ptgds^24,25 (Fig. 5b, c).

**Fig. 5: Oligodendrocytes are more vulnerable to acute mild traumatic brain injury.**

By performing UMAP, we found a distinct distribution among groups at the two-dimensional level, and the injury–1 h group had distinct features from the control and injury–24 h group. We found that MOL-4 clusters were predominantly included in the injury–1 h group (Fig. 5b, d). To investigate the effect of traumatic injury effect on the oligodendrocytes, we further evaluated the six expression patterns of genes across the groups using hdWGCNA (Fig. 5e–g) and noted that the gene expression patterns in the injury–1 h group were quite distinct from those in the control group and the injury–24 h group. In brief, genes whose expression was down-regulated only in the control group in Model 1 (M1), were enriched in the endocytosis, tight junction, autophagy, and AMPK signalling pathways. These genes were up-regulated in the injury–1 h group but down-regulated in the injury–24 h group. Conversely, the genes in Model 5 (M5), which were enriched in axon guidance, regulation of the actin cytoskeleton, GABAergic synapse, and glutamatergic synapse pathways, were down-regulated in the injury–1 h group but up-regulated in the injury–24 h group (Fig. 5h). Using Slingshot (2.12.0) and Monocle (2.30.1), we first confirmed that oligodendrocytes underwent significant changes 1 h post-injury by fitting group-related pseudo-time trajectory (Fig. 6a). We then modelled gene expression trends and identified two pseudo-time lineages, mapping the differentiation trajectory of oligodendrocytes. This analysis revealed a pathway from COPs to NFOLs and MFOLs, which then diversified into MOLs (MOL-4 and MOL-5) (Fig. 6b). We then performed a pseudotemporal ordering of gene dynamics across the groups according to the post-injury time and verified that the injury–1 h group presented distinct features at the endpoint of the pseudo-time transition (Fig. 6c). We also confirmed that MOL-4 was characteristic of the injury–1 h group by conducting a cell proportion analysis.

**Fig. 6: Oligodendrocytes located around the fourth ventricle exhibit widespread disconnection from neurons following acute mild traumatic injury.**

On the basis of the above evidence, we retrieved the pseudo-time-rated genes from COPs to MOL-4 and investigated their putative functions using KEGG analysis. The up-regulated genes were involved in endocytosis and the PI3K-Akt and AMPK signalling pathways, representing injury-induced early degeneration processes. By contrast, the down-regulated enrichment processes focused on glutamatergic, dopaminergic, and serotonergic synapses, adrenergic signalling in cardiomyocytes, and GABAergic synapses, indicating widespread disconnection of oligodendrocytes from neurons 1 h post-injury (Fig. 6d). By sub-setting the MOL-4 cluster and analysing the most down-regulated genes (Fig. 6e), we focused on the enriched pathways, such as extracellular matrix structural constituent and signalling receptor regulator activity (Fig. 6f). Interestingly, after further projecting the single-cell results of oligodendrocytes onto the region of our interest in the spatial slice, we found that MOL-4-related genes were predominantly expressed in the dorsal pons approaching the fourth ventricle and were spatially adjacent to neurons, including GABA(Pvalb)/Gly, GABA/Gly, Glu/GABA, Gly/GABA(Reln), Glu/Gly, Glu/Gly/GABA(Pvalb), glycinergic neurons, noradrenergic neurons, Net/Glu, and Serotonergic neurons (Fig. 6g). After considering short-range interactions, such as paracrine signalling, and somal or synaptic contact, we observed significantly attenuated connections and reduced cell adhesion, glutamate, or GABA signalling between MOL-4 oligodendrocytes and nearby neurons in the injury–1 h group, which aligned with the pseudo-time results (Fig. 6h).

Changes in the connections of noradrenergic neurons in the locus coeruleus are characteristic responses of the brainstem to traumatic brain injury

The nuclei of the brainstem are complex and functionally diverse, especially the dorsal pons, which are crucial for the autonomic nervous system and act as an integrating hub for diversified autonomic functions. Notably, the locus coeruleus plays a pivotal role in various functions, including stress responses, emotional regulation, attention, and arousal^6,26. Unfortunately, we still lack a comprehensive understanding of the changes in locus coeruleus-associated neurons following injury. Here, we primarily assessed locus coeruleus function, such as the tail-flick response to cold stimuli, by observing mouse behaviour. We found that the delay in the tail-flick response was prolonged 1 h post-injury but shortened 24 h post-injury, indicating a process of functional recovery (Fig. 7a). Notably, we performed the gene-set ranking using neuronal activity markers and reported that the injury caused a widespread reduction in neuronal activity, especially in the telencephalon, thalamus, and mesencephalon. However, we detected regions, such as the locus coeruleus and the intermediate reticular nucleus (IRt) with increased neuronal activity in the pons and medulla oblongata (Fig. 7b). To gain deeper insight into the functional changes in the locus coeruleus after injury, we extracted its barcode information and performed transcription factor (TF) and enrichment analyses. We found that 33 TFs were activated in the injury-24h group through TF certification and target gene identification. Essrg, Creb5, Atf2, Bmyc, and Trp53 were involved in up-regulating the noradrenaline pathway (Fig. 7c). Through unsupervised enrichment analysis, we identified injury-induced up-regulation in neuropeptide signalling and molecular transducer activity processes. Additionally, we observed an increase in noradrenaline-related blood pressure functions (Fig. 7d).

**Fig. 7: Changes in the locus coeruleus are a feature of the brainstem response to mild traumatic brain injury.**

snRNA-seq integrated with spatial transcriptome analysis revealed that most typical noradrenergic neurons were located within the locus coeruleus, whereas Net-Glu neurons were scattered around it. This distribution coincided with the ___location of MOL-4 oligodendrocytes, an area susceptible to injury. We then evaluated the expression of noradrenaline-synthesizing enzymes to assess the response of noradrenergic neurons to injury over time. We found that the pivotal enzymes Dbh and Net (Slc6a2) were gradually upregulated after the major traumatic brain injury (Figs. 7e, 8a). Moreover, along the pseudo-time curve, processes, such as neurotransmitter–receptor activity, serotonin binding, ionotropic glutamate receptor activity, and neurotransmitter–receptor activity involved in regulating postsynaptic membrane potential were up-regulated. This finding indicated that traumatic brain injury caused enhanced activity of noradrenergic neurons (Fig. 7f).

**Fig. 8: Changes in the connectivity between the locus coeruleus and multiple brain regions are a key characteristic of the brainstem’s response to mTBI.**

Next, we employed CellChat analysis with the ‘noradrenaline’ pathway to demonstrate the interactions between noradrenergic neurons and other neuron types, aiming to detect their contact preferences in typical and traumatic brain injury situations. Noradrenergic neurons usually connect with subtypes of Glutamatergic neurons, such as Glu-GABA-Net-choline, and glutamatergic neurons; subtypes of GABAergic neurons, such as GABAergic neurons and GABA-Gly-Glu; and Net-Glu, Gly-GABA (Reln), and small cluster populations, such as Cholinergic, Serotonergic, and Dopaminergic neurons (Fig. 8b). Spatially, these connected neurons are located in the sagittal middle of the brainstem, such as the IRt, the oral part of the PnO, the medullary dorsal vagal nucleus (MdV), and the nucleus of the solitary tract (Sol) according to the mouse atlas and verified by virus tracing, which is related to the autonomic nervous system function. We also detected the connected neurons in the thalamus and cerebral cortex, confirming that the noradrenergic neurons are involved in advanced neural regulation (Fig. 8c, d). One hour post-injury, the ‘noradrenaline’ interaction pathway was sensitive to the injury, and the interactions of noradrenergic neurons decreased especially the Glu-GABA-Net-choline neurons. However, after 24 h, the interactions of noradrenergic neurons were not only recovered but were also significantly enhanced, especially those of the Glu-GABA-Net-choline, glutamatergic, and Gly-GABA (Reln) neurons located in the IRt and Sol (Fig. 8c). Importantly, the changes in noradrenergic transmission could characterize the response of the brainstem to traumatic injury.

Discussion

The brainstem dominates numerous fundamental autonomic functions closely tied to vital signs, including arousal, sensation, and cardiorespiratory regulation. Its dense concentration of ascending and descending neural pathways, coupled with intricate interactions among nuclei, presents challenges for comprehensive research¹⁴. Previous case reports and imaging research have verified that the brainstem is vulnerable to TBI^27,28; however, a limited understanding exacerbates the difficulty of rescue and rehabilitation. In the present study, we performed snRNA-seq and ST-RNA-seq on brainstems under normal and injury conditions. By identifying genes involved in neurotransmitter synthesis and secretion, we identified seven main neuronal types. Interestingly, most neurons express multiple types of neuronal markers, contributing to the diversity of brainstem neurons and their complex functions. At the spatial level, adjacent neurons share gene expression patterns, which may facilitate the integration of biological information. We also found that mild traumatic injury had comprehensive effects on oligodendrocytes, especially during the acute stage, leading to widespread disconnection with neurons. By contrast, the injury contributed to the functional enhancement of the locus coeruleus, particularly stimulating the secretion of noradrenaline. Despite transient weakening of neuronal connections during the first hour post-injury, the noradrenergic neurons experienced a rebound with elevated connections 24 h after the injury.

Unlike previous reports that focused on the cortex^29,30,31, we found numerous Glycinergic neurons in the brainstem. Notably, these Glycinergic neurons are located primarily in the medulla oblongata at the tissue level of interest, especially in the spinal trigeminal nucleus (Sp5). Another notable finding is that our research highlights the diverse neuronal clusters among Glutamatergic, GABAergic, and Glycinergic neurons, in which subtypes co-express 2–4 neurotransmitter patterns (e.g., glutamatergic neurons that co-release GABA aminobutyric acid, noradrenaline, and acetylcholine; or GABAergic neurons that co-release glutamate and glycine). At the spatial level, the co-expression patterns appeared along the sagittal midline of the brainstem, particularly in the PnO, IRt, and MdV, which are involved in reticular formation. Moreover, spatially neighbouring neuronal clusters shared similar gene patterns. These results suggest the diversity and complexity of intercellular communication within the brainstem, particularly in the reticular formation. Interestingly, we also observed mTBI-related changes in the cell ratios of inhibitory neurons in the brainstem, while the number of excitatory neurons remained relatively constant. Numerous studies have illustrated that the projection of inhibitory neurons can be remodelled in regions such as the cortex and hippocampus^32,33,34. However, research in the brainstem remains limited. Our findings of a consistently increased ratio of glycinergic neurons and a decreased number of GABAergic neurons in the brainstem may play a role in the brainstem’s acute response to injury, warranting further investigation. (Figs. 5a and S5).

The recent spatial transcriptomics and single-cell studies of the cortex after mTBI have proposed that glial cells undergo a conversion to reactive states, particularly at the injury site. The emergence of reactive glial clusters has been shown to exhibit distinct patterns of reactivity, which likely play a crucial role in the response to injury and recovery³¹. Oligodendrocytes were the most abundant glial cells in our scRNA-seq study and were vulnerable to mTBI. Yoo et al.³⁵ have reported that on MRI imaging, the myelin volume fraction was significantly reduced in the brainstem of patients with mTBI, which could explain the post-concussion syndrome. However, the features that remain on imaging require further molecular verification. In our present study, we detected a marked decrease in oligodendrocytes within the injury-1 h group accompanied by opposing gene expression patterns. Among the oligodendrocyte clusters, we identified one mature oligodendrocyte group (MOL-4), notably enriched in the injury-1 h group. These types of cells exhibit widespread disconnection with numerous neuronal types. This heterogeneous group was predominantly located around the fourth ventricle, indicating a regional preference in the oligodendrocyte response to acute mTBI. Together, these findings suggest potential imaging correlates.

The locus coeruleus is a small segment adjacent to the fourth ventricle, that integrates the ascending and descending neural pathways and plays a vital role in regulating autonomic functions^36,37. We detected robust neuronal activity in the locus coeruleus, with vigorous noradrenaline secretion in response to mTBI 24 h post-injury. The cell interaction prediction also revealed enhanced communication among the noradrenergic neurons in the locus coeruleus. Interestingly, among the interactions, the noradrenergic neurons and the Glu-GABA-Net-choline interaction, which represent the communication between the locus coeruleus and the Sol, were particularly sensitive to injury. Vagus nerve stimulation (VNS) has emerged as a method to induce neuronal plasticity in the CNS with minimal side effects, especially those contributing to plasticity in the locus coeruleus^38,39. Previous studies have demonstrated indirect innervation from the vagus nerve through the Sol to the locus coeruleus. However, owing to the complex relationships between nuclei in the brainstem, the specific pathway involved is poorly understood⁴⁰. This study proposes direct communication between the locus coeruleus and the Sol, providing evidence to support further VNS treatments.

This study utilized snRNA-seq combined with ST RNA-seq to explore the characteristics of the brainstem under both normal and injury conditions. However, additional neuronal subtypes were not fully represented because of the limited sample size. We focused on one sagittal aspect of the tissue for our spatial analysis and plan to investigate additional aspects in future studies. Our current focus was on describing the brainstem features under different conditions. Future studies involving more extensive functional verification of the ascending and descending pathways are warranted.

In brief, the snRNA-seq and ST RNA-seq analyses revealed diverse types of neurons in the brainstem, characterized by complicated interactions and regional similarities in gene expression patterns. The disconnection of oligodendrocytes may be a distinctive feature of acute mTBI, whereas activation of the locus coeruleus represents a characteristic brainstem response to traumatic injury.

Methods

Animals

Specific pathogen-free male adult C57BL/6 mice aged approximately 6 weeks (17–20 g, n = 106) were obtained from the Beijing Vital River Laboratory Animal Technology. All experiments were approved by our local animal ethics committee (Beijing Neurosurgical Institute Experimental Animal Ethics Committee; approval no. 202103002). The experiments were designed and reported according to the Animal Research: Reporting of In Vivo Experiments guidelines. All the mice were housed in a controlled environment (five animals per cage; 12-h light/dark cycle; temperature, 22 ± 2 °C; humidity, 50–60%) with food and water ad libitum as reported in our previous study⁴¹. The mice were acclimated for 1 week before being divided without special selection into control (ctrl; n = 40), injury–1 h (n = 31), and injury–24 h (n = 35) groups.

Drop-weight impact model

For the mice in the brain-injury groups, the drop-weight impact model (DWI) was used to induce mild traumatic brain injury (mTBI) following and optimizing protocols from previous studies⁴². In brief, the mice were anaesthetized with 5% isoflurane (RWD Life Science) in 100% oxygen at a delivery rate of 1 L/min, and maintenance anaesthesia was then administered using 2% isoflurane (Sigma) inhalation. The animals were fixed in a prone position, and the skull was exposed using a midline incision. A spacer was then fixed on the surface of the bregma using dental cement. The impact was induced by the free fall of a weight (250 g) from a height of 2.6 cm onto the spacer. Mice that immediately fell into a coma or suffered from respiratory suppression met the standards for further analysis. After careful haemostasis, the skin incision was closed with 4-0 nylon suture thread, and the mice were transferred to their original controlled environment.

Tail-flick latency analysis

The tail-flick test measures an instinctive response to a stimulus and is widely used to evaluate brainstem function⁴³. In this procedure, we placed the mice into a cassette, leaving the tail outside, and approached it with a heat source. The tail-flick latency, which represents the time to tail-flick, was recorded. If the delay time exceeded 30 s, a value of 30 s was recorded for statistical analysis. Data were presented as means ± standard deviations and one-way ANOVA was used for statistical analysis. A p-value of <0.05 was considered statistically significant.

Immunofluorescence staining

Animals were fully anesthetized with 5% isoflurane and immediately sacrificed. The brains of the mice in each group were fixed in 4% paraformaldehyde for 24 h, followed by gradual dehydration and transparency procedures. After being embedded in paraffin wax, the tissues were cut into 5 μm sections. The procedure for immunofluorescence staining was performed as described in previous studies. After step-by-step dewaxing and rehydration, the following primary antibodies against the following antigens were used: Th (1:500; rabbit, ab137869; Abcam); Slc17a6 (1:500; rabbit, YN3008; Immunoway); Slc5a7 (1:500; rabbit, ab289635; Abcam); Dbh (1:500; rabbit, ab209487; Abcam); Reelin (1:100; rabbit, ab312310, Abcam); Glyt2 (1:100; rabbit, ab303556, Abcam); Vcan (1:500; rabbit, YT4874, Immunway); and Mbp(1:500; mouse, YM4830, Immunway). The tissues were incubated at 4 °C for 1 day, followed by incubation with the corresponding secondary antibodies conjugated to Alexa Fluor 488 or 594 (1:500; Abcam) for 1 h at room temperature. The samples were mounted using Vectashield containing DAPI for cell nuclei staining (ab104139; Abcam). Co-immunofluorescence was performed using a Treble-Fluorescence Immunohistochemical Mouse/Rabbit Kit according to the procedures specified by the manufacturer.

Fluorescent intensity of MBP and the number of positive neurons stained with Versican or Dbh were analysed using ImageJ software (version 1.53a). Data were presented as means ± standard deviations and statistical analysis was performed using one-way ANOVA with post hoc tests. Significance was defined as P < 0.05.

Virus injection and tracing

To determine the spatial association of the locus coeruleus and target nucleus projections, 1 μL of anterograde adeno-associated virus (AAV1) vector containing green fluorescent protein (GFP) was injected into the locus coeruleus of mice under 5% isoflurane anaesthesia (−5.45 mm anteroposterior (AP), +1.2 mm mediolateral (ML), −3.65 mm dorsoventral (DV)) using a stereotaxic instrument (RWD Life Science) in the control group⁴⁴. The injection speed was maintained at 0.1 μL/min, and the micropipette was retained for 10 min before retraction. Finally, a 10 μm sagittal tissue section was prepared, and fluorescently labelled neurons were identified using a fluorescence microscope (Axio Imager M2, Zeiss).

Western blot analysis

Animals were fully anesthetized with 5% isoflurane and promptly sacrificed. The whole brainstem tissue was collected from the three groups for Western blot (WB) analysis using gradient centrifugation as described previously⁴¹. Protein concentration was determined using the BCA method, and 30 μg of protein extract was separated by sodium dodecyl sulfate–polyacrylamide gel electrophoresis (SDS–PAGE). The separated proteins were then transferred onto polyvinylidene fluoride (PVDF) membranes. After incubation with 5% skim milk for 1 h at room temperature, the membranes were probed overnight at 4 °C with primary antibodies against Vcan (1:1000) and Mbp (1:1000). The PVDF membranes were then incubated with corresponding secondary antibodies (1:5000) for 1 hour at room temperature. The blot was visualized using chemiluminescence (Amersham Imager 600, GE). Relative band density was measured using ImageJ software (version 1.53a) and normalized to GAPDH band density.

Data were presented as means ± standard deviations and statistical analysis was performed using one-way ANOVA with post hoc analysis. Significance was defined as P < 0.05.

Flow cytometry analysis of cell ratios

Animals were anesthetized with 5% isoflurane and the total brainstem tissue was dissected and digested with enzymatic hydrolysate (Papain 5 mg/mL + DNase I 20 µg/mL). After incubation at 37 °C for 30 min, the enzyme reaction was terminated by adding serum-containing medium (DMEM + 10% FBS). The digested tissue was then filtered through a 70 µm cell strainer to remove undigested material and collect single-cell suspensions. The cells were fixed using BioLegend fixation buffer at room temperature for 20 min. The cell membrane was then permeabilized at room temperature for 15 min using BioLegend intracellular staining permeabilization wash buffer. The cells were then incubated with primary antibodies: Gad67 (1:100, mouse, ab26116, Abcam) and Glyt2 (1:100, rabbit, 23GB4145, Thermo Fisher) for 30 min, followed by incubation with the corresponding secondary antibodies for 30 min. Afterward, the cells were resuspended in 200 µL PBS. Flow cytometry analysis was performed using a Cytoflex S (Beckman Coulter, Inc., USA).

Cell ratio data were presented as means ± standard deviations and statistical analysis was performed using one-way ANOVA with post hoc analysis. Significance was defined as P < 0.05.

Single-nucleus RNA-sequencing

We utilized whole brainstem tissue from 20 male C57BL/6J mice (anesthetized with 5% isoflurane) for single-nucleus RNA-sequencing (snRNA-seq), with 5 merged in the control group, 5 merged in the injury−1 h group, and 10 merged in two batches in the injury−24 h group^45,46,47. All brainstem tissue from the midbrain to the foramen magnum was isolated.

Frozen tissue nuclear dissociation

Frozen brainstem samples were cut into <0.5 cm pieces and homogenized using a glass dounce tissue grinder (D8938, Sigma). The tissue was incubated and homogenized in ice-cold nuclei EZ lysis buffer. The nuclei were subsequently centrifuged at 500×g for 5 min at 4 °C. After centrifugation, the nuclear pellet was washed and resuspended in nuclear suspension buffer (NSB, consisting of 1 × PBS, 0.01% BSA, and 0.1% RNase inhibitor (2313 A, Clontech)), filtered through a 35 μm cell strainer (352235, Corning-Falcon) and counted. A final concentration of 1000 nuclei per μL was used for loading on a 10 × channel.

10 × library preparation and sequencing

Beads with a unique molecular identifier (UMI) and cell barcodes were loaded close to saturation so that each cell was paired with a bead in a gel beads-in-emulsion (GEM). After exposure to cell lysis buffer, polyadenylated RNA molecules hybridize to the beads. The beads were transferred to a single tube for reverse transcription. For cDNA synthesis, each cDNA molecule was tagged on the 3′ ends with UMI and a cell label indicating its cell of origin. Briefly, 10 × beads were then subjected to second-strand cDNA synthesis, adaptor ligation, and universal amplification. The sequencing libraries were prepared using randomly interrupted whole-transcriptome amplification products to enrich the 3′ ends of the transcripts linked with the cell barcode and UMI. All the remaining procedures, including the library construction, were performed according to the manufacturer’s standard protocol (Chromium Single Cell 3′ version 3.1). The sequencing libraries were quantified using a high-sensitivity DNA Chip (Agilent) on a Bioanalyzer 2100 and a Qubit high-sensitivity DNA assay (Thermo Fisher Scientific). The libraries were sequenced on a NovaSeq 6000 (Illumina) using 2 × 150 chemistry. Library construction and sequencing were performed by Shanghai Majorbio Bio-pharm Technology.

Single-nucleus RNA-sequence data processing

Reads were processed using the Cell Ranger version 7.1.0 pipeline with default and recommended parameters. FASTQs generated from Illumina sequencing output were aligned to the mouse genome, using the STAR algorithm⁴⁸. Next, gene–barcode matrices were generated for each individual sample by counting UMIs and filtering non-cell associated barcodes. Finally, we generate a gene-barcode matrix containing the barcoded cells and gene expression counts. The output was loaded as Seurat objects (version 5.1.0) in R (version 4.3.2) using the ‘CreateSeuratObject’ function⁴⁹. Low-quality cells in each object were discarded with the standards nFeature < 200 or >10,000, and the nuclei with >10% UMIs were mapped to mitochondrial genes. The data were normalized (normalizeData function in the Seurat package) to extract a subset of variable genes. Next, we integrated data from different samples after identifying ‘anchors’ between datasets using FindIntegrationAnchors and IntegrateData in the Seurat package. Principal component analysis (PCA) and ElbowPlot were subsequently performed to reduce the data to the top 20 PCA components. Clusters on a 2D map produced with uniform manifold approximation and projection (UMAP).

Identification of cell types and subtypes

The cells were clustered using graph-based clustering of the PCA-reduced data using the Louvain method after a shared nearest neighbour graph was computed. For each cluster, we used the Wilcoxon rank-sum test to find significantly deferentially expressed genes when the remaining clusters were compared. By setting the resolution to 0.7, we first separated the neuron and non-neuron clusters (e.g., oligodendrocytes, astrocytes, oligodendrocyte precursor cells (OPC), microglia, ependymal cells, endothelial cells, fibroblasts, and the choroid plexus). We then divided the neuron cluster into subsets and adjusted the resolution to 2 with the help of the Clustree package, identifying 40 sub-clusters. By selecting different markers, we finally pinpointed 7 types of neurons, including Cholinergic, Dopaminergic, GABAergic, Glutamatergic, Glycinergic, Noradrenergic, and Serotonergic neurons. The neuronal markers used were as follows^50,51:

Cholinergic neurons: Chat, Slc5a7;

Dopaminergic neurons: Th, Slc18a2;

GABAergic neurons: Gad1, Gad2;

Glutamatergic neurons: Slc17a6, Slc17a7, Slc17a8;

Glycinergic neurons: Slc6a5, Slc7a8;

Noradrenergic neurons: Slc6a2, Slc18a2; and

Serotonergic neurons: Slc6a4, Fev.

Spatial transcriptome sequence analysis

Sample selection and evaluation of RNA quality

We chose two brainstems in the control and injury-24 h groups, respectively, for further spatial transcriptome analysis. Each tissue sample was gently washed with cold PBS, filled with OCT (Sukura, 4583), and then sagittally cut into 10 μm sections (CM1950, Leica). The target region we selected was 830 μm lateral to the midline, spanning sagittally from the telencephalon to the medulla oblongata according to the Allen Mouse Brain Atlas (http://atlas.brain-map.org/)⁵². We collected the samples and optimized tissue experiments following the vision spatial tissue optimization panels (10× Genomics), including tissue staining, imaging, tissue permeabilization, fluorescent cDNA synthesis, RNA quality control (RIN > 7.0), and slide imaging.

Spatial transcriptomics using Visium Spatial

Gene expression (10× genomics)

The tissue slices of the samples were permeabilized for 30 min for gene expression. Libraries were prepared according to the Visium Spatial Gene Expression user guide and were loaded at 300 pM and sequenced on Xplus (Illumina) using a NovaSeq S4 Reagent Kit (200 cycles, catalogue no. 20027466, Illumina) at a sequencing depth of 20,000 mean read pairs/spot to determine the dilution for each sample. Libraries were generated from the cDNA and sequenced and the spatial barcodes were used to associate the reads with the tissue section images for spatial gene expression mapping. The sequencing and bioinformatics analysis were performed by Shanghai Majorbio Bio-pharm Technology.

ST RNA-seq data processing and identification of cell types and subtypes

Raw FASTQ files and histology images were processed by sample with the Space Ranger-1.3.1, which uses STAR for genome alignment, against the Cell Ranger hg38 reference genome⁵³. We evaluated the per-spot quality metrics using the Seurat (version 5.1.0) R Bioconductor package and excluded UMIs with fewer than 200 detected spots.

Gene expression from each voxel was normalized by the sctransform function in the Seurat package, which uses regularized negative binomial models to account for technical artefacts while preserving biological variance. The top 10 principal components in each slice were subsequently calculated and used to construct the KNN graph. The Louvain algorithm was used to cluster the voxels. We visualized the clusters on a 2D map produced with UMAP for each cluster. Wilcoxon rank-sum test was used to find significantly deferentially expressed genes when the remaining clusters were compared.

Integrating single-nucleus and Visium spatial gene expression data

Individual cell types were deconvoluted from the array of spots by using Robust Cell Type Decomposition (RCTD) with the annotated single-nucleus RNA-sequencing dataset from this study as a ref. ⁵⁴. In particular, RCTD identifies a set of characteristic genes that are uniquely expressed by each cell type, which assumes that the expression values of the spots follow a Poisson distribution. By combining the average expression values from single-nucleus data with proportion weights, RCTD constructs a hierarchical statistical model. This model yields the precise proportion of each cell type within each spot.

High-dimensional weighted gene co-expression network analysis (hdWGCNA)

The hdWGCNA is vital for revealing gene module associations and gene functional cooperativity⁵⁵. We performed hdWGCNA using total neuron clusters and oligodendrocyte clusters separately. The network was constructed by calculating the Pearson correlation matrix for all gene pairs and transforming it into an adjacency matrix. Soft-thresholding values of 7 and 4 were applied for the neuron cluster and oligodendrocyte cluster analyses, respectively. Modules were then defined by applying hierarchical clustering to the adjacency matrix and detecting distinct gene modules using the dynamic tree-cut algorithm. For each module, module eigengenes were calculated using the ModuleEigengenes function. These module eigengenes were correlated with external sample traits to identify significant associations (e.g., for the neuron clusters, we defined the cell type as traits; whereas, for the oligodendrocyte clusters, we defined the group as traits). Hub genes within each module are identified on the basis of their connectivity, and further, the functional enrichment analysis was performed to identify over-represented biological processes and pathways.

Differential gene expression analysis

We performed differential gene expression analysis both during the clustering step and to identify group differences related to injury time. Briefly, the normalized Seurat object was prepared, and the FindMarkers function was used on the basis of the ‘RNA’ assay. To identify cluster differential genes, we set the ‘min.pct’ as 0.25, and the threshold of logFC as 0.25. This ensured that the selected genes were expressed in at least 25% of the cells in the cluster and presented at least a 19% change in expression. We used a Wilcoxon test to identify group-related differentially expressed genes (DEGs). By setting ident.1 and ident.2, we evaluated the DEGs by comparing ident.1 with ident.2.

Gene Ontology (GO) and Kyoto Encyclopaedia of Genes and Genomes (KEGG) enrichment analyses

The GO and KEGG pathway enrichment analyses were conducted using the R package ClusterProfiler version 4.12.0. A mouse-specific database was used for annotation in enrichment analysis. Genes were selected using a false discovery rate (FDR) model, with a P value and adjusted P value both set to be <1. For the GO analysis, we focused on the biological process (BP) and molecular function (MF) categories.

Gene set enrichment analysis (GSEA)

GSEA is an unsupervised method with sensitivity to subtle changes. The DEGs were ranked using GSEA version 4.3.3, which utilized a Mouse_Gene_Symbol_Remapping_Human _Orthologs_MSigDB version 2023.2. Hs. chip to remap the gene symbols to their human orthologues. The gene set databases used were c5.go.mf.v2023.2.Hs.symbols.gmt for molecular function and c5.go.bp.v2023.2.Hs.symbols.gmt for biological processes. NES (normalized enrichment score) values were calculated, and a |NES| > 1 with NOM p < 0.05 was considered significant.

Inference cell–cell interactions

We used CellChat 2.1.2 to infer cell–cell communication between oligodendrocytes and neurons, and among neuron subtypes. The CellChat package incorporates a curated database of ligand–receptor interactions, covering a wide range of cell types and signalling pathways. This includes secreted signalling, ECM–receptor interactions, cell–cell contact, and non-protein signalling, such as neurotransmitter–receptor associations⁵⁶. For the scRNA-seq data, we extracted oligodendrocytes and neurons as a new Seurat object. The inference of cell–cell interactions included the identification of overexpressed genes, computation of the communication probability, and identification of significant cellular communications.

Trajectory analysis

To predict injury-time-related changes in oligodendrocytes and Noradrenergic neurons and to identify injury-related features, we performed Slingshot (2.10.0)⁵⁷ and Monocle 2 (2.30.1) analyses⁵⁸. Following dimensionality reduction, the hierarchical clustering of oligodendrocytes was performed on the UMAP coordinates. Slingshot was used to fit a minimum spanning tree (MST) to these clusters and determine the approximate trajectory structure. Then, the linear trajectory was smoothed using simultaneous principal curves to obtain the final trajectories and pseudo-time values. The groups were fitted to the pseudo-time to highlight injury-related cluster features. Finally, genes that changed along the trajectory were extracted to generate a heatmap, and enrichment analysis was performed. For the Noradrenergic neurons containing relatively few cells, we implemented Monocle 2, which is widely used and suitable for comprehensive analysis. Primarily, the ‘newCellDataSet’ function was utilized to create the object using cell type and group information. Then, genes that were differentially expressed between different groups were selected for the trajectory analysis using the ‘differentialGeneTest’ function. This was followed by dimensionality reduction with the ‘DDRTree’ method and cell ordering with the ‘orderCells’ function. Finally, genes involved in the pseudo-time branches were extracted for enrichment analysis.

Neuron activity score at the spatial level

We were particularly interested in the changes in neuron activity following injury, as the response of neurons in the brainstem remains unclear. Therefore, the ‘AddModuleScore’ function was applied to rank gene scores representing neuron activity on the basis of their spatial ___location in the control and injury groups⁵⁹. The hallmark genes we used were involved in neurotrophin signalling, synaptic transmission, and glycometabolism. These genes included Fos, Arc, Egr1, Syp, Bdnf, Npas4, Slc2a1, Slc2a3, Hk1, Pfkp, Slc1a, and Creb1. The spot colours in the ST RNA-seq data ultimately indicate the neuron activity scores.

Transcription factor analysis

Gene regulation was analysed using SCENIC (1.3.1) to investigate the transcriptional regulation features in the locus coeruleus after traumatic injury⁶⁰. Next, regulatory modules (regulons) were generated using co-expression and DNA motif analysis. The activity of the regulons was scored in each cell using the ‘AUCell’ function. Finally, we calculated the Spearman correlation coefficients between regulons. The regulons that had a correlation coefficient of >0.5 and were activated in at least 50% of the cells in each group were selected to show gene regulatory module features.

Statistics and reproducibility

All experimental studies were repeated at least three times per group, and the sample size is reported in the Source Data file.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Data availability

All relevant sequencing files used to produce the results in the present study are publicly available on SRA, and the BioProject accession number is PRJNA1173913. (http://www.ncbi.nlm.nih.gov/sra/PRJNA1173913). Source data are provided with this paper.

Code availability

The present study did not use any custom code, and all analyses were based on published code, which is properly referenced.

References

Ruder, L. & Arber, S. Brainstem circuits controlling action diversification. Annu. Rev. Neurosci. 42, 485–504 (2019).
Article CAS PubMed Google Scholar
Miller, J. R. et al. A subregion of the parabrachial nucleus partially mediates respiratory rate depression from intravenous remifentanil in young and adult rabbits. Anesthesiology 127, 502–514 (2017).
Article CAS PubMed Google Scholar
Krohn, F. et al. The integrated brain network that controls respiration. Elife 12, https://doi.org/10.7554/eLife.83654 (2023).
Benarroch, E. E. Brainstem integration of arousal, sleep, cardiovascular, and respiratory control. Neurology 91, 958–966 (2018).
Article PubMed Google Scholar
Paton, J. F. R. et al. Advancing respiratory-cardiovascular physiology with the working heart-brainstem preparation over 25 years. J. Physiol. 600, 2049–2075 (2022).
Article CAS PubMed Google Scholar
Lu, L. et al. Control of locomotor speed, arousal, and hippocampal theta rhythms by the nucleus incertus. Nat. Commun. 11, 262 (2020).
Article ADS CAS PubMed PubMed Central Google Scholar
Kaur, S. & Saper, C. B. Neural circuitry underlying waking up to hypercapnia. Front. Neurosci. 13, 401 (2019).
Article PubMed PubMed Central Google Scholar
Nardone, S. et al. A spatially-resolved transcriptional atlas of the murine dorsal pons at single-cell resolution. Nat. Commun. 15, 1966 (2024).
Article ADS CAS PubMed PubMed Central Google Scholar
Benghanem, S. et al. Brainstem dysfunction in critically ill patients. Crit. Care 24, 5 (2020).
Article PubMed PubMed Central Google Scholar
Karakis, I. Brainstem mapping. J. Clin. Neurophysiol. 30, 597–603 (2013).
Article PubMed Google Scholar
Marklund, N. et al. Treatments and rehabilitation in the acute and chronic state of traumatic brain injury. J. Intern. Med. 285, 608–623 (2019).
Article CAS PubMed PubMed Central Google Scholar
Flerlage, W. J. et al. Effects of repetitive mild traumatic brain injury on corticotropin-releasing factor modulation of lateral habenula excitability and motivated behavior. J Neurotrauma. 42, 832–850 (2024).
Capizzi, A., Woo, J. & Verduzco-Gutierrez, M. Traumatic brain injury: an overview of epidemiology, pathophysiology, and medical management. Med Clin. North Am. 104, 213–238 (2020).
Article PubMed Google Scholar
Callaway, C. C. M. & Kosofsky, B. E. Autonomic dysfunction following mild traumatic brain injury. Curr. Opin. Neurol. 32, 802–807 (2019).
Article CAS PubMed Google Scholar
Kim, E. et al. Reduced brainstem volume after mild traumatic brain injury. Am. J. Phys. Med. Rehabil. 100, 473–482 (2021).
Article PubMed Google Scholar
Peng, H. et al. Morphological diversity of single neurons in molecularly defined cell types. Nature 598, 174–181 (2021).
Article ADS CAS PubMed PubMed Central Google Scholar
Velmeshev, D. et al. Single-cell analysis of prenatal and postnatal human cortical development. Science 382, eadf0834 (2023).
Article CAS PubMed PubMed Central Google Scholar
Huang, M. et al. SNAP-25, but not SNAP-23, is essential for photoreceptor development, survival, and function in mice. Commun. Biol. 7, 34 (2024).
Article CAS PubMed PubMed Central Google Scholar
Kiss, T. et al. Single-cell RNA sequencing identifies senescent cerebromicrovascular endothelial cells in the aged mouse brain. Geroscience 42, 429–444 (2020).
Article PubMed PubMed Central Google Scholar
Li, Q. et al. Developmental heterogeneity of microglia and brain myeloid cells revealed by deep single-cell RNA sequencing. Neuron 101, 207–223.e210 (2019).
Article CAS PubMed Google Scholar
Allen, W. E., Blosser, T. R., Sullivan, Z. A., Dulac, C. & Zhuang, X. Molecular and spatial signatures of mouse brain aging at single-cell resolution. Cell 186, 194–208.e118 (2023).
Article CAS PubMed Google Scholar
Lim, L., Mi, D., Llorca, A. & Marin, O. Development and functional diversification of cortical interneurons. Neuron 100, 294–313 (2018).
Article CAS PubMed PubMed Central Google Scholar
Gouwens, N. W. et al. Integrated morphoelectric and transcriptomic classification of cortical GABAergic cells. Cell 183, 935–953 e919 (2020).
Article CAS PubMed PubMed Central Google Scholar
Marques, S. et al. Oligodendrocyte heterogeneity in the mouse juvenile and adult central nervous system. Science 352, 1326–1329 (2016).
Article ADS CAS PubMed PubMed Central Google Scholar
Emery, B. & Wood, T. L. Regulators of oligodendrocyte differentiation. Cold Spring Harb. Perspect. Biol. 16, a041358 (2024).
Article PubMed Google Scholar
Breton-Provencher, V., Drummond, G. T. & Sur, M. Locus coeruleus norepinephrine in learned behavior: anatomical modularity and spatiotemporal integration in targets. Front. Neural Circuits 15, 638007 (2021).
Article CAS PubMed PubMed Central Google Scholar
Humble, S. S. et al. Prognosis of diffuse axonal injury with traumatic brain injury. J. Trauma Acute Care Surg. 85, 155–159 (2018).
Article PubMed PubMed Central Google Scholar
Moe, H. K. et al. The influence of traumatic axonal injury in thalamus and brainstem on level of consciousness at scene or admission: a clinical magnetic resonance imaging study. J. Neurotrauma 35, 975–984 (2018).
Article PubMed PubMed Central Google Scholar
Herring, C. A. et al. Human prefrontal cortex gene regulatory dynamics from gestation to adulthood at single-cell resolution. Cell 185, 4428–4447.e4428 (2022).
Article CAS PubMed Google Scholar
Trevino, A. E. et al. Chromatin and gene-regulatory dynamics of the developing human cerebral cortex at single-cell resolution. Cell 184, 5053–5069.e5023 (2021).
Article CAS PubMed Google Scholar
Koupourtidou, C. et al. Shared inflammatory glial cell signature after stab wound injury, revealed by spatial, temporal, and cell-type-specific profiling of the murine cerebral cortex. Nat. Commun. 15, 2866 (2024).
Article ADS CAS PubMed PubMed Central Google Scholar
Frankowski, J. C. et al. Brain-wide reconstruction of inhibitory circuits after traumatic brain injury. Nat. Commun. 13, 3417 (2022).
Article ADS CAS PubMed PubMed Central Google Scholar
Pereira, M., Birtele, M. & Rylander Ottosson, D. Direct reprogramming into interneurons: potential for brain repair. Cell Mol. Life Sci. 76, 3953–3967 (2019).
Article CAS PubMed PubMed Central Google Scholar
Kumru, H., Vidal, J., Kofler, M., Portell, E. & Valls-Sole, J. Alterations in excitatory and inhibitory brainstem interneuronal circuits after severe spinal cord injury. J. Neurotrauma 27, 721–728 (2010).
Article PubMed Google Scholar
Yoo, R. E. et al. Myelin content in mild traumatic brain injury patients with post-concussion syndrome: quantitative assessment with a multidynamic multiecho sequence. Korean J. Radiol. 23, 226–236 (2022).
Article PubMed PubMed Central Google Scholar
Poe, G. R. et al. Locus coeruleus: a new look at the blue spot. Nat. Rev. Neurosci. 21, 644–659 (2020).
Article CAS PubMed PubMed Central Google Scholar
Fortin, S. M. et al. The locus coeruleus contributes to the anorectic, nausea, and autonomic physiological effects of glucagon-like peptide-1. Sci. Adv. 9, eadh0980 (2023).
Article CAS PubMed PubMed Central Google Scholar
Berger, A. et al. How is the norepinephrine system involved in the antiepileptic effects of vagus nerve stimulation?. Front. Neurosci. 15, 790943 (2021).
Article PubMed PubMed Central Google Scholar
Hachem, L. D., Wong, S. M. & Ibrahim, G. M. The vagus afferent network: emerging role in translational connectomics. Neurosurg. Focus 45, E2 (2018).
Article PubMed Google Scholar
Farrand, A., Jacquemet, V., Verner, R., Owens, M. & Beaumont, E. Vagus nerve stimulation parameters evoke differential neuronal responses in the locus coeruleus. Physiol. Rep. 11, e15633 (2023).
Article CAS PubMed PubMed Central Google Scholar
Zhuang, Y. et al. Megf10-related engulfment of excitatory postsynapses by astrocytes following severe brain injury. CNS Neurosci. Ther. 29, 2873–2883 (2023).
Article CAS PubMed PubMed Central Google Scholar
Lu, S. et al. Decoupling the mutual promotion of inflammation and oxidative stress mitigates cognitive decline and depression-like behavior in rmTBI mice by promoting myelin renewal and neuronal survival. Biomed. Pharmacother. 173, 116419 (2024).
Article CAS PubMed Google Scholar
Borrelli, K. N. et al. Sex differences in behavioral and brainstem transcriptomic neuroadaptations following neonatal opioid exposure in outbred mice. eNeuro. 8, ENEURO.0143–21 (2021).
Takeuchi, T. et al. Locus coeruleus and dopaminergic consolidation of everyday memory. Nature 537, 357–362 (2016).
Article ADS CAS PubMed PubMed Central Google Scholar
Wang, H. et al. Antiandrogen treatment induces stromal cell reprogramming to promote castration resistance in prostate cancer. Cancer Cell 41, 1345–1362.e1349 (2023).
Article CAS PubMed Google Scholar
Niu, N. et al. Tumor cell-intrinsic epigenetic dysregulation shapes cancer-associated fibroblasts heterogeneity to metabolically support pancreatic cancer. Cancer Cell 42, 869–884.e869 (2024).
Article CAS PubMed Google Scholar
Wang, Z. et al. Leucine-tRNA-synthase-2-expressing B cells contribute to colorectal cancer immunoevasion. Immunity 55, 1067–1081 e1068 (2022).
Article CAS PubMed Google Scholar
Le, D. T. et al. Mismatch repair deficiency predicts response of solid tumors to PD-1 blockade. Science 357, 409–413 (2017).
Article ADS CAS PubMed PubMed Central Google Scholar
Stuart, T. et al. Comprehensive integration of single-cell data. Cell 177, 1888–1902.e1821 (2019).
Article CAS PubMed PubMed Central Google Scholar
Wang, K., Cai, B., Song, Y., Chen, Y. & Zhang, X. Somatosensory neuron types and their neural networks as revealed via single-cell transcriptomics. Trends Neurosci. 46, 654–666 (2023).
Article CAS PubMed Google Scholar
Zhang, M. et al. Molecularly defined and spatially resolved cell atlas of the whole mouse brain. Nature 624, 343–354 (2023).
Article ADS CAS PubMed PubMed Central Google Scholar
Wang, Q. et al. The allen mouse brain common coordinate framework: a 3D reference atlas. Cell 181, 936–953.e920 (2020).
Article CAS PubMed PubMed Central Google Scholar
Hafemeister, C. & Satija, R. Normalization and variance stabilization of single-cell RNA-seq data using regularized negative binomial regression. Genome Biol. 20, 296 (2019).
Article CAS PubMed PubMed Central Google Scholar
Cable, D. M. et al. Robust decomposition of cell type mixtures in spatial transcriptomics. Nat. Biotechnol. 40, 517–526 (2022).
Article CAS PubMed Google Scholar
Morabito, S., Reese, F., Rahimzadeh, N., Miyoshi, E. & Swarup, V. hdWGCNA identifies co-expression networks in high-dimensional transcriptomics data. Cell Rep. Methods 3, 100498 (2023).
Article CAS PubMed PubMed Central Google Scholar
Jin, S. et al. Inference and analysis of cell–cell communication using CellChat. Nat. Commun. 12, 1088 (2021).
Article ADS CAS PubMed PubMed Central Google Scholar
Bi, K. et al. Tumor and immune reprogramming during immunotherapy in advanced renal cell carcinoma. Cancer Cell 39, 649–661.e645 (2021).
Article CAS PubMed PubMed Central Google Scholar
Qiu, X. et al. Reversed graph embedding resolves complex single-cell trajectories. Nat. Methods 14, 979–982 (2017).
Article CAS PubMed PubMed Central Google Scholar
Liu, J., Shi, Y. & Zhang, Y. Multi-omics identification of an immunogenic cell death-related signature for clear cell renal cell carcinoma in the context of 3P medicine and based on a 101-combination machine learning computational framework. EPMA J. 14, 275–305 (2023).
Article PubMed PubMed Central Google Scholar
Bravo Gonzalez-Blas, C. et al. SCENIC+: single-cell multiomic inference of enhancers and gene regulatory networks. Nat. Methods 20, 1355–1367 (2023).
Article CAS PubMed PubMed Central Google Scholar

Download references

Acknowledgements

This study was supported by the Beijing Natural Science Foundation (Grant No. J230008) and the National Key R&D Programme of China (Grant No. 2023YFC2410605) for Guoyi Gao. The National Natural Science Foundation of China (Grant Nos. 82401612 for Yuan Zhuang and 82302792 for Xiang Wu). The authors thank Xiaolong He, Guoyu Chu, and Yuling Jiang from Shanghai Majorbio Bio-Pharm Technology Co., Ltd. for their technical assistance in the scRNA-seq and spatial transcriptomic analyses.

Author information

Authors and Affiliations

Department of Neurosurgery, Beijing Tiantan Hospital, Capital Medical University, Beijing, China
Yuan Zhuang, Xixian Liao, Ming Li, Yu Yan, Chuanhang He, Runfa Tian & Guoyi Gao
Beijing Key Laboratory of Central Nervous System Injury, Beijing Neurosurgical Institute, Capital Medical University, Beijing, China
Fei Niu & Guoyi Gao
Department of Neurosurgery, Shanghai General Hospital, School of Medicine, Shanghai Jiao Tong University, Shanghai, China
Xiang Wu

Authors

Yuan Zhuang
View author publications
Search author on:PubMed Google Scholar
Xixian Liao
View author publications
Search author on:PubMed Google Scholar
Fei Niu
View author publications
Search author on:PubMed Google Scholar
Ming Li
View author publications
Search author on:PubMed Google Scholar
Yu Yan
View author publications
Search author on:PubMed Google Scholar
Chuanhang He
View author publications
Search author on:PubMed Google Scholar
Xiang Wu
View author publications
Search author on:PubMed Google Scholar
Runfa Tian
View author publications
Search author on:PubMed Google Scholar
Guoyi Gao
View author publications
Search author on:PubMed Google Scholar

Contributions

G.Y.G. and W.X. designed the study; Y.Z., M.L. and Y.Y. conducted the single-nucleus and spatial transcriptomic data analysis; Y.Z., X.X.L. and C.H.H. performed the experimental study; Y.Z. made manuscript preparation; G.Y.G., F.N. and R.F.T. made manuscript review. All authors approved the final version of the manuscript.

Corresponding authors

Correspondence to Xiang Wu, Runfa Tian or Guoyi Gao.

Ethics declarations

Competing interests

The authors declare no competing interests.

Peer review

Peer review information

Nature Communications thanks Kevin Wang and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. A peer review file is available.

Additional information

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Supplementary Information

Reporting Summary

Transparent Peer Review file

Source data

Source Data

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

About this article

Cite this article

Zhuang, Y., Liao, X., Niu, F. et al. Single-nucleus and spatial signatures of the brainstem in normal brain and mild traumatic brain injury in male mice. Nat Commun 16, 5082 (2025). https://doi.org/10.1038/s41467-025-59856-2

Download citation

Received: 06 September 2024
Accepted: 07 May 2025
Published: 31 May 2025
DOI: https://doi.org/10.1038/s41467-025-59856-2