A structural connectivity atlas of limbic brainstem nuclei

Levinson, Simon; Miller, Michelle; Iftekhar, Ahmed; Justo, Monica; Arriola, Daniel; Wei, Wenxin; Hazany, Saman; Avecillas-Chasin, Josue M.; Kuhn, Taylor P.; Horn, Andreas; Bari, Ausaf A.

doi:10.3389/fnimg.2022.1009399

ORIGINAL RESEARCH article

Front. Neuroimaging, 12 January 2023
Sec. Brain Imaging Methods
Volume 1 - 2022 | https://doi.org/10.3389/fnimg.2022.1009399

A structural connectivity atlas of limbic brainstem nuclei

Simon Levinson^1,2

Michelle Miller¹

Ahmed Iftekhar¹

Monica Justo¹

Daniel Arriola¹

Wenxin Wei¹

Saman Hazany³

Josue M. Avecillas-Chasin⁴

Taylor P. Kuhn⁵

Andreas Horn^6,7,8

Ausaf A. Bari¹^*

¹Department of Neurosurgery, David Geffen School of Medicine at the University of California, Los Angeles, Los Angeles, CA, United States
²Stanford Department of Neurosurgery, Stanford University, Palo Alto CA, United States
³Department of Radiology, VA Greater Los Angeles Healthcare System, David Geffen School of Medicine at UCLA, Los Angeles, CA, United States
⁴Department of Neurosurgery, University of Nebraska Medical Center, Omaha, NE, United States
⁵Department of Psychiatry and Biobehavioral Sciences, University of California, Los Angeles, Los Angeles, CA, United States
⁶Movement Disorder and Neuromodulation Unit, Department of Neurology, Charité–Universitätsmedizin Berlin, Corporate Member of Freie Universität Berlin and Humboldt–Universität zu Berlin, Berlin, Germany
⁷Department of Neurology, Center for Brain Circuit Therapeutics, Harvard Medical School, Brigham and Women's Hospital, Boston, MA, United States
⁸Massachusetts General Hospital Neurosurgery and Center for Neurotechnology and Neurorecovery (CNTR) at MGH Neurology Massachusetts General Hospital, Harvard Medical School, Boston, MA, United States

Background: Understanding the structural connectivity of key brainstem nuclei with limbic cortical regions is essential to the development of therapeutic neuromodulation for depression, chronic pain, addiction, anxiety and movement disorders. Several brainstem nuclei have been identified as the primary central nervous system (CNS) source of important monoaminergic ascending fibers including the noradrenergic locus coeruleus, serotonergic dorsal raphe nucleus, and dopaminergic ventral tegmental area. However, due to practical challenges to their study, there is limited data regarding their in vivo anatomic connectivity in humans.

Objective: To evaluate the structural connectivity of the following brainstem nuclei with limbic cortical areas: locus coeruleus, ventral tegmental area, periaqueductal grey, dorsal raphe nucleus, and nucleus tractus solitarius. Additionally, to develop a group average atlas of these limbic brainstem structures to facilitate future analyses.

Methods: Each nucleus was manually masked from 197 Human Connectome Project (HCP) structural MRI images using FSL software. Probabilistic tractography was performed using FSL's FMRIB Diffusion Toolbox. Connectivity with limbic cortical regions was calculated and compared between brainstem nuclei. Results were aggregated to produce a freely available MNI structural atlas of limbic brainstem structures.

Results: A general trend was observed for a high probability of connectivity to the amygdala, hippocampus and DLPFC with relatively lower connectivity to the orbitofrontal cortex, NAc, hippocampus and insula. The locus coeruleus and nucleus tractus solitarius demonstrated significantly greater connectivity to the DLPFC than amygdala while the periaqueductal grey, dorsal raphe nucleus, and ventral tegmental area did not demonstrate a significant difference between these two structures.

Conclusion: Monoaminergic and other modulatory nuclei in the brainstem project widely to cortical limbic regions. We describe the structural connectivity across the several key brainstem nuclei theorized to influence emotion, reward, and cognitive functions. An increased understanding of the anatomic basis of the brainstem's role in emotion and other reward-related processing will support targeted neuromodulatary therapies aimed at alleviating the symptoms of neuropsychiatric disorders.

Highlights

– The brainstem plays a key role in the processing of emotional stimuli and is intricately linked with the limbic system.

– Anatomic data for these connections is limited in humans.

– We describe the structural connectivity of five brainstem nuclei (locus coeruleus, ventral tegmental area, periaqueductal grey, dorsal raphe nucleus, and nucleus tractus solitarius) in relation to limbic circuits.

– Our results present a comprehensive delineation of the brainstem-limbic structural connectivity of these nuclei and are compiled into a freely available tractographic atlas.

– Applications include future targeting of these structures for neuropsychiatric conditions.

Introduction

The brainstem, comprised of the medulla oblongata, pons, and midbrain, comprises ~3% of the mass of the brain and contains about 2% of the neurons in the central nervous system. Yet, what it lacks in size it makes up for in complexity and a disproportional influence on processes ranging from autonomic functions to arousal and consciousness (Azevedo et al., 2009; Herculano-Houzel, 2009). The brainstem provides vital autonomic regulation and homeostatic maintenance. It also serves as a conduit for all fibers linking the cerebral cortex and cerebellum to the spinal cord. Furthermore, it functions as a major afferent sensory system, receiving input from visceral fibers and cranial nerve nuclei which it then filters and transmits (often through several intermediary nuclei), to higher cortical centers.

The brainstem also plays a key role in emotional processing. Three brainstem networks have been identified that are thought to contribute to limbic processing: (1) the ascending sensory network consisting of the spinothalamic tracts, medial forebrain bundle, nucleus of the tractus solitarius (NTS), parabrachial nuclear complex and thalamic nuclei, (2) the descending motor network consisting of the periaqueductal grey (PAG), caudal raphe nucleus, and locus coeruleus (LC), and (3) the modulatory network with the serotonergic dorsal raphe (DRN), noradrenergic LC, and dopaminergic ventral tegmental area (VTA) (Angeles Fernández-Gil et al., 2010; Venkatraman et al., 2017). The ascending, descending and modulatory brainstem networks allow for the progressive integration and processing of information as signals move rostrally through the brainstem, thalamus and then to the cortex, but also carry information in reverse, with cortical regions regulating the action-response relationships of phylogenetically older structures (Tucker et al., 2002).

The anatomic and structural basis for the brainstem's role in limbic processing is theorized to involve several key nuclei which are the sole or major source of potent monoamine neurotransmitters for the higher cerebral cortex: the LC (norepinephrine), DRN (serotonin), and VTA (dopamine). Additionally, the PAG and NTS serve as inputs or centers of modulation to these monoamine neural networks (Figure 1). Several monoamine neurotransmitters have, individually or in combination, been implicated in disease states including Parkinson's disease, major depressive disorder, or addiction, and are essential for physiological activities including arousal, sleep/wake cycles, perception of pain, affect, and goal directed behavior. Existing evidence obtained largely from animal studies (Table 1) describes how these brainstem nodes interact with cortical limbic structures to convey body state and homeostatic information to result in behaviors such as heightened alertness, arousal from sleep, fear, and defense measures (Omar et al., 2008; Ulrich-Lai and Herman, 2009). Lesion studies have also demonstrated the emotional role of the brainstem utilizing models or subjects with disease states such as Alzheimer's Disease, brainstem-cerebellar pathology, and brainstem infarcts (van Zandvoort et al., 2003; Omar et al., 2008; Mariën and D'aes, 2015; Dutt et al., 2021). However, despite their importance, it has been challenging to study brainstem nuclei in humans in vivo because of difficulties in accurately defining the nuclei on imaging (Edlow et al., 2012; Song et al., 2014; Tang et al., 2018). Therefore, limited anatomic data exists for these brainstem-limbic relationship in humans.

FIGURE 1

Figure 1. Schematic representation of key brainstem-cortex circuits. DLPFC, dorsolateral-prefrontal cortex; rACC, rostral anterior cingulate cortex; INSULA, insular cortex; AMY, amygdala; HIPPO, hippocampus; VTA, ventral tegmental area; PAG, periaqueductal grey; DRN, dorsal raphe nucleus; LC, locus coeruleus; NTS, nucleus tractus solitarius; Pb, parabrachial nucleus; NE, norepinephrine; DA, dopamine; 5-HT, serotonin.

TABLE 1

Table 1. Limbic brainstem nuclei key characteristics.

There are two general approaches to defining regions of interest (ROI) on MRI, manual and automated segmentation. Automated segmentation has been demonstrated to reliably delineate cortical structures (Fischl et al., 2002; Woolrich et al., 2009). It has also been used successfully to define certain regions of the brainstem (particularly for studies interested in masking the medulla, pons and midbrain separately) (Patenaude et al., 2011; Iglesias et al., 2015; Wang et al., 2016; Sander et al., 2019), yet it remains technically difficult to define most individual brainstem nuclei via this method. While there are several automated techniques, most recent attempts have been to mask brainstem nuclei using convolutional neural network-based segmentation, a deep learning image recognition technique, to delineate the substantia nigra (Berre et al., 2019) and LC (Dünnwald et al., 2020). However, this method relies on neuromelanin-MRI scans, limiting its application to non-pigmented regions. Furthermore, other techniques such as voxel intensity-based algorithms still require some manual delineation and thresholding and can be complicated by homogenously intense regions (Ogisu et al., 2013; Zhang et al., 2016; Berre et al., 2019; Dünnwald et al., 2020). Therefore, although automated processing technologies are being developed, manual segmentation still has significant advantages and is considered the gold standard.

Manual segmentation is labor intensive, often resulting in studies with small sample sizes and high inter-rater variability. In addition, utilizing a manually-defined mask as an atlas for new subjects can be more computationally intensive than automated segmentation techniques (Aljabar et al., 2009; Berre et al., 2019). Yet, a well-trained individual, given proper anatomic knowledge and tools, can produce reliable results and manually delineated atlases are still considered by many to be the gold standard (Aljabar et al., 2009; Morey et al., 2009; Iglesias et al., 2015).

In this study, we manually segmented five brainstem nuclei to perform probabilistic tractography to selected limbic targets in 197 human subjects. We hypothesized that autonomic and monoamine brainstem nuclei would demonstrate structural connectivity to cortical limbic regions, namely the amygdala, insula and hippocampus. Additionally, since there were no existing atlases of the structural connectivity of these nuclei, we aimed to produce an accurate anatomical atlas for use in future research.

We overcame some of the challenges of performing tractography on the brainstem by using a rigorous anatomic definition scheme with a combination of voxel measurements and anatomic landmarks to reduce variability while defining each nucleus. Additionally, we were able to achieve a scale of nearly 200 subjects, similar to many automated based approaches, thereby increasing the power of the study and reducing the effects of outliers or isolated errors. Lastly, we obtained diffusion MRI (dMRI) scans from the Human Connectome Project which provided a large dataset with scans acquired in a high resolution allowing for increased accuracy (Glasser et al., 2013; Van Essen et al., 2013). Overall, our results present a comprehensive delineation of the brainstem-limbic structural connectivity of these nuclei and are compiled into a freely available tractographic atlas.

Methods

Subjects

Data were obtained from the publicly available WU-Minn HCP 1,200 Subjects data release repository (Glasser et al., 2013; Van Essen et al., 2013). The scanning protocol was approved by Human Research Protection Office (HRPO), Washington University (IRB# 201 204 036). No human subject experimental procedures were undertaken at the authors' home institutions. The participants included in the HCP 1,200 Subjects data release provided written informed consent as approved by the Washington University IRB. From this repository, 200 total non-twin subjects were randomly selected. The analysis was limited to these subjects based on available computational resources. Three subjects were excluded due to lack of required imaging data files. The remaining 197 subjects were included in our analyses and a description of their demographic characteristics is provided in Table 2.

TABLE 2

Table 2. Cohort.

MRI acquisition

The data were acquired on a modified Siemens 3T Skyra scanner with a customized protocol²¹. The T1-weighted MRI has an isotropic spatial resolution of 0.7 mm, and the dMRI data have an isotropic spatial resolution of 1.25 mm. The multi-shell dMRI data were collected over 270 gradient directions distributed over three b-values (1,000, 2,000, 3,000 s/mm²). For each subject, the multi-shell dMRI data were collected with both L/R and R/L phase encodings using the same gradient table, which were then merged into a single copy of multi-shell dMRI data after the correction of distortions with the HCP Preprocessing Pipeline. Average T1w and T2w images were then aligned to MNI space (with 0.7 mm resolution), with 6 degrees of freedom which also aligns with AC-PC line and the inter-hemispheric plane but which does not alter the original size and shape of the brain. Acquisition time was 32 min for T1w scans which the majority of the atlas used here was based on [FOV = 224 mm, matrix = 320, 256 sagittal slices in a single slab, TR=24,00 ms, TE = 2.14 ms, TI = 1,000 ms, FA=8°, Bandwidth (BW) = 210 Hz per pixel] (Glasser et al., 2013; Sotiropoulos et al., 2013).

Masking of seed structures and anatomic boundaries

The areas of interest in this study consisted of brainstem structures with widespread projections. Due to the wide inter-subject anatomic variation and the small nature of many of these structures, all seed masks were created manually for each subject on the original Human Connectome Project (HCP) structural MRI scans using FSLeyes software. The seed masks generated in this fashion included the Locus Coeruleus (LC), Nucleus Tractus Solitarius (NTS), Periaqueductal Grey (PAG), Dorsal Raphe Nucleus (DRN) and Ventral Tegmental Area (VTA).

The anatomic boundaries used for each set of masks are described below and shown in Figure 2 and Anatomic Supplement 1. Multiple histologic and radiographic resources were used to carefully define each region of interest (Desikan et al., 2006; Naidich et al., 2009; Edlow et al., 2012; Mai et al., 2016; Vanderah, 2019). Individual raters were first trained on a sample data set and their accuracy was assessed relative to a template mask. Raters had to achieve an error rate of <5% to begin masking subjects. To ensure a high degree of internal consistency, no more than two individuals were responsible for creating the masks of each nucleus. We utilized a two expert review method where an expert neuroradiologist and neurosurgeon evaluated each mask based on the below standardized anatomic boundaries. This provided internal validity. After review between the two expert raters, there was agreement that all masks were within their tolerance as described by the anatomic protocols in Figure 2. Between 30 and 60 min were spent per subjects creating, editing, and reviewing the 5 seed masks.

FIGURE 2

Figure 2. Anatomic boundaries. Diagram of how anatomic boundaries of nuclei were masked. First radiographic²³ (left pane) and histologic atlas²⁴ (middle pane) were carefully examined for each nucleus. Using a combination of anatomic landmarks and voxel measurements (not shown), the resulting mask (right pane) was created in FSL. Locus Coreleus (A), Dorsal Raphe Nucleus (B), Solitary nucleus (C), Ventral tenemental area (D), Periaquaductal Grey (E).

After initial masking was completed, the fslmaths boxv command was used to dilate each mask by a factor of 5 followed by an erosion of 5 to ensure edges were smoothed and gaps were filled.

Locus coeruleus

The LC was defined in both the caudal-rostral and medial-lateral planes. The most caudal point was the midpoint of the 4^th ventricle and the most rostral was the rostral pons. Laterally, the LC was defined to be 3 mm lateral to the midline at the anterior-lateral angle of the 4^th ventricle (Figure 2A). Average volume of the seed mask was 94.5 mm³.

Dorsal raphe nucleus

Beginning at the caudal midbrain, the DRN was masked caudally until its termination in the mid pons. The cerebral aqueduct and fourth ventricle were used as guides with the rostral portions measuring 12 voxels wide tapering to 6 voxels at the most caudal aspect (Figure 2B). Average volume of the seed mask was 153.8 mm³.

Nucleus tractus solitarius

The caudal aspect of the NTS was defined to begin 3 mm rostral to the obex. It then proceeded rostrally in a “V” shape. Standardized voxel measurements were used to ensure a consistent shape across individual. The NTS was masked caudally until the middle cerebellar peduncles were clearly visible on a horizontal section. Since the NTS also tracks slightly more anteriorly as it progresses rostrally, the anterior-caudal boundary was defined in the horizontal plane as the point where a line would transect from the root of cranial nerve VIII, with the point halfway between the midpoint and lateral aspect of the 4^th ventricle on either side of the pons (Figure 2C). Average volume of the seed mask was 207.3 mm³.

Ventral tegmental area

The most inferior transverse section of the VTA was defined at the section of the midbrain where the interpeduncular fossa opens to the interpeduncular cistern. The VTA boundary on inferior sections was the medial border of the Substantia Nigra. The lateral and medial borders were directly adjacent to the interpeduncular fossa. As the VTA progresses rostrally it becomes a contiguous structure with its medial boundaries joining in the midline bordering the medial aspect of the interpeduncular fossa and extending in a posterior direction to the midpoint of the medial edge of the red nucleus. The most superior section of the VTA was defined at the level of the mammillary bodies. Standardized voxel measurements were used to ensure a consistent shape across individual subjects (Figure 2D). Average volume of the seed mask was 214.1 mm³.

Periaqueductal grey

From a mid-coronal slice, the aqueduct was identified on the sagittal perspective. On the axial perspective, the lateral ventricles were traced inferiorly until it became the cerebral aqueduct which delineated the superior margin of the mask. A diamond-shaped border surrounding the aqueduct was demarcated as the PAG on the axial plane. The inferior border of the mask was determined to be the location where the mammillary bodies became fully defined, which correlated to the opening of the fourth ventricle. Standardized voxel measurements were used to ensure a consistent shape across individual subjects (Figure 2E). Average volume of the seed mask was 172.2 mm³.

Probabilistic tractography

Probabilistic tractography was carried out using FSL's FMRIB Diffusion Toolbox (probtrackx) with modified Euler streaming (Woolrich et al., 2009; Jenkinson et al., 2012). We used the diffusion tensor model that was fitted on processed diffusion data. ROI were delineated using standard FSL parameters (specified at http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/) to derive fractional anisotropy maps for each ROI. The bedpost command was used to generate an estimate of crossing fiber orientation at the level of individual voxels. Then, utilizing a so called “ball and stick” model, we selected 3 sticks and left all other options as default. Importantly, at each voxel, all possible fiber orientations were considered taking into account standard assumptions about uncertainty when computing whether a given set of voxels will be considered to be “in-line” so as to count as a streamline. For more detail on the underlying algorithms please see these references (Behrens et al., 2003, 2007; Woolrich et al., 2009).

Target masks were generated using the Harvard-Oxford subcortical atlas (Desikan et al., 2006). Target masks included the amygdala (AMY), dorsolateral prefrontal cortex (DLPFC), hippocampus, insula, Nucleus accumbens (NAc), orbitofrontal cortex (OFC), and rostral anterior cingulate cortex (rACC). Seed masks also served as target masks once produced such that the number of targets increased overtime as new seeds were created. All tractography was performed between each (right and left) seed and the ipsilateral target. Each target mask was also a termination mask such that tractography was terminated once a streamline entered the target. Additionally, Freesurfer standard lookup tables were used to generate Ipsilateral white matter masks which were used as waypoints. The ventricles and cerebellum masks were similarly generated with Freesufer and used as exclusion masks. We used the “onewaycondition,” curvature 0.2, 2,000 samples, steplength = 0.5, fibthresh = 0.01, distthresh = 1 and sampvox = 0.0. This resulted in 14 or more seed_to_target output files representing a voxelwise map of the number of seed samples from each seed to target.

To calculate the probability of connectivity (POC) between each brainstem seed voxel to each of the 7 cortical and to the other 4 brainstem nuclei targets, we ran the FSL proj_thresh subroutine with a threshold of 1,250 on each probtrackx output. For each voxel in the seed mask with a value above the threshold, proj_thresh calculates the number of samples reaching each of the target masks as a proportion of the total number of samples reaching any of the target masks. The fsl waytotal and waynorm commands were used to normalize the tractography results by seed and target size. This yielded a separate map of each seed mask for each target with each voxel having a value between 0 and 1 representing the POC of that voxel to the given target. Thus, there were 2 maps (one for each hemisphere) per seed and per target for each subject. To produce an overall POC from each seed mask to target, probabilities were averaged across all voxels in each map. We next created a population connectivity map across all 197 subjects. Each of the previously created proj_thresh maps was registered to MNI 1 mm standard space, thresholded at a level of 0.1 and binarized. These maps were then added across all 197 subjects such that each voxel value now represented the number of subjects with connectivity to the target.

Parallel data processing

FSL software was implemented in a distributed fashion using Amazon Web Services (AWS, http://aws.amazon.com) EC2 instances running in parallel. Each AWS EC2 instance was an r4.large clone of an Amazon Machine Image (AMI) running Ubuntu 14.04 with FSL software version 5.0.10. FSL bedpostx directories for each subject and the probtrackx output files were stored on an Amazon S3 bucket.

Statistical analysis

Statistical analysis was carried out using the R software package (http://www.r-project.org/) and Prism 8 software from Graphpad (https://www.graphpad.com/).

To analyze the structural data, the subject specific output from the proj_thresh files stored on Amazon S3 were downloaded to a new EC2 instance running Ubuntu 14.04.1. FSLmaths was then used to compute the mean connectivity to each target using fslstats. The means were then imported to RStudio (version 1.3.959) running R (version 3.6.3) and the means to each target region were compared with a single factor analysis of variance (ANOVA). A Tukey HSD test was then run to determine statistical significance of the variance in the means.

Overall connectivity measurements were obtained by first taking the average of right and left connectivity for each subject specific seed-target relationship and then computing the mean across all subjects.

MNI atlas

Separately, each seed mask was warped to MNI space using the FSL applywarp command and the HCP subject specific nonlinear acpc_dc2standard file. The MNI152 NLIN 2009b T1 0.5 mm brain was used as reference (ICBM 152 Nonlinear atlases, 2009). Once warped, each mask was averaged across all 197 subjects and set the threshold to two standard deviations from the mean to exclude extreme values.

Results

Summative images in MNI space are publicly available at: https://www.uclahealth.org/neurosurgery/research-areas. Additionally, the atlas comes preinstalled on the widely used neurostimulation software, Lead-DBS (https://www.lead-dbs.org/helpsupport/knowledge-base/atlasesresources/atlases). All other data is available upon request of the corresponding author.

Complete results for the average POC for each seed to target relationship is reported in Table 3. Figures 3–7 details the average streamline paths, provides a graphic representation of mean connectivity, and portrays the anatomic boundaries of each seed mask averaged over all subjects in MNI space. The three highest POCs are reported here for each seed. A summative three-dimensional atlas generated with Lead-DBS utilizing an MNI-152 0.5 mm template brain is also showed (Figure 8).

TABLE 3

Table 3. Average probabilities of connectivity; format: Mean [95% CI interval].

FIGURE 3

Figure 3. Dorsal raphe nucleus structural results. (A) MNI space structural connectivity results visual representation averaged over all 197 subjects. Brighter yellow on heat map indicates a high number of samples passing through a given point that will eventually reach a target map (brighter yellow = more samples). Dark green: DLPFC; pink: OFC; brown: AMY; blue: HIPPO; purple: insula; orange: NAc; light green: rACC. (B) Mean connectivity results with dashed line showing mean and 95% CI, each point on graph shows result from individual subject. (C) Anatomic MNI mask of seed region.

Within each seed region, ANOVA with multiple comparisons (ANOVA-MC) and post-hoc t-test decompositions demonstrated significant (adjusted p<0.05) difference in the POC between each seed and target structure except as noted below. The only other exceptions to this were if the total POC was <0.03 for both structures.

The dorsal raphe nucleus demonstrated average POC to each target as follows: AMY 0.204 [95% CI 0.195, 0.214], VTA 0.193 [95% CI 0.177, 0.208], DLPFC 0.189 [95% CI 0.179, 0.200]. However, ANOVA-MC did not demonstrate a significant difference between these structures (Figure 3).

The locus coeruleus demonstrated the greatest POC to the DLPFC 0.349 [95% CI 0.333, 0.364], followed by the AMY 0.273 [95% CI 0.262, 0.284] and Hippocampus 0.150 [95% CI 0.143, 0.157] (Figure 4).

FIGURE 4

Figure 4. Locus coeruleus structural results. (A) MNI space structural connectivity results visual representation averaged over all 197 subjects. Brighter yellow on heat map indicates a high number of samples passing through a given point that will eventually reach a target map (brighter yellow = more samples). Dark green: DLPFC; pink: OFC; brown: AMY; blue: HIPPO; purple: insula; orange: NAc; light green: rACC. (B) Mean connectivity results with dashed line showing mean and 95% CI, each point on graph shows result from individual subject. (C) Anatomic MNI mask of seed region.

The nucleus tractus solitarius showed the greatest POC to the LC 0.260 [95% CI 0.244, 0.275], followed by the DLPFC 0.257 [95% CI 0.244, 0.270] and the AMY 0.184 [95% CI 0.176, 0.192]. There was no significant difference between NTS-LC and NTS-DLPFC POC (Figure 5).

FIGURE 5

Figure 5. Nucleus tractus solitarius structural results. (A) MNI space structural connectivity results visual representation averaged over all 197 subjects. Brighter yellow on heat map indicates a high number of samples passing through a given point that will eventually reach a target map (brighter yellow = more samples). Dark green: DLPFC; pink: OFC; brown: AMY; blue: HIPPO; purple: insula; orange: NAc; light green: rACC. (B) Mean connectivity results with dashed line showing mean and 95% CI, each point on graph shows result from individual subject. (C) Anatomic MNI mask of seed region.

For the periaqueductal grey the highest POC was to the AMY 0.273 [95% CI 0.262, 0.283], followed by DLPFC 0.271 [95% CI 0.257, 0.285], however this difference was not significant (P > 0.99). The third greatest POC was to the hippocampus at 0.157, significantly less than either the AMY (p < 0.001) or the DLPFC (p < 0.001) [95% CI 0.151, 0.164] (Figure 6).

FIGURE 6

Figure 6. Periaqueductal grey structural connectivity. (A) MNI space structural connectivity results visual representation averaged over all 197 subjects. Brighter yellow on heat map indicates a high number of samples passing through a given point that will eventually reach a target map (brighter yellow = more samples). Dark green: DLPFC; pink: OFC; brown: AMY; blue: HIPPO; purple: insula; orange: NAc; light green: rACC. (B) Mean connectivity results with dashed line showing mean and 95% CI, each point on graph shows result from individual subject. (C) Anatomic MNI mask of seed region.

The ventral tegmental area showed the greatest POC to the DLPFC 0.252 [95% CI 0.236, 0.268], AMY 0.221 [95% CI 0.210, 0.232], and the hippocampus 0.174 [95% CI 0.166, 0.182]. Of note, there was no significant difference between VTA-DLPFC and VTA-AMY POC (p = 0.26) (Figure 7).

FIGURE 7

Figure 7. Ventral tegmental area structural connectivity. (A) MNI space structural connectivity results visual representation averaged over all 197 subjects. Brighter yellow on heat map indicates a high number of samples passing through a given point that will eventually reach a target map (brighter yellow = more samples). Dark green: DLPFC; pink: OFC; brown: AMY; blue: HIPPO; purple: insula; orange: NAc; light green: rACC; (B) Mean connectivity results with dashed line showing mean and 95% CI, each point on graph shows result from individual subject. (C) Anatomic MNI mask of seed region.

FIGURE 8

Figure 8. MNI 3D Atlas. Green: nucleus tractus solitaris; yellow: locus coeruleus; pink: periaqueductal grey; blue: ventral tegmental area; orange: dorsal raphe nucleus. For this figured, Lead-DBS was used to visualize the same volumetric MNI seed masks as shown in Figures 3–7 which are here overlaid on the MNI template brain to provide a 3-D representation of the relationship between these brainstem nuclei.

The extent to which each nucleus demonstrated laterality of connectivity was also calculated. The mean connectivity for each seed to target relationship was computed for the left and right and the difference between left and right was divided by the total connectivity to give a relative value for laterality. Left was given a negative value and right a positive value. Interestingly, there was no significant difference (p > 0.05, two way ANOVA with multiple comparison) between the mean laterality for each seed nuclei regardless of which target was selected (Figure 9).

FIGURE 9

Figure 9. Laterality. Arbitrary measure of laterality by seed nuclei. (Left) is set as negative and right as positive. More extreme values indicate greater difference in mean connectivity between (left) and (right).

The demographic characteristics of the data set was homogenous by design (Table 2). The t-test was performed between each seed-target POC and compared between males (n = 96) and females (n = 101). P > 0.05 for all relations except for the PAG-insula connection with male POC = 0.135 and female POC = 0.154 (P = 0.025).

Discussion

In this study, we carried out manual in vivo segmentation of defined brainstem regions of interest in a large cohort of subjects and performed probabilistic tractography between these regions and several limbic target regions of interest in the study of emotion and reward. The seeds selected were the NTS, LC, DRN, VTA, and PAG. These were chosen because of their clinical significance, the importance and wide-ranging effects of their respective neurotransmitters, and their established role in emotion and reward processing. We evaluated the connectivity between each of these seeds individually to a set of limbic processing centers to test the hypothesis of whether there would be a high POC between brainstem nuclei and these regions. Lastly, we aggregated these results into an MNI atlas and have made this publicly available to all interested researchers. This data represents a unique comprehensive in vivo analysis of key brainstem nuclei in a large human population and is also one of the first such studies to examine the relationship of these structures utilizing the Human Connectome Project.

Overall, our connectivity analysis of these five brainstem nuclei supports the role of the brainstem in emotional processing and confirmed our hypothesis that these brainstem nuclei demonstrate structural connectivity with known limbic regions. There is strong evidence that most mammals have analogous experiences to human emotion, even with far less evolved higher cortical centers, indicating that important processing occurs at levels below the cortex (Craig, 2003; Omar et al., 2008; Angeles Fernández-Gil et al., 2010; Holstege, 2016; Venkatraman et al., 2017). In order to further understand the basis of emotion and reward, and intervene in cases of disorders of emotion and reward processing, it is critical to understand phylogenetically older structures such as the brainstem and recognize them as critical factors in emotional processing.

Emotions can be thought of as mental representations associated with distinct sensory states with the evolutionary goal of producing relevant behavioral responses in an organism (Venkatraman et al., 2017). It is therefore logical that the brainstem is a key modulator of this system as it is able to receive afferent visceral and somatic sensory information, begin filtering and processing these signals, and transmit them to higher cortical centers. Existing literature suggests that many of these processes are related to the monoamine neuropeptides, norepinephrine, serotonin, and dopamine. These neurotransmitters are produced by these brainstem nuclei which then project to higher cortical centers (Holstege, 1992; Plutchik, 2001; Edlow et al., 2012; Venkatraman et al., 2017). Networks and nodes in the brainstem are increasingly considered key aspects of the conscious experience of emotion and therefore are of great clinical significance to neuropsychiatric disorders (Edlow et al., 2012, 2015; Venkatraman et al., 2017).

In our analysis, each seed region demonstrated its own pattern of connectivity largely consistent with existing anatomical information from other imaging, postmortem and animal studies (Naidich et al., 2009; Edlow et al., 2015, 2019; Bianciardi et al., 2016; Mai et al., 2016; Venkatraman et al., 2017; Tang et al., 2018). A general trend was observed for high probabilities of connectivity with the amygdala, DLPFC, and hippocampus, with relatively lower connectivity with the OFC, NAc, insula, and rACC, suggesting that while there is heterogeneity among these nuclei, there is also a distinct pattern of brainstem connectivity with both frontal and temporal lobe structures.

Dorsal raphe nucleus

The neurons which transmit serotonin and project to the cerebral cortex are mainly clustered around the dorsal and rostral aspects of the raphe nucleus (Hornung and De Tribolet, 1995). The results here pertain only to the rostral portion as the largely inferiorly projecting caudal half was selectively excluded in the anatomic boundaries for this analysis (Figure 2B). The serotonin system is widely studied and is known to modulate fear and anxiety and other social behaviors (Zangrossi et al., 2001; Moskowitz et al., 2003; Arbib and Fellous, 2004; Gobrogge et al., 2017). Treatments for depression and anxiety uses medications selective serotonin reuptake inhibitors (SSRIs) targeted at this system (Adell, 2015; Pillai et al., 2019). In animals, serotonergic projection to the amygdala and hippocampus has been associated with anxiety and retrieval of fear memories, and this relationship is under study in humans as a possible mechanistic explanation for major depressive disorder (Lowry et al., 2008; Dayan and Huys, 2009; Ohmura et al., 2010; Weinstein et al., 2015; Brakowski et al., 2017; Anand et al., 2018).

Based on histological data, the DRN's projection patterns are complex and include temporal lobe structures such as the hippocampus and amygdala, the PAG, LC, and frontal and insular cortices as well (Hornung, 2003). Functional connectivity studies have shown alterations in DRN connectivity to the frontal and cingulate cortices in subjects with depression (Weinstein et al., 2015; Anand et al., 2018). Given the DRN's significance in depression, there have been studies analyzing DRN-amygdala structural connectivity (Schmaal et al., 2016; Pillai et al., 2019). While, the DRN has been included in whole-brain brainstem connectivity studies which demonstrated a connectivity profile congruent with the histological evidence described above (Bianciardi et al., 2016), to our knowledge, this is the first study to focus specifically on the structural connectivity of the DRN to the limbic system in humans.

We found significant cortical connectivity of the DRN to the amygdala (0.20), DLPFC (0.19), hippocampus (0.12), insula (0.08) and OFC (0.07), as well as to brainstem structures such as the NTS (0.12) and VTA (0.19). There was also relatively greater variation between subjects for the DRN than for any other structure analyzed, congruent with prior work on the serotonin system that has shown significant difference across individuals (Meneses and Liy-Salmeron, 2012; Gold, 2015). Interestingly, given the prior histologic evidence for LC-DRN and PAG-DRN connectivity (Hornung, 2003; Groves and Brown, 2005; Baker and Lui, 2020), we found a relatively low POC between these structures (0.02 and 0.05, respectively), which may reflect a limitation of our method for resolving extremely short-range connectivity within the brainstem (Figure 3). Taken together, our results provide important structural connectivity information, confirming significant connectivity between the DRN, amygdala, and hippocampus.