Inhibitory Neurons in Nucleus Tractus Solitarius Are Involved in Decrease of Heart Rate Variability and Development of Depression-Like Behaviors in Temporal Lobe Epilepsy

Abstract

Background

Diminished heart rate variability (HRV) has been observed in epilepsy, especially in epilepsy with depressive disorders. However, the underlying mechanism remains elusive.

Methods

We studied HRV, spontaneous recurrent seizures, and depression-like behaviors in different phases of pilocarpine-induced temporal lobe epilepsy (TLE) in mice. Single-cell RNA sequencing analysis was used to identify various nerve cell subsets in TLE mice with and without depression. Differentially expressed gene (DEG) analysis was performed in epilepsy, depression, and HRV central control–related brain areas.

Results

We found decreased HRV parameters in TLE mice, and alterations were positively correlated with the severity of depression-like behaviors. The severity of depression-like behaviors was correlated with the frequency of spontaneous recurrent seizure. Characteristic expression of mitochondria-related genes was significantly elevated in mice with depression in glial cells, and the enrichment analysis of those DEGs showed an enriched GABAergic synapse pathway in the HRV central control–related brain area. Furthermore, inhibitory neurons in the nucleus tractus solitarius, which is an HRV central control–related brain area, were specifically expressed in TLE mice combined with depression compared with those in mice without depression. A significantly enriched long-term depression pathway in DEGs from inhibitory neurons was found.

Conclusions

Our study reported correlations between HRV and epilepsy–depression comorbidity in different phases of TLE. More importantly, we found that HRV central control-related inhibitory neurons are involved in the development of depression in TLE, providing new insights into epilepsy comorbid with depression.

INTRODUCTION

Depression is one of the most common psychiatric comorbidities in people with epilepsy (Fiest et al., 2013), particularly among patients with drug‐resistant temporal lobe epilepsy (TLE) (Tellez-Zenteno et al., 2007; Kanner et al., 2012; Kanner, 2014). Correlations between the development of seizures and depression‐like behaviors have been demonstrated in various animal models of epilepsy, suggesting a relationship between depression and epilepsy (Kanner et al., 2012; Hester and Danzer, 2014; D’Alessio et al., 2019). The underlying pathogenic mechanisms remain unknown, but alterations of the hypothalamus–pituitary–adrenal axis, neurotransmitter disturbances, immunologic disturbances, and genetic factors (Heiman et al., 2010; Lang and Borgwardt, 2013; Hester and Danzer, 2014; Kanner, 2014; Lacey et al., 2014; Wulsin et al., 2016) have been frequently observed in major depression and likely constitute a common pathway that may also be disturbed in the combination of epilepsy and depression.

Decreased heart rate variability (HRV), as a result of stressful events perpetuating the rhythms of sympathetic nerve–dominant states (Ernst, 2017), has been observed in epilepsy (Lotufo et al., 2012; Romigi et al., 2016; Arbune et al., 2020; Dono et al., 2020; Powell et al., 2021). Furthermore, previous work suggested that abnormal HRV parameters were associated with depressive disorders (Sgoifo et al., 2015; Borrione et al., 2018), with a plausible mechanistic link with autonomic nervous system (ANS) dysfunction (Shaffer et al., 2014). We have found a correlation between HRV alterations and epilepsy depression status in previous clinical work (Aimaier et al., 2022). However, the underlying mechanism of diminished HRV in epilepsy comorbid with depression is not yet clear.

In this study, we hypothesized that seizure activity in chronic TLE is a chronic endogenous stress that may lead to alterations in HRV and depression-like behaviors through a common pathway. Following this hypothesis, we explored the variation tendency of depression-like behaviors and HRV parameters in different phases of the chronic TLE mouse model. Single-cell RNA sequencing (scRNA-seq) was performed in epilepsy, depression, and HRV central control–related brain areas in TLE mice with and without depression. We found that the alterations in HRV parameters were correlated with the severity of depression-like behaviors in TLE mice and found upregulated inhibitory neurons in the HRV central control–related brain area nucleus tractus solitarius (NTS) in TLE mice with depression compared with those without depression, which provides a new target for the development of depression and HRV alterations in epilepsy.

MATERIALS AND METHODS

Experimental Design

We studied HRV, spontaneous recurrent seizure (SRS), and depression-like behaviors 1 week before status epilepticus (SE) induction (baseline stage) and every 4 weeks after SE induction (Figure 1A). The electrophysiological signal recording pin was implanted into male adult C57BL/6J mice (6–8 weeks old) 2 weeks before SE induction, and baseline stage recording was started 1 week after recovery from the surgery. Mouse brain tissue collection was performed for scRNA-seq 12 weeks after SE induction (Figure 1A).

Mice

All experimental procedures followed the guidelines of the National Institutes of Health and were approved by the Animal Care and Use Committee at the Institute of Neuroscience, Chinese Academy of Sciences. All experiments were performed using young mice (only male mice were used, >6 weeks at the time of electrode implantation surgery). Mice were housed in a 12-hour-light/-dark cycle (lights on at 7 am) with food and water available ad libitum. Mice with implants for electroencephalogram (EEG)/electrocardiograph (ECG)/electromyography (EMG) recording were individually housed. All of the C57BL/6J mice were purchased from institute-approved vendors (Shanghai Silaike or LingChang Experiment Animal Co., Shanghai, China). A total of 53 adult male mice were used in the experiment, and 26 of them died during pilocarpine-injected SE induction (supplementary Table 1). No mice died during the latent period, but 2 of them died in SRS attack (supplementary Table 1). The remaining 25 mice were observed, and behavioral and electrophysiological data were recorded throughout the experimental cycle (14 weeks).

Surgical Procedures

All surgical procedures were conducted under general anesthesia using continuous isoflurane (5% for induction; 1.5%–2% for maintenance). Depth of anesthesia was continuously monitored and adjusted when necessary. Following induction of anesthesia, the mice were placed on a stereotaxic frame with a heating pad. After anaesthetization, the skin was incised to expose the skull. For electrophysiological signal recordings, a reference screw was inserted into the skull on top of the left cerebellum. EEG recordings were made from 2 screws on top of the left and right cortex. Two EMG electrodes were inserted into the neck musculature. The mice were rotated onto their right side and the most obvious point of heartbeat on the left side was selected. A small area (approximately 2 cm²) around the point was shaved, and the long ECG wire was extended down to the shaved point through the subcutaneous region. The ECG wire was sutured to the surrounding muscle tissue. Insulated leads from the EEG/EMG/ECG electrodes were soldered to a pin header, which was secured to the skull using dental cement.

Establishment of the Pilocarpine-Induced SE Mouse Model

The pilocarpine-induced SE mouse model was established as described previously (Arshad and Naegele, 2020). Briefly, mice received IP injections of scopolamine methyl bromide (1 mg/kg) to reduce peripheral muscarinic effects. Thirty minutes later, the muscarinic agonist pilocarpine (300–400 mg/kg; the specific dose of pilocarpine to individual mice is shown in supplementary Table 1) was IP injected to induce SE. Seizures started 10 to 20 minutes after the pilocarpine injection, and their pathological behavior was visually monitored using a scale adapted from Racine’s standard criteria (Racine, 1972). After 3 to 5 stage 4 or 5 Racine scale seizures (supplementary Table 2), the mouse typically developed SE, defined as a continuous state of seizures (Arshad and Naegele, 2020). We noted the time of SE and terminated seizures at 30 minutes to 1 hour after SE onset by injecting diazepam solution (0.01–0.02 mL, IP); otherwise, they were excluded from the experiment. We administered 1 mL of 5% dextrose per individual mouse (IP; 1 mL) to provide an energy source and hydration that could increase the survival rate (Kim and Cho, 2018). As the higher dosage of pilocarpine significantly increased the mortality rate, the success rate of validated SRS occurrence was not increased (supplementary Table 1). Based on the above correlations among the dose of pilocarpine, mortality, and success rate of SRS, we recommend a concentration between 320 mg/kg and 360 mg/kg.

Video Monitoring SRS

Pilocarpine administration to mice results in SE, and after a latency period, SRS occurs. The mice that survived the first 2 weeks after SE induction were monitored with a video surveillance system. The monitoring period lasted 10 weeks (24 h/day from 2 to 12 weeks after SE induction) (Figure 1A). Only seizures with a score ≥4 on the Racine scale (supplementary Table 2) were counted as valid SRS. Mice with more than 1 valid SRS attack were defined as an effective chronic TLE mouse model.

Physiological Signal Recordings

Physiological signals EEG/EMG/ECG were recorded in a free-roaming state in home cages placed in sound-attenuating boxes (Figure 1B). EEG/EMG/ECG electrodes were connected to flexible recording cables via a mini-connector. Recordings started after 20 to 30 minutes of habituation. The signals were recorded with a TDT amplifier (sampling rate, 610 Hz). Spectral analysis was carried out using fast Fourier transform, and brain states were classified into wakefulness and nonrapid eye movement states.

Depression-Like Behavior Tests

Sucrose Preference Test (SPT)—

This test assesses the depression-like behavior of anhedonia based on the innate preference of rodents toward sweets (Pucilowski et al., 1993). Before the test, 2 identical bottles, 1 with pure water and 1 with water containing 2% sucrose, were placed in the cage for at least 48 hours to adapt to the 2 sources of drink. After adaptation, the mice were deprived of water for 24 hours. On the day of the test, the 2 identical bottles were placed in the cage for 2 hours, and at 1 hour after the onset of the test the locations of the 2 bottles were exchanged (Figure 1F). After the test, the liquid intake was calculated as follows: sucrose consumption rate (SucroRate) = sucrose consumption/ (sucrose consumption + water consumption). A low SucroRate was indicative of the state of anhedonia. Mice that were inborn insensitive to sweets (SucroRate was <50% when the SPT was repeated twice before SE induction) were excluded.

Tail Suspension Test (TST)—

This test assesses the depression-like behavior of despair based on the immobility time (IMT) (Can et al., 2012). The mouse’s tail was suspended from a plastic rod with tape, and the mouse was placed in a head-down position. We controlled the distance from the ground to the head of the mouse to approximately 30 cm (Figure 1G). The experiment lasted 6 minutes. During the 6-minute period of the test, each mouse was adapted for 1 minute after being suspended, and the amount of IMT was recorded during the remaining 5 minutes. A longer IMT was indicative of the state of despair.

To avoid the influence of seizures on the outcome of the behavioral test, the SPT/TST was performed after verifying that no seizures had developed for at least 6 hours before the tests.

HRV Analysis

Calculation of HRV parameters was carried out with an open-source MATLAB toolbox for analyzing HRV (MarcusVollmer-HRV-1b). Considering the possible effects of wakefulness and sleep stages on HRV parameters mentioned before (Herzig et al., 2017), the indices were compared in the same stage. Consecutive 1-minute period ECG data without artifacts and arousals from the wakefulness and nonrapid eye movement sleep stages were analyzed.

The common HRV time-domain measures were calculated: SDNN (standard deviation of all RR intervals), RMSSD (root mean square of the difference of adjacent RR intervals), and probability of RR intervals >50 milliseconds. SDNN represents a global measure of HRV and provides information about all HRV components. RMSSD is considered a powerful measure of high-frequency power (HF, 0.15–0.4 Hz) variations in short-term recording, as it provides a useful evaluation of high-frequency power and vagal tone (Task Force of the European Society of Cardiology and the North American Society of Pacing and Electrophysiology, 1996). Probability of RR intervals >50 milliseconds is the percentage of consecutive intervals that differ by more than 50 milliseconds, which reflects fast, high-frequency variability changes (Georgieva-Tsaneva et al., 2020).

scRNA-seq Analysis and Data Processing

C57BL/6J adult male mice with SRS were divided into groups with and without depression according to whether they developed depression-like behaviors after SE induction. Brain slices were prepared from the epilepsy control–related brain area hippocampus (HPC) (Liu et al., 2008), depression control-related area prefrontal cortex (PFC) (Pizzagalli and Roberts, 2022) and ventral tegmental area (VTA) (Douma and de Kloet, 2020), and HRV central control–related area NTS (Benarroch, 1993). The brain slices were dissected and dissociated using the Papain Dissociation System (Worthington, Lakewood, NJ, USA; Cat#LK003153) according to the manufacturer’s instructions. The cell suspension was filtered with a 40-mm filter (Thermo Fisher Scientific, Cat# 352340) and resuspended in DMEM (Thermo Fisher Scientific, Cat# 12800017) containing 10% fetal bovine serum (Thermo Fisher Scientific, Cat# C838T52). Single cells were captured using 10X Chromium (10X Genomics, Pleasanton, CA, USA), and libraries were prepared according to the manufacturer’s instructions (Chromium Single Cell 30 Library & Gel Bead Kit v2, 10× Genomics).

Sequencing was performed on an Illumina NovaSeq 6000 (Novogene, Beijing, China). Sequencing data were processed with 10× Cell Ranger (10× Genomics), and the sequencing depth of different conditions and batches was normalized. RunTSNE and FindClusters functions were used for visualization and clustering. The cell markers of different clusters were identified using the FindConservedMarkers function. The neuron-responsive genes in each cell cluster were identified using the FindMarkers function.

Quantification and Statistical Analysis

Statistical analysis was performed using GraphPad Prism 8 (GraphPad Software, La Jolla, CA, USA). No randomization algorithm was used, although mice were randomly assigned to experimental conditions.

Data were analyzed blinded to the experimental condition. Data are presented as the mean ± SD. The Gaussian distribution of the data was assessed using the Shapiro‒Wilk normality test or KS normality test. Data matrices were built with the obtained behavioral or physiological data and partial least squares regression using SOLO 8.6 chemometric software from Eigenvector Research (Manson, WA, USA). We utilized R packages to filter out the data of the matrix to transcriptome and normalize the data. The percentage of the gene numbers, cell counts, and mitochondria sequencing count were calculated. Per Seurat’s recommendation for quality control, genes expressed in fewer than 3 cells and cells containing fewer than 50 genes were discarded. To remove potential cell doublets and low-quality cells, cells with a percentage of mitochondrial genes >5% were filtered out (Ilicic et al., 2016). The standards for identifying differentially expressed gene (DEG) were a fold change ≥2 and a false discovery rate (FDR) <0.05.

RESULTS

SRS in Pilocarpine-Induced SE in C57BL/6J Mice

The development of SRS in 25 SE-induced animals was evaluated by continuous video monitoring. Overall, we observed a total of 42 SRS attacks from 7 animals (Figure 1C), all of which had a seizure scale score ≥4 (supplementary Table 2). All the SRSs were “convulsive,” accompanied by unilateral or bilateral forelimb clonus and balance loss after rearing. We did not observe significant “convulsive” SRS in the remaining 18 animals during the recording periods (Figure 1C). Overall, all the remarkable SRSs developed in the chronic period (at least 4 weeks after SE induction), with a mean latency to the first SRS attack of 85.3 ± 21.2 days. The overall frequency and mean duration of SRS attacks increased from week 4 to week 12 of the recording period (Figure 1D). SRS occurs mostly during the sleep-wake transition in mice (Figure 1E).

Depression-Like Behaviors in Mice With SRS

The SucroRate of the SPT, an indicator of anhedonia, was essentially the same in the mice with and without SRS at baseline and 4 weeks after SE induction (0.83 ± 0.06 vs 0.81 ± 0.11). However, with the occurrence of SRS, the mice with SRS showed a significant decrease in SucroRate (baseline 0.83 ± 0.06 vs 12 weeks after 0.56 ± 0.23) and a significant difference from the mice without SRS (P < .001) (Figure 1F). The IMT of the TST, an indicator of despair, was largely the same in the mice with and without SRS at baseline (126.6 ± 59.1 vs 131.7 ± 64.9 seconds). However, with the occurrence of SRS, the IMT of mice with SRS showed a gradually increasing trend, while the mice without SRS showed no significant trend at different stages (Figure 1G).

To further analyze the correlation between the development of depression-like behaviors and SRS in pilocarpine-induced chronic TLE mice, all the SucroRate and IMT results were pooled together to draw scatter plots, and mice with and without SRS could be clearly separated (Figure 2A). The results of the SucroRate and IMT were combined as the severity of depression-like behaviors and positively correlated with the frequency of SRS episodes (Spearman correlation, R²=0.84, P = .01; Figure 2B). These results provide evidence for the hypothesis that seizure activity in chronic TLE might contribute to the development of depression.

Alterations in HRV Parameters Correlate With Severity of Depression-Like Behaviors

For the analysis of HRV parameters in the different phases of pilocarpine-induced SE in mice with and without SRS, the differences in HRV parameters between baseline and 12 weeks after SE induction were calculated. We found an increase in heart rate (baseline 531.6 ± 38.5 bmp vs 12 weeks after 616.7 ± 100.8 bmp) and a decrease in time-domain parameters (SDNN baseline 10.1 ± 1.7 milliseconds vs 12 weeks after 7.7 ± 6.3 milliseconds; RMSSD baseline 5.9 ± 1.7 milliseconds vs 12 weeks after 5.3 ± 3.5 milliseconds) in mice with SRS. No significant difference in mice without SRS was found. Partial least squares discriminant analysis was employed to evaluate the characterization, and classification of the HRV parameters in mice with and without SRS and showed a different distribution between the 2 groups (Figure 2C). We calculated the difference in HRV parameters between baseline and 12 weeks after SE induction as the alteration of different phases and found that the degree of HRV alterations was positively correlated with the severity of depression in mice with SRS (Figure 2D). The results above suggest that validated seizure attacks lead to both alterations in HRV parameters and the development of depression simultaneously.

Basic Data Description for scRNA-seq

The pilocarpine-induced chronic TLE mice were divided into groups with and without depression according to the variation in depression-like behaviors in the progression of SRS. scRNA-seq was performed in glial cells and neurons from the TLE control-related HPC, depression control-related PFC and VTA, and HRV central control-related NTS (Figure 3A). To discover the altered regulation of gene expression, we combined HPC with NTS and PFC with VTA in with- and without-depression groups of mice and obtained 4 groups of samples for the analysis (Figure 3B). The harvested nerve cells from different groups were sequenced on a 10× Genomics platform. After application of quality control filters (Figure 3B), we obtained a total of 34 193 cells and 20 373 transcriptomes from the above 4 samples. The cell identities were inferred by reliably identifying specific nerve cell markers (supplementary Table 3) (Deczkowska et al., 2018; Morarach et al., 2021). We visually inspected the nerve cell type of computed clusters, projecting the data onto 2 dimensions by t-distributed stochastic neighbor embedding (t-SNE) and obtained 7 clusters of cell types, including neurons, microglia, oligodendrocytes, astrocytes, endothelial cells, mural cells, and oligodendrocyte precursor cells (Figure 3C). To characterize the different types of nerve cells in the above 4 samples, we applied violin plots to describe the expression of specific markers (Figure 3D), and each group demonstrated a similar color to the corresponding cell clusters in Figure 3C.

Characteristic Expression of Mitochondria-Related Genes Is Elevated in Glial Cells in TLE Mice With Depression

DEG analysis was performed in glial cells between chronic TLE mice with and without depression. Types of glial cells from related brain regions are identified by the expression of specific markers (Deczkowska et al., 2018). We found that in both the HPC+NTS and PFC+VTA brain areas, there was a large number of upregulated DEGs related to mitochondrial function (Figure 4A shows the PFC+VTA brain area). Furthermore, we analyzed the proportion of mitochondrial-related genes among glial cell DEGs and found that mitochondrial-related genes accounted for the majority of upregulated DEGs in glial cells (Figure 4B). The nonmitochondrial-related DEGs were irregularly scattered in other organelles (Figure 4C). We further analyzed the overlapping and uniquely expressed DEGs among microglia, oligodendrocytes, and astrocytes (Figure 4D) and found significantly uniquely expressed DEGs in microglia and astrocytes (Figure 4E).

Specific Upregulation of Inhibitory Neurons in NTS of TLE Mice With Depression

DEG analysis was performed in neurons between chronic TLE mice with and without depression. Excitatory and inhibitory neurons from related brain regions were identified by the expression of specific markers (Figure 5A) (Morarach et al., 2021). DEGs, such as Hsp90aa1, Pcdh11x, Vwa5a, and Ttr, which have been linked to depression in previous studies, have been found (Crow, 2011; Xiang et al., 2018; Szabo et al., 2021; Yu et al., 2021). In total, 5110 neurons in the above 4 groups were obtained and classified into 12 subtypes (Figure 5B) based on the expression of specific markers (supplementary Table 4). Furthermore, we analyzed the distribution of neuron subtypes in the with- and without-depression groups (Figure 5C) and found a specifically expressed subtype in the HPC+NTS area in the depression group (Figure 5C; correspondence Cluster 4 in Figure 5B). To localize and characterize the clusters obtained in Figure 5B, the expression of mRNAs marking glutamatergic (vesicular glutamate transporter [VGLUT1: slc17a7 and VGLUT2: slc17a6]) and GABAergic (vesicular GABA transporter [slc32a1]) (Deczkowska et al., 2018) (supplementary Table 5) was described (Figure 5D), and the specifically upregulated cluster 4 in Figure 5C was identified as a class of inhibitory neuron subtype from the NTS (Figure 5D).

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) Enrichment Analysis of DEGs

To further explore the biological functions of the DEGs in glial cells and neurons, GO and KEGG enrichment analyses were conducted. Biological process of GO enrichment analysis in glial cells showed that changed genes were mainly implicated in regulation of metabolic process (Figure 6A). More importantly, we found significantly enriched GABAergic synapse pathways in KEGG enrichment analysis in astrocytes and oligodendrocytes in the HPC+NTS brain area (Figure 6B) but not in microglia or the PFC+VTA brain area (supplementary Figure 1).

GO enrichment analysis was also performed in DEGs from excitatory and inhibitory neurons. The biological process of GO enrichment analysis in excitatory neuron DEGs showed that the genes were involved in the regulation of metabolic processes and developmental processes (supplementary Figure 2); however, the inhibitory neuron DEGs were mostly connected with the response to endogenous stimulus (Figure 6C). Furthermore, we found a significantly enriched long-term depression pathway in DEGs from inhibitory neurons in KEGG enrichment analysis (Figure 6D). No depression-related significantly enriched pathway was found in excitatory neurons (supplementary Table 6), suggesting that the development of depression in chronic TLE was mostly associated with inhibitory neurons.

DISCUSSION

In a recent clinical study, we found that autonomic dysfunctions reflected by decreased HRV were significant in epilepsy with depressive disorders (Aimaier et al., 2022). However, the exact relationship between diminished HRV and epilepsy–depression comorbidity is not yet clear. We have now extended the phenomenological analysis to animal behavior tests and scRNA-seq analysis, revealing that severe HRV and neuron subtype alterations occur in chronic TLE mice with depression.

Epilepsy can be interpreted as a disruption of the normal balance of excitation and inhibition in the brain, and this imbalance can be induced by changes in all aspects of brain function, ranging from genetic mutations and altered neuronal subtype signaling to abnormalities in a wide range of neural circuits (Stafstrom and Carmant, 2015). A bidirectional relationship between depressive disorders and epilepsy has been suggested by several population-based studies and is supported by animal model experiments (Adelöw et al., 2012; Hesdorffer et al., 2012) However, the underlying etiology remains elusive. Sympathetic nerve system–dominant autonomic dysfunctions have been found in epilepsy (Evrengül et al., 2005; Dütsch et al., 2006). HRV is thought to be regulated by the central autonomic network (Thayer and Lane, 2009; Sgoifo et al., 2015), and NTS inhibitory GABAergic neurons play a crucial role in the regulation of HRV (Sgoifo et al., 2015). HRV parameters can be used to noninvasively assess the dynamic homeostasis of sympathetic and parasympathetic nerves (Sgoifo et al., 2015). Previous studies of depression have also found that the decrease in HRV was correlated with depression (Sgoifo et al., 2015; Borrione et al., 2018; Zhang et al., 2021) and suggested that the ANS is 1 of the main neural pathways activated by chronic endogenous and exogenous stress. In conditions associated with chronic stress, such as depression and epilepsy, the sympathetic nervous system can be continuously activated while parasympathetic output is suppressed, and impaired autonomic homeostasis can be reflected by decreased HRV (Won and Kim, 2016). Sustained sympathetic nerve activation is associated with cardiovascular diseases (Julius, 1990; Esler and Kaye, 2000), metabolic diseases (Seravalle and Grassi, 2016), and mood disorders (Won and Kim, 2016; Borrione et al., 2018). The pilocarpine-induced mouse model of chronic epilepsy is characterized by epilepsy-induced psychotic behavior and is recommended as an animal model for studying psychiatric comorbidity in epilepsy (Müller et al., 2009; Wulsin et al., 2018). Therefore, the dynamic changes in depression-like behaviors and HRV parameters in chronic TLE mice at different phases were evaluated in this study. Significant alterations in HRV parameters in chronic TLE mice were found, suggesting a shift to sympathetic output dominance in the ANS. The alterations were positively correlated with the severity of depression, suggesting a strong relationship between depressive disorders and sympathetic nerve dominant autonomic dysfunctions in epilepsy.

scRNA-seq analysis further revealed the mechanism of depressive disorders and autonomic dysfunctions in epilepsy. scRNA-seq was performed in epilepsy-, depression-, and HRV control–related areas in chronic TLE mice with and without depression. We found extensively upregulated mitochondria-related genes in glial cells in the depression group, and the DEGs in glial cells were mostly associated with the regulation of metabolic processes. The function of mitochondria and disrupted homeostasis of metabolism have been considered to play an important role in the development of depression (Rezin et al., 2009; Bansal and Kuhad, 2016; Feuerstein et al., 2016; Petschner et al., 2018).

More importantly, we found upregulated inhibitory neurons in the NTS in the depression group compared with the without-depression group. The findings of significant enrichment of the GABAergic synapse pathway in glial cells from the HPC+NTS brain area further confirmed the result of specifically upregulated inhibitory neurons in the NTS. GABAergic inhibitory neurons in the NTS have been found to regulate the output of HRV (Barbas et al., 2003; Saha, 2005; Thayer and Lane, 2009; Sgoifo et al., 2015). In addition, long-term activity of GABAergic neurons emerged in the NTS in association with epileptogenesis, which might contribute to an increased risk of cardiorespiratory autonomic dysfunction and sudden death in pilocarpine-induced TLE (Derera et al., 2017, 2019). The pathway enrichment analysis in inhibitory neurons showed that the DEGs were implicated in long-term depression, suggesting a possible relationship of inhibitory neurons with depression in chronic TLE. Our findings in this study further clarify the changes in neurons from the NTS in the TLE mouse model and suggest that the inhibitory neurons of the NTS might be the common regulatory target of HRV and depressive disorders in epilepsy.

The limitations of the study are as follows. First, the current research results cannot prove the exact relationship between the changes in NTS inhibitory neuron activity and the alterations in HRV and depression-like behaviors in epilepsy. To further confirm our hypothesis, calcium imaging and photogenetic regulation need to be performed to verify the activity of inhibitory neurons in different phases of epilepsy and their specific role in epileptic autonomic nerves and depressive behaviors. Second, pilocarpine-induced TLE mice served as an animal model for TLE in humans, in which genetic background and lifestyle interplay play a pivotal role. It remains a hurdle to accurately distinguish whether the transcriptomic differences observed in our study are caused by the genetic background. In our study, we could also not distinguish the causal changes from the downstream effects. Third, we mainly presented descriptive results in the scRNA-seq analysis. As the main purpose of the study was to serve as a resource for the characterization of cell type-specific changes in chronic TLE mice with depression, only the comparisons of neurons and glial cells in 4 related brain areas in TLE mice with and without depression were included in the study. Building on the present results, future work might include functional studies of potentially interesting molecules. Finally, only male mice were included in our study, which might differ from female mice. To determine the difference in gene expression profiling in different sexes, future experiments are needed.

CONCLUSIONS

In summary, autonomic dysfunctions represented by decreased HRV parameters were positively correlated with the severity of depression in TLE. Significantly upregulated mitochondria-related genes in glial cells in the depression group were found, which have been considered to play a pivotal role in the development of depression. More importantly, the HRV central control-related neuron subtype (inhibitory neurons of the NTS) was specifically expressed in the depression group compared with the without-depression group. The enrichment of specifically expressed inhibitory neurons was associated with the long-term depression pathway. These findings suggest that specifically expressed inhibitory neurons in the NTS are involved in development of depression and decrease in HRV in TLE. Further mechanistic studies are required to substantiate the pathological importance of our findings by exploring the activity of NTS inhibitory neurons in different phases of epilepsy.

Acknowledgments

This work was supported by a grant from the National Natural Science Foundation of China (grant no. 82171444).

Author Contributions

Guliqiemu Aimaier: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Writing—original draft, Software. Kun Qian: Conceptualization, Data curation, Investigation, Methodology, Resources, Project administration. Huateng Cao: Data curation, Formal analysis, Software. Weifeng Peng: Conceptualization and methodology. Zhe Zhang: Investigation, Methodology, Data curation, Software. Jianhua Ma: Investigation, Methodology, Data curation. Jing Ding: Conceptualization, Methodology, Validation, Resources, Writing—review and editing, Supervision, Project administration. Xin Wang: Conceptualization, Resources, Supervision, Project administration, Funding acquisition.

Ethics Statement

All experimental procedures were approved by the Animal Care Committee of Shanghai Institutes for Biological Sciences, Chinese Academy of Sciences (Shanghai, China).

Interest Statement

The authors declare no conflict of interest.

Data Availability

The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

References

Author notes