Abstract
Alcohol withdrawal syndrome (AWS) is characterized by neuronal hyperexcitability, autonomic dysregulation, and severe negative emotion. The nucleus tractus solitarius (NTS) likely plays a prominent role in the neurological processes underlying these symptoms as it is the main viscerosensory nucleus in the brain. The NTS receives visceral interoceptive inputs, influences autonomic outputs, and has strong connections to the limbic system and hypothalamic-pituitary-adrenal axis to maintain homeostasis. Our prior analysis of single neuronal gene expression data from the NTS shows that neurons exist in heterogeneous transcriptional states that form distinct functional subphenotypes. Our working model conjectures that the allostasis secondary to alcohol dependence causes peripheral and central biological network decompensation in acute abstinence resulting in neurovisceral feedback to the NTS that substantially contributes to the observed AWS. We collected single noradrenergic and glucagon-like peptide-1 (GLP-1) neurons and microglia from rat NTS and measured a subset of their transcriptome as pooled samples in an alcohol withdrawal time series. Inflammatory subphenotypes predominate at certain time points, and GLP-1 subphenotypes demonstrated hyperexcitability post-withdrawal. We hypothesize such inflammatory and anxiogenic signaling contributes to alcohol dependence via negative reinforcement. Targets to mitigate such dysregulation and treat dependence can be identified from this dataset.
Introduction
Alcohol withdrawal syndrome (AWS) is characterized by adverse physical and emotional symptoms. Physical symptoms are driven by autonomic dysregulation, γ-aminobutyric acid (GABA) hypoactivity, and increased glutamatergic signaling leading to dysphoria, nausea, diaphoresis, tachycardia, hypertension, seizures, and delirium tremens (). Fear and anxiety are the principle emotional symptoms. The negative reinforcement hypothesis of substance dependence postulates that these negative physical and emotional symptoms experienced in withdrawal motivate alcohol dependence (, ; ; ). We conjecture that peripheral network decompensation is a central facet of this model, and that neurovisceral feedback via the vagus nerve conveying peripheral information to the central nervous system contributes substantially to the severity of symptoms experienced () (Figures 1A,B). Investigation into the underlying mechanisms producing these symptoms may provide insight into targets that mitigate acute and protracted AWS severity and prevent relapse following abstinence. Such treatments may provide clinical utility for other substances of abuse with severe withdrawal syndromes such as opioids.
FIGURE 1
Neuroinflammatory processes have emerged as an important contributor to the severity of AWS symptoms, especially in the amygdala (
The nucleus tractus solitarius (NTS) is another brain region that contributes to the symptoms of AWS and is implicated in alcohol dependence (
Indeed, NE+ and GLP-1+ neurons in the NTS have been shown to contribute to symptoms of AWS and alcohol intake in withdrawal (
Here, we measured how the functional states of single-neuron samples containing neuronal phenotypes enriched with NE neurons or GLP-1 neurons and microglia in the NTS change over the course of alcohol withdrawal. Single-cell approaches allow for the identification of cellular subphenotypes—morphologically indistinguishable cells anatomically localized that use the same primary neurotransmitter yet have distinct transcriptomic profiles. Our previous work has demonstrated the heterogeneity of single-cells, and the functional importance of subphenotypes that may be missed in tissue-level approaches (
Results
We used laser capture microdissection (LCM) to gather single cells from rat NTS in control, chronic ethanol (EtOH), 8- h wd, 32 h wd, or 176 h wd treatments (Supplementary Figure 1A) 0.10 cells were pooled to comprise a sample that underwent microfluidic reverse transcription quantitative polymerase chain reaction (RT-qPCR) to measure a subset of the transcriptome in these samples comprising 10 individually selected single-cells (Supplementary Figure 1B) (
The expression of the neurotransmitter precursor genes Th and Gcg for NE and GLP-1 neuron enriched samples, respectively, is plotted across treatment time points in Figures 2A,B. We find that expression of these transcripts is inversely correlated—at time points in which Th expression is relatively high in NE neuron enriched samples, Gcg expression is relatively low in GLP-1 neuron enriched samples, and vice versa (Figures 2A,B). We speculate this may indicate a push-pull dynamic mechanism in the genetic regulation of neurotransmission by these neurons. Notably, Gcg expression was induced in the three withdrawal time points though moderately at the 32 h wd time point, in which Th expression was induced, consistent with the aforementioned push-pull dynamic. Additionally, Gcg was low in control and EtOH treatments which may suggest that GLP-1 neurotransmission is pathologically elevated during the withdrawal process. These hypothesis-generating observations require verification. At the 176 h wd time point, bimodal distribution of Th expression in NE neuron enriched samples and trimodal Gcg expression in GLP-1 neuron enriched samples is observed (Figures 2A,B). High and low Th-expressing NE neuron enriched samples and Gcg-expressing GLP-1 neuron enriched samples from this time point were separated into heat maps organized by Euclidian distance clustering of gene expression (Figures 2C,D). Th expression in NE neurons and Gcg expression in GLP-1 neurons is moderately predictive of cellular subphenotypes that loosely organize co-expression gene clusters. However, Th-expression and Gcg-expression alone, which is to say neurotransmitter expression, is not the best determinate of cellular subphenotypes—a finding we have observed previously in other neuronal nuclei (
FIGURE 2

Neurotransmitter Expression in Neurons. **p < 0.003 *p < 0.05, #p < 0.1; two-tailed heteroscedastic t-test. Bars show standard error. (A) Tyrosine hydroxylase expression (Th) in norepinephrine (NE) neuron enriched samples. Bimodal Th expression at the 176 h wd time point is explored in the heat map in panel (C). (B) Preproglucagon (Gcg) expression in GLP-1 neuron enriched samples. Trimodal Gcg expression at the 176 h wd time point is explored in the heat map in panel (D). (C) Heat map exploring high Th-expressing neurons and low Th-expressing neurons as subclusters of NE neuron enriched samples in the 176 h wd treatment. Single-cell gene expression is shown as z-scores of –ΔΔCt values. Neurotransmitter expression levels does not determine high-template co-expression gene clusters. (D) Heat map exploring high, middle, and low Gcg-expressing neurons as subclusters of GLP-1 enriched neuronal samples in the 176 h wd treatment. Single-cell gene expression is shown as z-scores of –ΔΔCt values. Neurotransmitter expression levels does is more predictive of co-expression gene clusters than Th.
Heat maps establishing well-defined data-driven cellular subphenotypes in neuronal samples and microglia were generated using Euclidean distance clustering of z-scores of the –ΔΔCt values for each sample and gene in the dataset (Figures 3–5). Co-expression gene clusters are labeled with numbers and cellular subphenotype groupings are labeled with letters. GLP-1 neuron enriched samples had the same cellular subphenotypes with the same gene clusters across all treatments while NE neuron enriched samples and microglia had two co-expression configurations comprising the identified subphenotypes. The proportion of cells constituting a subphenotype in addition to gene cluster expression levels shifted with the treatment. GABAAR subunit genes clustered together in every configuration, and their expression was largely indicative of subphenotype groupings.
FIGURE 3

Heat Map of NE Neurons. Heat map displays cellular subphenotypes within NE neurons enriched samples through alcohol withdrawal time series. Rows represent 10-cell pooled samples with cellular subphenotype clusters labeled with uppercase letters. Columns represent the z-score of –ΔΔCt gene expression values on a –1 to 1 color scale for that gene in that sample. Gene clusters are labeled by number.
FIGURE 4

Heat Map of GLP-1 Neurons. Heat map displays cellular subphenotypes within GLP-1 neuron enriched samples through alcohol withdrawal time series. Rows represent 10-cell pooled samples with cellular subphenotype clusters labeled with uppercase letters. Columns represent the z-score of –ΔΔCt gene expression values on a –1 to 1 color scale for that gene in that sample. Gene clusters are labeled by number.
FIGURE 5

Heat Map of Microglia. Heat map displays microglia cellular subphenotypes through alcohol withdrawal time series. Rows represent 10-cell pooled samples with cellular subphenotype clusters labeled with uppercase letters. Columns represent the z-score of –ΔΔCt gene expression values on a –1 to 1 color scale for that gene in that sample. Gene clusters are labeled by number.
Two prominent subphenotypes, C and D, emerged in NE neuron enriched samples in withdrawal time points. Subphenotype C highly expressed gene cluster 4 which is rich in inflammatory ligands and receptors including Crh, Il1b, and Ptgs2. This was labeled the “inflammatory” gene cluster (see below). Subphenotype C suppressed gene cluster 5, which includes Th and GABAAR subunits, and 6. Subphenotype D, conversely, had the opposite expression pattern in these gene clusters. At 8 h wd, subphenotype C was predominant, but subphenotype D made up a higher proportion of cells at 32 and 176 h wd. Another subphenotype, E, emerged in the 32 h wd treatment that had moderately high expression of gene clusters 4 and 5. Subphenotypes D and E were further split into D1, D2, E1, and E2 based on medium or high expression, respectively, of gene cluster 6 which includes Cd200, cFos, and Mif.
Noradrenergic (NE) neuron enriched samples from control and EtOH treatments shared gene clusters distinct from the withdrawal treatments. GABAAR subunits were co-expressed consistent with all sampled cell types. GABAAR subunit genes showed high expression in subphenotype B in control and EtOH treatments, moderate expression in subphenotype D and E in 8 and 32 h wd treatments, and returned to high expression in 176 h wd subphenotype D. In subphenotype E, only found at the 32 h wd time point, all assayed genes were at least moderately expressed which may suggest that the regulatory mechanisms of gene expression that control this transcriptomic profile are activated at this phase of the withdrawal process. This observation requires further mechanistic study.
Co-expression genes clusters in GLP-1 neuron enriched samples were consistent throughout the time series. Subphenotype A highly expressed “inflammatory” gene cluster 1 rich in cytokine and chemokine ligands and receptors including Crh, Il1b, and Ptgs2, while suppressing “GABAAR” gene cluster 2. Subphenotype B had the opposite pattern. GLP-1 co-expression clusters 1 and 2 were surprisingly similar to NE co-expression clusters 4 and 5, again suggesting the mechanisms of regulatory constraint are shared between these phenotypes. Interestingly, GLP-1 subphenotype B emerged only in the EtOH treatment. At 8 h wd, subphenotypes A and B suppressed expression of their high-expressing gene clusters, 1 and 2, respectively, compared to control. At 32 and 176 h wd, subphenotype A gene cluster 1 was more highly expressed than in the control condition. Subphenotype B gene demonstrated a steady decrease in expression of gene cluster 2, and by the 176 h wd time had only moderate expression of “GABAAR” gene cluster 2. Concurrently, expression of gene cluster 3 for subphenotype B consistently increases throughout the withdrawal process and by 176 h wd is the most prominently upregulated gene cluster. Tnf did not group into any gene cluster and is isolated in the GLP-1 neuron enriched samples heat map to display this clearly (Figure 4).
Microglia shared co-expression clusters in subphenotypes A, B, and C for control, 32 h wd, and 176 h wd treatments. EtOH and 8 h wd treatments shared co-expression clusters in subphenotypes D and E. Subphenotype A was the exclusive expression pattern for the 32 h wd treatment which upregulated “inflammatory” gene cluster 1 including Crh, Il1b, and Ptgs2. This subphenotype A is most similar to the so-called M1 phenotype (
Next, single-cell samples were combined within their subphenotypes to yield an average value of expression for each gene within that subphenotype. This data is displayed in the cellular diagrams of Figures 6–8. Each box represents an assayed gene. Its color indicates the average z-score of –ΔΔCt expression values. The location of this box in the cellular cartoon corresponds to the protein function of that gene. These diagrams provide a higher-level display of the functional state of the subphenotype at that time point and display the transcriptional dynamics occurring in each subphenotype in a readable way. In brief, Figure 6 shows a clear upregulation of GABAAR genes at 176 h wd as compared to 8 h wd in Group C and D of NE neuron enriched samples. Additionally, CD200 expression is one of the primary distinguishers of Group D1 vs. Group D2. Figure 7 displaying GLP-1 neuron enriched samples shows that GABAAR gene expression at 176 h wd is decreased in both subphenotypes. At the 32 h wd time point, Cxcl10, Cxcr1, Cxcr2, and Cxcr3 expression distinguish Group A1 most prominently from Group A2. Microglia displayed in Figure 8 showed the most Tnf expression in groups D and E, subphenotypes only identified at the EtOH and 8 h wd treatments. Cx3cr1, a mircoglia gene prominently involved in neuronal adhesion, showed increased expression at the 176 h wd time point in all subphenotypes (Wolf et al., 2013). GLP-1 neuron enriched samples Group B had the most increased Cx3cl1 expression possibly indicated this subphenotype interacts most with microglia at this time point.
FIGURE 6

Subphenotype Gene Expression in Norepinephrine Neuron Enriched Samples. Cellular diagrams display boxes representing relative gene expression (average z-score of –ΔΔCt values) of subphenotypes shown in prior heatmaps. Legend with gray boxes in lower right labels which boxes correspond to which gene. Box color represents expression (blue is low expression and yellow is high expression). The location of the box represents the localization or function of the protein product from that gene transcript. Green numbers indicate subgroups within subphenotypes. Groups A and B shown in Supplementary Figure 5.
FIGURE 7

Subphenotype Gene Expression in GLP-1 Neuron Enriched Samples. Cellular diagrams display boxes representing relative gene expression (average z-score of –ΔΔCt values) of subphenotypes shown in prior heatmaps. Legend with gray boxes in on right of figure labels which boxes correspond to which gene. Box color represents expression (blue is low expression and yellow is high expression). The location of the box represents the localization or function of the protein product from that gene transcript. Green numbers indicate subgroups within subphenotypes. Group C shown in Supplementary Figure 5.
FIGURE 8

Subphenotype Gene Expression in Norepinephrine Neuron Enriched Samples. Cellular diagrams display boxes representing relative gene expression (average z-score of –ΔΔCt values) of subphenotypes shown in prior heatmaps. Legend with gray boxes in center of figure labels which boxes correspond to which gene. Box color represents expression (blue is low expression and yellow is high expression). The location of the box represents the localization or function of the protein product from that gene transcript.
Discussion
Nucleus tractus solitarius (NTS) neurons regulate emotion, autonomic homeostasis, and stress responses. Multiple neuronal nuclei, ligands, receptors, and signaling dynamics are involved in these complex functions including NE, GLP-1, CRH, and GABA (Figures 1A,B). Moreover, local glial-neuronal paracrine signaling via inflammatory cytokines like tumor necrosis factor-alpha (TNF-α) also play a role. We microdissected single Th + neurons, Th- neurons, and microglia from the rat NTS as 10-cell pooled samples using LCM and measured their expression of 96 gene transcripts in an alcohol withdrawal time series (Supplementary Figure 1). Time points were chosen based on rat alcohol metabolism and withdrawal symptomatology (
Further analysis of Th and Gcg expression showed an inverse relationship with respect to time point with Gcg expression demonstrating elevated expression levels only during withdrawal (Figures 2A,B). However, expression of these neurotransmitter precursor genes did not organize the other genes assayed into distinct subphenotypes correlated to their expression levels (Figures 2C,D). A data-driven approach to cellular subphenotype organization identified stark subphenotypes unique to each cell type likely with discrete functions (Figures 3–5). Strikingly, these subphenotypes shared similarities in their expression of their inflammatory gene clusters (Supplementary Table 4). Single-cell pooling may contribute to the observed heterogeneity of transcriptomic subphenotypes though single-neuron datasets also demonstrate high heterogeneity (
Gene cluster 4 in NE neurons enriched samples, gene cluster 1 in GLP-1 neuron enriched samples, and gene cluster 1 in microglia samples constituted these “inflammatory” clusters (Supplementary Table 4). 18 genes were shared across all of these co-expression clusters and only 5 genes were unique to a single cluster suggesting similar mechanisms across cell types that regulate their expression. In NE enriched neuronal samples, subphenotype C highly expressed this inflammatory cluster while subphenotype E had moderate inflammatory co-expression cluster elevation. At 8 h wd, NE subphenotype C was 62.5% of the samples (5/8) and at 32 h wd C and E combined to 62.2% of the samples (23/37). By 176 h wd, subphenotype C was only 29% of NE neuron enriched samples (5/17). This may suggest that this subphenotype of NE neurons experiences a marked increase in inflammation during acute AWS, but that this subphenotype is not involved in protracted withdrawal symptoms such as low-grade anxiety (
In GLP-1 neuron enriched samples, subphenotype A highly expressed the “inflammatory” gene cluster (gene cluster 1). The pattern of expression in this inflammatory subphenotype of GLP-1 neuron enriched samples (A) is similar to the inflammatory subphenotypes of NE neuron enriched samples (C and E). In samples enriched with GLP-1 neurons, subphenotype A makes up 33.3% of control samples (3/9), 0% of EtOH samples, 62.5% (5/8) of 8 h wd samples, and 55.2% (16/29) of 32 h wd samples. By 176 h wd, this inflammatory subphenotype has decreased back near control levels: 35.7% (5/14) of GLP-1 neurons.
Surprisingly, microglia demonstrated a similar pattern. Microglia subphenotype A also highly expressed the inflammatory gene cluster (cluster 1). High gene expression in this cluster is indicative of M1 microglia phenotypes as this cluster includes the M1 markers Il1β, Il6, Nos1, Ptgs2, and TLRs 1,4, and 5 (
The 176 h wd time point is meant to measure long term changes in gene expression that occur in protracted withdrawal. At this time point, some similarities across cell types were observed in the subphenotypes that highly expressed GABAR subunits as was observed at other time points. NE neuron enriched sample cluster 5, GLP-1 neuron enriched sample cluster 2 and microglia cluster 2 contained the majority of the GABAR subunit genes, and the makeup of this “GABAR” co-expression cluster was not as consistent as the inflammatory cluster across cell types—16 genes are shared across all cell types and 16 genes are unique to a single cell type within its respective GABAR cluster (Supplementary Table 4). GLP-1 neuron enriched samples in subphenotype B upregulates this co-expression cluster in the control treatment, but the relative level of expression of this cluster decreases throughout the time series within this subphenotype (Figure 4). At the 176 h wd time point, this GABAR cluster is only moderately expressed which may suggest long term changes to this neuronal subphenotype following alcohol dependence and withdrawal. The decrease in expression of inhibitory GABAR gene transcripts, along with the concurrent upregulation of co-expression cluster 3, which contains Gcg, may suggest that this GLP-1 enriched sample subphenotype increases its GLP-1 neurotransmission in protracted withdrawal. Literature indicates that GLP-1 signaling from the NTS to the amygdala and other nuclei is anxiogenic (
Microglia also showed elevated GABAR expression at the 176 h wd time point, but the pattern of increased GABAR expression was unexpected. Control microglia in subphenotype C show moderate expression of both cluster 1 (inflammatory) and cluster 2 (GABAAR) (Figure 5). Expression of both clusters increase at the 176 h wd time point. This may suggest elevated inflammation, but not by distinct M1 phenotype microglia (subphenotype A), and also elevated GABAR expression. These observations are best visualized in Figure 8. Of note, there are many genes in microglia cluster 2 that are not GABAR subunits. Moreover, microglial Tnf expression was significantly elevated in control, EtOH and 8 h wd treatments compared to 176 h wd independent of subphenotype (Supplementary Table 3). Indeed, Tnf expression by microglia did not fit neatly into a gene cluster. Cluster C has some cells that demonstrate high Tnf expression in both control and 176 h wd, where Cluster A showed a decrease in Tnf expression between these two time points. Cluster B, conversely, increased its expression of Tnf from control to 176 h wd. This apparent absence of a pattern in microglia Tnf expression suggests that in microglia this gene that is central to neuroinflammation is constrained by a mechanism that is independent of the other genes measured in this study. Further, the decrease in overall microglia Tnf expression at 176 h wd as measured by an average of –ΔΔCt values and two-tailed heteroscedastic t-tests may be misleading. A single-cell analysis reveals that overall expression may not be the best indicator of inflammation. Rather, shifts in subphenotype proportion, and the number of cells showing a moderately increased Tnf expression, as seen in subphenotype C, may have more of a physiologic impact than total gene expression levels.
Cell diagrams in Figures 6–8 average the expression of a gene within a subphenotype designated by color and display that color in a location on the diagram that corresponds to the protein function. This method of data presentation allows for analysis of receptor-ligand interactions within and between subphenotypes. For example, Figure 6 displaying NE neuron enriched samples shows that at 32 h wd, subphenotype C experiences an increase in expression of ligand-receptor pair Ccl-Ccr and Cxcl10-Cxcr. This may indicate that CCL-CCR and CXCL10-CXCR signaling is elevated at this timepoint in AWS. Figure 6 also provides clarity in subphenotype D upregulation of Mapk1 at 176 h wd which may suggest long term transcription is altered during protracted withdrawal in this subphenotype. Moreover, transcription factor genes cFos, Junb, NfkB, and Stat3 have increased expression in subphenotype D2 which is consistent with this subset of NE neuron enriched samples experiencing long-term changes in transcription following alcohol withdrawal. Microglia in subphenotype C upregulate IL1a, IL1b, and IL1r1 at 176 h wd in subphenotype C as compared to control, while subphenotype B downregulate these genes at 176 h wd compared to control (Figure 8). This dynamic may suggest that subphenotype B provides an anti-inflammatory function that is most active in protracted withdrawal. Similarly, it may suggest that microglia subphenotype C, identified here as a microglia subset that can function in a multitude of processes whether inflammatory or anti-inflammatory based on their lack of a clear co-expression module pattern in control, is pushed toward an inflammatory state in protracted withdrawal.
This dataset has allowed the identification of cellular subphenotypes and their gene expression dynamics in alcohol withdrawal through time. The fusion of single cells into 10-cell pools can result in some obfuscation of phenotypic dynamics; The details of which can be resolved at a higher resolution. Nevertheless, this analysis has revealed valuable observations in both neurotransmission signaling and local paracrine signaling processes that aid in hypothesis-generation while relating to what is observed clinically in the context of what is already established about such neurotransmission. The dataset is unique in that microfluid RT-qPCR, a method lower in throughput but more reliable than RNA-seq (
The major weakness of this study is the number of animals assayed. Ten rats total were assayed and single cells were collected from a single animal from some conditions (control neurons, chronic ethanol, 8 h wd, 32 h wd microglia). This design was due to both cost, and previous studies from our group that consistently demonstrate that single-cells within an animal have as much transcriptional heterogeneity, or variance, as between animals (
We have collected the data, validated the accuracy of the dataset, and identified cellular subphenotypes and their major signaling dynamics. However, signaling dynamics measured in our dataset can be further investigated and may identify clinical targets to treat acute or protracted AWS and potentially alcohol dependence itself. Future studies analyzing these signaling dynamics with the addition of female rats that also include other brain cell types such as astrocytes and endothelial cells are needed to further understand the underlying pathophysiology of AWS and dependence.
Lastly, these findings are consistent with our hypothesis that neuroinflammation in the visceral-emotional neuraxis contributes to antireward which motivates alcohol, and opioid, dependence (Figures 1A,B and Supplementary Figure 7) (
Materials and Methods
Animals
Approval of protocols was given by Institutional Animal Care and Use Committee of Thomas Jefferson University. The study was carried out in compliance with ARRIVE guidelines and in accordance with all relevant guidelines and regulations. Ten young, male, Sprague Dawley rats (35–45 grams) ordered from Harlan Laboratory were housed individually in the Thomas Jefferson University Alcohol Research Center Animal Core Facility. Standard chow and water were given until rats weighed 120 grams. Rats were then fed an alcohol-free, maltose-dextrin substituted, Lieber-DeCarli liquid diet (bioServe, Frenchtown, NJ) for three days (
Rapid Decapitation, Fast Staining Protocol, Laser Capture Microdissection
Dissected brainstems were frozen in Optimal Cutting Temperature (O.C.T.) following rapid decapitation for cryostat sectioning and stored at -80°C for nucleic acid preservation. An in-house rapid immunofluorescent staining protocol developed to preserve nucleic acid integrity was used to visualize cell types for single-cell LCM as explained elsewhere (Supplementary Figure 1B) (
Single Cell Sampling and High-Throughput RT-qPCR
3230 single brain cells, 950 Th + neurons, 1030 Th− neurons, and 1250 microglia, were collected from the NTS using LCM. Cells were grouped into 10-cell pools comprising 323 total samples analyzed. This pooling of cells increases the number of samples analyzed by the microfluidic RT-qPCR platform. cDNA from mRNA transcripts was generated by reverse transcription (SuperScriptTM VILOTM cDNA Synthesis Kit; ThermoFisher). TaqMan PreAmp Master Mix was used for pre-amplification of cDNA (22 cycles) with forward and reverse PCR primers (96 pairs). The Biomark microfluidic qPCR platform (Fluidigm©) was used to measure expression levels of 96 genes. Four batches of probe-based qPCR measured the previously amplified 96 cDNA transcripts. Supplementary Table 1 lists primers used. Primer amplicon validation was performed on agarose gel electrophoresis. Following strict quality control protocols, a total of 229 10-cell pooled samples (70 NE neuron samples, 65 GLP-1 neuron samples, and 94 microglial samples) and 65 gene transcripts were used for data analysis.
The four microfluidic RT-qPCR batches run for this study were assessed for intra- and inter-batch experimental quality (Supplementary Figures 3, 4). Technical replicates assessing intra-batch quality demonstrated high similarity with r values listed (Supplementary Figure 3). Inter-batch replicates demonstrated high batch similarity, though batch 4 sample 40 showed contamination (Supplementary Figure 4). A dilution series using standard rat brain RNA was also included in each batch for quantitative analysis (Supplementary Table 2); However, the data normalization method explained below calculated relative expression and was used for all analysis in this study.
Data Normalization
A two-step median-centering –ΔΔCt method was used for expression level normalization was explained elsewhere (
Publisher’s Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Statements
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.
Ethics statement
The animal study was reviewed and approved by Institutional Animal Care and Use Committee of Thomas Jefferson University.
Author contributions
SO’S performed microfluidic qPCR, data analysis, figure generation, and writing of manuscript. DM-C collected single-cell samples under the guidance of JP. JP also designed the experiments. RV and JS were involved with figure design and editing. All authors discussed the results and commented on the manuscript.
Funding
The work presented here was funded through NIH HLB U01 HL133360 awarded to JS and RV and NIDA R21 DA036372 awarded to JS and Elisabeth Van Bockstaele and T32 AA-007463 awarded to Jan Hoek and supporting SO’S.
Acknowledgments
SO’S would like to acknowledge Jan Hoek for his support with the T32 AA 007463.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnsys.2021.739790/full#supplementary-material
Supplementary Figure 1Experimental design and single-cell selection. (A) Rat triplets were randomly assigned one of five treatments studied. (B) Single-cell selection and measure of transcription.
Supplementary Figure 2Density plots of normalized data across all genes. (A) Three separate rats provided single-cell samples for 32 h wd neurons. Gene expression variance between animals was similar or less than gene expression variance within an animal. (All 32 h wd neurons N = 3093, variance = 4.16; Animal 6 N = 1743, variance = 2.88; Animal 9 N = 730, variance = 5.76; Animal 10 N = 620, variance = 5.58). (B) (All 32 h wd NE neurons N = 1815, variance = 3.96; Animal 6 N = 1112, variance = 3.04; Animal 9 N = 391, variance = 4.24, Animal 10 N = 312, variance = 6.63). (C) (All 32 h wd GLP-1 neurons N = 1278, variance = 4.44; Animal 6 N = 631, variance = 2.60, Animal 9 N = 339, variance = 7.39, Animal 10 N = 308, variance = 4.52). (D) Two separate rats provided single-cell samples for control microglia. Gene expression variance between animals was similar or less than gene expression variance within an animal (All Control Microglia N = 1451, variance = 3.97; Animal 2 N = 958, variance = 3.40; Animal = 4 N = 493, variance = 5.01).
Supplementary Figure 3Technical replicate plots of raw Ct values.
Supplementary Figure 4Inter-batch replicate plots of raw Ct values.
Supplementary Figure 5Inter-batch replicate plots of raw Ct values.
Supplementary Figure 6Suphenotype gene expression in NE and GLP-1 neurons. Cellular cartoons display boxes representing relative gene expression (average z-score of −ΔΔCt values) of subphenotypes shown in prior heatmaps. Legend on right labels which boxes correspond to which gene and the color that represents expression (blue is low expression and yellow is high expression). The location of the box represents the localization or function of the protein product from that gene transcript. Legend is shown in gray boxes with labels.
Supplementary Figure 7Schematic of opponent-process model of addiction. Originally published in
References
1
AchantaS.VermaA.SrivastavaA.NilakantanH.HoekJ. B.VadigepalliR.et al (2018). Single cell gene expression analysis identifies chronic alcohol-mediated shift in hepatocyte molecular states after partial hepatectomy.Gene Expr.1997–119. 10.3727/105221618X15361728786767
2
AiminoM.CokerC.SilbermanY. (2018). Acute ethanol modulation of neurocircuit function in the nucleus of the tractus solitarius.Brain Res. Bull.1385–11. 10.1016/j.brainresbull.2017.07.019
3
ArmstrongD. M.RossC. A.PickelV. M.JohT. H.ReisD. J. (1982). Distribution of dopamine-, noradrenaline-, and adrenaline-containing cell bodies in the rat medulla oblongata: demonstrated by the immunocytochemical localization of catecholamine biosynthetic enzymes.J. Comp. Neurol.212173–187. 10.1002/cne.902120207
4
BakerT. B.PiperM. E.McCarthyD. E.MajeskieM. R.FioreM. C. (2004). Addiction Motivation Reformulated: an Affective Processing Model of Negative Reinforcement.Psychol. Rev.11133–51. 10.1037/0033-295X.111.1.33
5
BakkenT. E.HodgeR. D.MillerJ. A.YaoZ.NguyenT. N.AevermannB.et al (2018). Single-nucleus and single-cell transcriptomes compared in matched cortical cell types.PLoS One13:e0209648. 10.1371/journal.pone.0209648
6
BärK. J.BoettgerM. K.BoettgerS.GrotelüschenM.NeubauerR.JochumT.et al (2006). Reduced baroreflex sensitivity in acute alcohol withdrawal syndrome and in abstained alcoholics.Drug Alcohol. Depend.8566–74. 10.1016/j.drugalcdep.2006.03.014
7
BreeseG. R.KnappD. J. (2016). Persistent adaptation by chronic alcohol is facilitated by neuroimmune activation linked to stress and CRF.Alcohol529–23. 10.1016/j.alcohol.2016.01.005
8
CarbiaC.LannoyS.MaurageP.López-CanedaE.O’RiordanK. J.DinanT. G.et al (2021). A biological framework for emotional dysregulation in alcohol misuse: from gut to brain.Mol. Psychiatry261098–1118. 10.1038/s41380-020-00970-6
9
CovarrubiasM. Y.KhanR. L.VadigepalliR.HoekJ. B.SchwaberJ. S. (2005). Chronic alcohol exposure alters transcription broadly in a key integrative brain nucleus for homeostasis: the nucleus tractus solitarius.Physiol. Genomics2445–58. 10.1152/physiolgenomics.00184.2005
10
de TimaryP.LeclercqS.StärkelP.DelzenneN. (2015). A dysbiotic subpopulation of alcohol-dependent subjects.Gut Microbes6388–391. 10.1080/19490976.2015.1107696
11
DeCiccoD.ZhuH.BrureauA.SchwaberJ. S.VadigepalliR. (2015). MicroRNA network changes in the brain stem underlie the development of hypertension.Physiol. Genomics47388–399. 10.1152/physiolgenomics.00047.2015
12
FerryB.RoozendaalB.McGaughJ. L. (1999). Role of norepinephrine in mediating stress hormone regulation of long-term memory storage: a critical involvement of the amygdala.Biol. Psychiatry461140–1152. 10.1016/s0006-3223(99)00157-2
13
FreemanK.BrureauA.VadigepalliR.StaehleM. M.BrureauM. M.GonyeG.et al (2012a). Temporal changes in innate immune signals in a rat model of alcohol withdrawal in emotional and cardiorespiratory homeostatic nuclei.J. Neuroinflamm.9:97. 10.1186/1742-2094-9-97
14
FreemanK.StaehleM. M.GümüşZ. H.VadigepalliR.GonyeG. E.NicholsC. N.et al (2012b). Rapid temporal changes in the expression of a set of neuromodulatory genes during alcohol withdrawal in the dorsal vagal complex: molecular evidence of homeostatic disturbance.Alcohol. Clin. Exp. Res.361688–1700. 10.1111/j.1530-0277.2012.01791.x
15
FreemanK.StaehleM. M.VadigepalliR.GonyeG. E.OgunnaikeB. A.HoekJ. B.et al (2013). Coordinated dynamic gene expression changes in the central nucleus of the amygdala during alcohol withdrawal.Alcohol. Clin. Exp. Res.37E88–E100. 10.1111/j.1530-0277.2012.01910.x
16
GeislerR. F.HunterB. E.WalkerD. W. (1978). Ethanol dependence in the rat: temporal changes in neuroexcitability following withdrawal.Psychopharmacology56287–292. 10.1007/BF00432851
17
GrodinE.BujarskiB.BurnetteE.NietoS.LimA.LinJ.et al (2021). Ibudilast, a neuroimmune modulator, reduces heavy drinking and alcohol cue-elicited neural activation: a randomized trial.Transl. Psychiatry11:355. 10.1038/s41398-021-01478-5
18
GuG.RolandB.TomaselliK.DolmanC. S.LoweC.HeiligJ. S.et al (2013). Glucagon-like peptide-1 in the rat brain: distribution of expression and functional implication.J. Comp. Neurol.5212235–2261.
19
HanV. K.HynesM. A.JinC.TowleA. C.LauderJ. M.LundP. K.et al (1986). Cellular localization of proglucagon/glucagon-like peptide I messenger RNAs in rat brain.J. Neurosci. Res.1697–107. 10.1002/jnr.490160110
20
HermanJ. P. (2012). Neural pathways of stress integration relevance to alcohol abuse.Curr. Rev.34441–447.
21
HoltM. K.TrappS. (2016). The physiological role of the brain GLP-1 system in stress.Cogent Biol.2:1229086. 10.1080/23312025.2016.1229086
22
HunterB. E.RileyJ. N.WalkerD. W.FreundG. (1975). Ethanol dependence in the rat: a parametric analysis.Pharmacol. Biochem. Behav.3619–629. 10.1016/0091-3057(75)90183-5
23
JerlhagE. (2020). Alcohol-mediated behaviours and the gut-brain axis; with focus on glucagon-like peptide-1.Brain Res.1727:146562. 10.1016/j.brainres.2019.146562
24
KaliaM.FuxeK.GoldsteinM. (1985). Rat medulla oblongata. II. Dopaminergic, noradrenergic (A1 and A2) and adrenergic neurons, nerve fibers, and presumptive terminal processes.J. Comp. Neurol.233308–332. 10.1002/cne.902330303
25
KingA. C.ErricoA. L.ParsonsO. A.LovalloW. R. (1991). Blood Pressure Dysregulation Associated with Alcohol Withdrawal.Alcohol. Clin. Exp. Res.15478–482.
26
KoobG. F. (2009). Brain stress systems in the amygdala and addiction.Brain Res.129361–75.
27
KoobG. F. (2014). Neurocircuitry of alcohol addiction: synthesis from animal models.Handb. Clin. Neurol.12533–54. 10.1016/b978-0-444-62619-6.00003-3
28
KoobG. F.Le MoalM. (2001). Drug Addiction, Dysregulation of Reward, and Allostasis.Neuropsychopharmacology2497–129.
29
KoobG. F.Le MoalM. (2008a). Addiction and the Brain Antireward System.Annu. Rev. Psychol.5929–53. 10.1146/annurev.psych.59.103006.093548
30
KoobG. F.Le MoalM. (2008b). Review. Neurobiological mechanisms for opponent motivational processes in addiction.Philos. Trans. R. Soc. Lond. B. Biol. Sci.3633113–3123. 10.1098/rstb.2008.0094
31
KostenT. R.O’ConnorP. G. (2003). Management of Drug and Alcohol Withdrawal.N. Engl. J. Med.3481786–1795.
32
LarsenP. J.Tang-ChristensenM.HolstJ. J.OrskovC. (1997). Distribution of glucagon-like peptide-1 and other preproglucagon-derived peptides in the rat hypothalamus and brainstem.Neuroscience77257–270. 10.1016/s0306-4522(96)00434-4
33
LieberC. S.DecarliL. M. (1974). An experimental model of alcohol feeding and liver injury in the baboon.J. Med. Primatol.3153–163. 10.1159/000459999
34
LieberC. S.DecarliL. M. (1994). Animal models of chronic ethanol toxicity.Methods Enzymol.233585–594. 10.1016/S0076-6879(94)33061-1
35
MaceyD. J.SchulteisG.HeinrichsS. C.KoobG. F. (1996). Time-dependent quantifiable withdrawal from ethanol in the rat: Effect of method of dependence induction.Alcohol13163–170. 10.1016/0741-8329(95)02030-6
36
ManiscalcoJ. W.RinamanL. (2018). Vagal Interoceptive Modulation of Motivated Behavior.Physiology33151–167. 10.1152/physiol.00036.2017
37
MartinezF. O.GordonS. (2014). The M1 and M2 paradigm of macrophage activation: time for reassessment.F1000Prime Rep.6:13. 10.12703/P6-13
38
McBrideW. J.KimpelM. W.SchultzJ. A.McClintickJ. N.EdenbergH. J.BellR. L.et al (2010). Changes in gene expression in regions of the extended amygdala of alcohol-preferring rats after binge-like alcohol drinking.Alcohol44171–183. 10.1016/j.alcohol.2009.12.001
39
McDonaldM.HoekJ.OgunnaikeB.SchwaberJ. (2008). Behavioral and neurobiological changes within a period of heightened susceptibility to voluntary alcohol withdrawal.FASEB J. Federat. Am. Soc. Exp. Biol.22, 946.7–946.7. 10.1096/fasebj.22.1_supplement.946.7
40
MeckelK. R.KiralyD. D. (2019). A potential role for the gut microbiome in substance use disorders.Psychopharmacology vol.2361513–1530.
41
MurrayP. J.AllenJ. E.BiswasS. K.FisherE. A.GilroyD. W.GoerdtS.et al (2014). Macrophage Activation and Polarization: nomenclature and Experimental Guidelines.Immunity4114–20. 10.1016/j.immuni.2014.06.008
42
O’SullivanS. J.SchwaberJ. S. (2021). Similarities in alcohol and opioid withdrawal syndromes suggest common negative reinforcement mechanisms involving the interoceptive antireward pathway.Neurosci. Biobehav. Rev.125355–364. 10.1016/j.neubiorev.2021.02.033
43
O’SullivanS. J.ReyesB. A. S.VadigepalliR.Van BockstaeleE. J.SchwaberJ. S. (2020). Combining laser capture microdissection and microfluidic qpcr to analyze transcriptional profiles of single cells: a systems biology approach to opioid dependence.J. Vis. Exp.2020:e60612. 10.3791/60612
44
O’SullivanS. J.MalahiasE.ParkJ.SrivastavaA.ReyesB. A. S.GorkyJ.et al (2019). Single-Cell Glia and Neuron Gene Expression in the Central Amygdala in Opioid Withdrawal Suggests Inflammation With Correlated Gut Dysbiosis.Front. Neurosci.13:665. 10.3389/fnins.2019.00665
45
ParkJ.AnthonyB.KateK.AlexandriaS.SonaliG.BabatundeO.et al (2014). Inputs drive cell phenotype variability.Genome Res.24930–941. 10.1101/gr.161802.113
46
ParkJ.HaisunZ.SeanO. S.OgunnaikeB. A.WeaverD. R.SchwaberJ. S.et al (2016). Single-Cell Transcriptional Analysis Reveals Novel Neuronal Phenotypes and Interaction Networks Involved in the Central Circadian Clock.Front. Neurosci.10:481. 10.3389/fnins.2016.00481
47
PhelpsE. A.LeDouxJ. E. (2005). Contributions of the Amygdala to Emotion Processing: from Animal Models to Human Behavior.Neuron48175–187.
48
RayL. A.BujarskiS.ShoptawS.RocheD. J.HeinzerlingK.MiottoK.et al (2017). Development of the Neuroimmune Modulator Ibudilast for the Treatment of Alcoholism: a Randomized, Placebo-Controlled, Human Laboratory Trial.Neuropsychopharmacology421776–1788. 10.1038/npp.2017.10
49
ReinerP.VincentS. (1986). The distribution of tyrosine hydroxylase, dopamine-beta-hydroxylase, and phenylethanolamine-N-methyltransferase immunoreactive neurons in the feline medulla oblongata.J. Comp. Neurol.248518–531.
50
RetsonT. A.SterlingR. C.Van BockstaeleE. J. (2015). Alcohol-induced dysregulation of stress-related circuitry: The search for novel targets and implications for interventions across the sexes.Prog. Neuropsychopharmacol. Biol. Psychiatry10.1016/j.pnpbp.2015.05.009
51
RinamanL. (1999). Interoceptive stress activates glucagon-like peptide-1 neurons that project to the hypothalamus.Am. J. Physiol.277R582–R590.
52
RobertoM.KirsonD.KhomS. (2020). The Role of the Central Amygdala in Alcohol Dependence.Cold Spring Harb. Perspect. Med.11:a039339. 10.1101/cshperspect.a039339
53
SchäfersM.SorkinL. (2008). Effect of cytokines on neuronal excitability.Neurosci. Lett.437188–193. 10.1016/j.neulet.2008.03.052
54
SEQC/MAQC-III Consortium. (2014). A comprehensive assessment of RNA-seq accuracy, reproducibility and information content by the Sequencing Quality Control Consortium.Nat. Biotechnol.32903–914. 10.1038/nbt.2957
55
SkosnikP. D.Cortes-BrionesJ. A. (2016). Targeting the ecology within: the role of the gut–brain axis and human microbiota in drug addiction.Med. Hypotheses9377–80. 10.1016/j.mehy.2016.05.021
56
VallöfD.VestlundJ.JerlhagE. (2019). Glucagon-like peptide-1 receptors within the nucleus of the solitary tract regulate alcohol-mediated behaviors in rodents.Neuropharmacology149124–132. 10.1016/j.neuropharm.2019.02.020
57
WakiH.GouraudS. S.MaedaM.PatonJ. F. R. (2010). Evidence of specific inflammatory condition in nucleus tractus solitarii of spontaneously hypertensive rats.Exp. Physiol.95595–600.
58
WalkerD. W.HunterB. E.RileyJ. (1975). A behavioral and electrophysiological analysis of ethanol dependence in the rat.Adv. Exp. Med. Biol.59353–372. 10.1007/978-1-4757-0632-1_25
59
WangY.ChiuW.ChengC.LeeC.ChihW. H. (2021). A. Examination of neuroinflammatory cytokine interleukin-1 beta expression in the medial prefrontal cortex, amygdala, and hippocampus for the paradoxical effects of reward and aversion induced by morphine.Neurosci. Lett.760:136076. 10.1016/j.neulet.2021.136076
60
WhitmanB. A.KnappD. J.WernerD. F.CrewsF. T.BreeseG. R. (2013). The Cytokine mRNA Increase Induced by Withdrawal from Chronic Ethanol in the Sterile Environment of Brain is Mediated by CRF and HMGB1 Release.Alcohol. Clin. Exp. Res.372086–2097. 10.1111/acer.12189
61
WilliamsC. L.MenD.ClaytonE. C.GoldP. E. (1998). Norepinephrine release in the amygdala after systemic injection of epinephrine or escapable footshock: contribution of the nucleus of the solitary tract.Behav. Neurosci.1121414–1422. 10.1037/0735-7044.112.6.1414
62
WilsonJ. S.KorstenM. A.LieberC. S. (1986). The combined effects of protein deficiency and chronic ethanol administration on rat ethanol metabolism.Hepatology6823–829.
63
WolfY.YonaS.KimK. W.JungS. (2013). Microglia, seen from the CX3CR1 angle.Front. Cell. Neurosci.7:26. 10.3389/fncel.2013.00026
64
YangL.WangM.GuoY. Y.SunT.LiY. J.YangQ.et al (2016). Systemic inflammation induces anxiety disorder through CXCL12/CXCR4 pathway.Brain. Behav. Immun.56352–362. 10.1016/j.bbi.2016.03.001
65
ZhengH.CaiL.RinamanL. (2015). Distribution of glucagon-like peptide 1-immunopositive neurons in human caudal medulla.Brain Struct. Funct.2201213–1219. 10.1007/s00429-014-0714-z
Summary
Keywords
alcohol withdrawal, neuroinflammation, RT-PCR, subphenotypes, single-cell heterogeneity, microglia
Citation
O’Sullivan SJ, McIntosh-Clarke D, Park J, Vadigepalli R and Schwaber JS (2021) Single Cell Scale Neuronal and Glial Gene Expression and Putative Cell Phenotypes and Networks in the Nucleus Tractus Solitarius in an Alcohol Withdrawal Time Series. Front. Syst. Neurosci. 15:739790. doi: 10.3389/fnsys.2021.739790
Received
12 July 2021
Accepted
22 October 2021
Published
19 November 2021
Volume
15 - 2021
Edited by
James W. Grau, Texas A&M University, United States
Reviewed by
Thomas Alexander Lutz, University of Zurich, Switzerland; Morgane Thomsen, Region Hovedstad Psychiatry, Denmark
Updates

Check for updates
Copyright
© 2021 O’Sullivan, McIntosh-Clarke, Park, Vadigepalli and Schwaber.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: James S. Schwaber, James.Schwaber@jefferson.edu
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.