Loneliness is linked to specific subregional alterations in hippocampus-default network covariation.

Social interaction complexity makes humans unique. But in times of social deprivation, this strength risks exposure of important vulnerabilities. Human social neuroscience studies have placed a premium on the default network (DN). In contrast, hippocampus (HC) subfields have been intensely studied in rodents and monkeys. To bridge these two literatures, we here quantified how DN subregions systematically covary with specific HC subfields in the context of subjective social isolation (i.e., loneliness). By codecomposition using structural brain scans of ∼40,000 UK Biobank participants, loneliness was specially linked to midline subregions in the uncovered DN patterns. These association cortex patterns coincided with concomitant HC patterns implicating especially CA1 and molecular layer. These patterns also showed a strong affiliation with the fornix white matter tract and the nucleus accumbens. In addition, separable signatures of structural HC-DN covariation had distinct associations with the genetic predisposition for loneliness at the population level.NEW & NOTEWORTHY The hippocampus and default network have been implicated in rich social interaction. Yet, these allocortical and neocortical neural systems have been interrogated in mostly separate literatures. Here, we conjointly investigate the hippocampus and default network at a subregion level, by capitalizing structural brain scans from ∼40,000 participants. We thus reveal unique insights on the nature of the "lonely brain" by estimating the regimes of covariation between the hippocampus and default network at population scale.

INTRODUCTION

Social isolation can be detrimental. Accumulating evidence suggests social disconnection as a major risk factor for morbidity and mortality (1, 2). Indeed, it has been found that ∼10–20% of adults who live alone report feeling lonely (3). Yet, it is becoming increasingly clear that experiencing feelings of social isolation, namely loneliness, has diverging biological correlates compared with an objective lack of social contact with others (4). Conversely, individuals who are well surrounded by others can still experience a feeling of inadequate social connection (5). Therefore, it is widely acknowledged that there can be a divide between subjective satisfaction with their social relationships and people’s objective level of social support.

Lonely individuals are particularly characterized by consistent cognitive biases related to processing social cues in the environment. Additionally, loneliness was reported to entail a state of hypervigilance for threats in one’s social environment (5). According to some authors, an aversive feeling of loneliness serves as a biological warning signal that alerts individuals to improve their social relationships (6). Yet, the impact of loneliness extends well beyond the realm of social interaction. For example, previous studies showed that lonely individuals tend to prioritize nonsocial rewards over social rewards (7) and have reduced cognitive control (8), reduced immune response to viruses (6, 9), heightened stress response (10, 11), poorer mental health (1), increased susceptibility to major psychiatric disorders (1), as well as higher risks for alcohol abuse (12), committing suicide (1, 3), cardiovascular disease (13, 14), cognitive decline, and Alzheimer’s disease (15–17). The association with Alzheimer’s disease also motivated our investigation of the hippocampus (HC) and default network (DN) because the hippocampus (18) and DN (19) have long been recognized as a primary neural pathway implicated in the pathophysiology of this major neurodegenerative disease.

Furthermore, while disruption of memory capacity is one of the hallmarks of Alzheimer’s disease (20), the DN and HC have both been implicated in environment-independent processing. This invokes cognitive processes such as episodic memory and mental scene construction, which also show alterations in lonely individuals (21–23). That these two neural systems should relate to loneliness may thus be unsurprising. Indeed, past literature on loneliness has insisted on the importance of rumination on self-focused memories (24), poor executive control (8, 25), and negative bias in the perception of social cues (cf. above). In fact, the aspects of loneliness that rely upon internally generated dimensions of cognition have been argued to especially relate to the DN (26). In rodents and monkeys, the HC has been found to support highly synchronized and spontaneous neural activity bursts during states of rest, so-called sharp-wave ripples (SWRs). These rapid electrophysiological bursts originate in specific anatomical areas of the hippocampus, such as CA1, known to be closely affiliated with DN regions, especially with the medial prefrontal cortex (mPFC) (27). For example, in monkey brains, hippocampal CA1 neurons have been shown to send direct axonal output connections to the mPFC through the fornix white matter pathway (28). This fact could in part explain the recent finding that the fornix is the major white matter tract that is most strongly linked to loneliness (29). The nucleus accumbens (NAc), which is another major target site of hippocampal CA1 neurons and the fornix (28, 30), has likewise been linked to the craving for social connection (31). Alteration of the NAc in loneliness, and thus the reward circuitry of the brain, may therefore be a product of this common anatomical pathway. In light of these details, it has been concluded that since “the pyramidal neurons of the CA1 region provide the only hippocampal output to cortical targets, this activity must have functional significance. We just have to figure out what it is” (32). To take a few steps in this direction, our study aimed to zoom in on the HC-DN correspondence in the context of subjective social isolation.

Indeed, recent evidence has shown that DN gray matter (GM) morphology is strongly associated with loneliness (29). Analogously, the hippocampus has been discussed to be particularly affected by social isolation in various animal species (33–36). Moreover, researchers have been able to perform these invasive studies on the animal hippocampus at a subregion and even single-cell level. Such fine-grained resolution has been difficult to achieve with respect to DN subregions in the living human brain. Since the HC and DN subregions are closely related, studying their structural divergence in lonely individuals provides a window that can offer refined understanding of the association between loneliness and its brain basis. By joint analysis of the HC and DN, the existing in-depth knowledge of the animal HC may help elucidate the nature of potentially human-specific DN subregions. Moreover, in addition to high-resolution structural brain imaging and social isolation information, the availability of genetic data allowed us to also investigate the differing contributions of genetic influences on structural brain patterns related to loneliness.

Previous brain imaging studies aimed at brain parcellation have typically investigated either the default network alone (e.g., Ref. 37) or the hippocampus alone (e.g., Ref. 38). Additionally, there is still insufficient work relating animal studies on anatomically defined subregions of the hippocampus to what can be reliably measured in the human medial temporal lobe with brain imaging techniques. Advances in the automatic segmentation of the hippocampus using ex vivo brain imaging (39, 40) now allow the reliable assessments of microanatomically defined hippocampus subregions in a way that scales to the ∼40,000 UK Biobank Imaging cohort. This enables deeper analyses of the principled interrelationships between the evolutionarily more conserved hippocampus of the allocortex and the default network of the recently expanded neocortex. By leveraging a framework for high-dimensional codecomposition at a fine-grained subregion resolution, we here investigate the structural deviations of the HC-DN covariation signatures that characterize loneliness.

MATERIALS AND METHODS

Population Data Resource

The UK Biobank is a prospective epidemiology resource that offers extensive behavioral and demographic assessments, medical and cognitive measures, as well as biological samples in a cohort of ∼500,000 participants recruited from across Great Britain (https://www.ukbiobank.ac.uk/). This openly accessible population data set aims to provide brain imaging for ∼100,000 individuals, planned for completion in 2022. The present study was based on the recent data release from February/March 2020. To improve comparability and reproducibility, our study built on the uniform data preprocessing pipelines designed and carried out by FMRIB, Oxford University, Oxford, UK (41). Our study involved data from 38,701 participants with brain imaging measures and expert-curated image-derived phenotypes of gray matter morphology (T1-weighted MRI) from 48% men and 52% women, aged 40–69 yr when recruited [mean age 55 yr, standard deviation (SD) 7.5 yr]. The present analyses were conducted under UK Biobank application number 25163. All participants provided written informed consent. Further information on the consent procedure can be found elsewhere (http://biobank.ctsu.ox.ac.uk/crystal/field.cgi?id=200).

Loneliness Target Phenotype

Regarding the loneliness target phenotype, we used the yes/no answer from UK Biobank participants to the question “Do you often feel lonely?” (data field 2020) as a subjective indicator of the quality of social interactions. Measures of social relationship quality represent a widely accepted and commonly investigated component of social embeddedness (42, 43). Loneliness is more commonly viewed as a subjective feeling of being alone, regardless of social encounter frequency (44). Conceptually similar and highly correlated scales (45) are also contained in other standard measurement tools of social embeddedness, such as the Revised UCLA Loneliness Scale (43) and the Interpersonal Support Evaluation List (42). In general, a variety of studies showed single-item measures of social isolation traits to be reliable and valid (46–48). The demographic differences between lonely and nonlonely individuals in the UK Biobank have been reported elsewhere (29).

Brain Imaging and Preprocessing Procedures

Magnetic resonance imaging scanners (3-T Siemens Skyra) were matched at several dedicated data collection sites with the same acquisition protocols and standard Siemens 32-channel radiofrequency receiver head coils. To protect the anonymity of the study participants, brain imaging data were defaced and any sensitive metainformation was removed. Automated processing and quality control pipelines were deployed (41, 49). To improve homogeneity of the imaging data, noise was removed by means of 190 sensitivity features. This approach allowed for the reliable identification and exclusion of problematic brain scans, such as due to excessive head motion.

The structural MRI data were acquired as high-resolution T1-weighted images of brain anatomy with a three-dimensional (3-D) magnetization-prepared rapid gradient echo (MPRAGE) sequence at 1-mm isotropic resolution. Preprocessing included gradient distortion correction (GDC), field of view reduction using the Brain Extraction Tool (50), and FLIRT (51, 52), as well as nonlinear registration to MNI152 standard space at 1-mm resolution using FNIRT (53). To avoid unnecessary interpolation, all image transformations were estimated, combined, and applied by a single interpolation step. Tissue type segmentation into cerebrospinal fluid, gray matter, and white matter was applied by using FAST [FMRIB’s Automated Segmentation Tool (54)] to generate full bias-field-corrected images. SIENAX (55), in turn, was used to derive volumetric measures normalized for head size.

For the default network, volume extraction was anatomically guided by the Schaefer–Yeo reference atlas (56). Among the total of 400 parcels, 91 subregion definitions are provided as belonging to the default network among the seven canonical networks. For the hippocampus, 38 volume measures were extracted with the automatic FreeSurfer subsegmentation (40). The allocortical volumetric segmentation draws on a probabilistic hippocampus atlas with ultrahigh resolution at ∼0.1 mm isotropic. In particular, this tool from the FreeSurfer 7.0 suite gives special attention to surrounding anatomical structures to refine the hippocampus subregion segmentation in each participant.

As a preliminary data-cleaning step, building on previous UK Biobank research (29, 57), interindividual variations in brain region volumes that could be explained by nuisance variables of no interest were regressed out: body mass index, head size, head motion during task-related brain scans, head motion during task-unrelated brain scans, head position and receiver coil in the scanner (x, y, and z), position of scanner table, as well as the data acquisition site, in addition to age, age2, sex, sex × age, and sex × age2. The cleaned volumetric measures from the 91 DN subregions in the neocortex and the 38 HC subregions in the allocortex served as the basis for all subsequent analysis steps.

Structural Variation Relationships of the Fornix and Nucleus Accumbens

In the context of HC-DN correspondence, the fornix white matter tract plays a central role for several important reasons. The fornix constitutes the major output pathway of the hippocampus, carrying information unidirectionally toward the neocortical mantle (28). Our previous research has shown that the fornix is the major white matter tract whose microstructural variation is by far most predictable based on DN subregion volumes alone (58). This fiber pathway has also recently been shown to be linked to loneliness (29). Additionally, on its way, the fornix carries axons from hippocampus neurons to the nucleus accumbens (28, 30)—a core node of the reward circuitry, which has been implicated in the sequelae of loneliness, such as substance abuse (2).

For these reasons, we have conducted preliminary regression analyses for the fornix [UK Biobank (UKB) data fields 25095, 25094, and 25061] and for the nucleus accumbens (UKB data fields 25024 and 25023). The ancillary analyses explored their structural relationships with the subregions composing our DN atlas and those composing our HC atlas (cf. above). The model specification for these regression analyses was as follows: where β denotes the slope parameters corresponding to the brain volumes x for all subregion volumes of either the DN atlas or the HC atlas (z scored across participants) and y denotes the mean fractional anisotropy (FA) of the fornix or the volume of the nucleus accumbens (z scored across participants) for the UK Biobank participants.

Analysis of Covariation between Hippocampus Subregions and Default Network Subregions

As the central step of the analytical workflow, we sought dominant regimes of structural correspondence—signatures or “modes” of population covariation that provide insights into how structural variation among the segregated HC can explain structural variation among the segregated DN. Canonical correlation analysis (CCA) was a natural choice of method to interrogate such a multivariate interrelation between two high-dimensional variable sets (59–61). A first variable set X was constructed from the DN subregion volumes (number of participants × 91 DN parcels matrix). A second variable set Y was constructed from the HC subregion volumes (number of participants × 38 HC parcels matrix): where n denotes the number of observations or participants, p is the number of DN subregions, and q is the number of HC subregions. Each column of the two data matrices was z scored to zero mean (i.e., centering) and unit variance (i.e., rescaling). CCA addresses the problem of maximizing the linear correlation between low-rank projections from the two variable sets or data matrices. The two sets of linear combinations of the original variables are obtained as follows: where V and U denote the respective contributions of X and Y, L and L denote the respective latent “modes” of joint variation based on patterns derived from X and patterns derived from Y, l is the lth column of L, and l is the lth column of L. We define modes as general principles of population variation in our target neural circuits that can be reliably extracted in brain structure at the population level. The goal of our CCA application was to find pairs of latent vectors l and l with maximal correlation in the derived latent embedding. In an iterative process, the data matrices X and Y were decomposed into L components, where L denotes the number of modes given the model specification. In other words, CCA involves finding the canonical vectors u and v that maximize the (symmetric) relationship between a linear combination of DN volume expressions (X) and a linear combination of HC volume expressions (Y). CCA thus identifies the two concomitant projections Xv and Yu. These yielded the optimized co-occurrence between patterns of subregion variation inside the segregated DN and patterns of subregion variation inside the segregated HC across participants.

In other words, each estimated cross-correlation signature identified a constellation of within-DN volumetric variation and a constellation of within-HC volumetric variation that go hand in hand with each other. The set of k orthogonal modes of population covariation are mutually uncorrelated by construction (60). They are also naturally ordered from the most important to the least important HC-DN covariation mode based on the amount of variance explained between the allocortical and neocortical variable sets. The first and strongest mode therefore explained the largest fraction of joint variation between combinations of HC subregions and combinations of DN subregions. Each ensuing cross-correlation signature captured a fraction of structural variation that is not explained by one of the k − 1 other modes. The variable sets were entered into CCA after a confound-removal procedure based on previous UK Biobank research (cf. above).

Group Difference Analysis

For the derived population modes of HC-DN covariation, we then performed a rigorous group contrast analysis in the context of loneliness. We aimed to identify which anatomical subregions show statistically defensible deviation in the lonely group compared with the control group. We carried out a principled test for whether the solution vector obtained from the CCA (i.e., canonical vectors, cf. above) in the lonely group is systematically different from the solution vector in the control group.

More specifically, following previous UK Biobank research (57), we carried out a bootstrap difference test of the CCA solution from the lonely versus nonlonely groups (62). In 100 bootstrap iterations, we randomly pulled participant samples with replacement to build an alternative data set (with the same sample size) that we could have gotten. We subsequently performed CCA in parallel by fitting one separate model to each of the two groups. In each resampling iteration, this approach thus carried out a separate estimation of the doubly multivariate correspondence between HC subregions and DN subregions in each particular group. The two distinct CCA solutions from each iteration were then matched mode by mode regarding sign invariance and mode order. Canonical vectors of a given mode that carried opposite signs were aligned by multiplying one with −1. The order of the CCA modes was aligned based on pairwise Pearson’s correlation coefficient between the canonical vectors from each estimated CCA model. After mode matching, we directly estimated the resample-to-resample effects of model parameter estimates by elementwise subtraction of the corresponding canonical vectors of a given mode k between the two groups. We finally recorded these difference estimates for each vector entry (each corresponding to the degree of deviation in 1 particular anatomical subregion). The subregion-wise differences were ultimately aggregated across the 100 bootstrap data sets to obtain a nonparametric distribution of group contrast estimates.

We thus propagated the variability attributable to participant sampling into the computed uncertainty estimates of group differences in the UK Biobank population cohort. Statistically relevant alterations of anatomical subregions in loneliness were determined by whether the two-sided confidence interval included zero or not according to the 10/90% bootstrap-derived distribution of difference estimates (57) in an approach that is faithful to our multivariate analytical strategy and research question at hand. This nonparametric approach directly quantified the statistical uncertainty of how loneliness is manifested in specific subregions of the HC-DN axis.

Analysis of How Individual Expressions of Hippocampus-DN Covariation Are Linked to the Genetic Predisposition for Loneliness

Polygenic risk score (PRS) is a genome-wide analysis technique that has been shown to successfully quantify individuals’ genetic predisposition for a variety of phenotypes. The approach has become especially potent for complex phenotypes that implicate tens of thousands of single-nucleotide polymorphisms (SNPs) with individually small effect sizes, such as major psychiatric diseases (63–67). PRSs have also come to be a sharp tool for heritability analyses because of the advent of large population data sets such as the UK Biobank (68, 69). Such data resources have allowed the investigation of the relationship between SNP variation and interindividual differences in a particular phenotype, which includes neuroimaging-derived phenotypes (70, 71). For the purpose of the present study, we have constructed PRS models for the loneliness trait. The subject-specific risk scores were then regressed onto our expressions of HC-DN modes (i.e., canonical variates). The integrative imaging-genetics approach allowed disentangling which mode expressions showed reliable relationships to the heritability of loneliness as attributable to thousands of genetic variants.

As is common for PRS analysis workflows, summary statistics from previously conducted genome-wide association studies (GWASs) on our target phenotypes were used as the backdrop to determine how several hundred thousand SNPs are associated with the loneliness trait. The summary statistics for loneliness were obtained from a GWAS that was conducted as part of the Psychiatric Genomics Consortium. Quality control was implemented by excluding SNPs with a minor allele frequency of <1%, as well as excluding SNPs with imputation information score of <0.8. Mismatching, duplicate, and ambiguous SNPs were also disregarded from further analysis. Quality control on the base data also involved excluding individuals with a difference between reported sex and that indicated by their sex chromosomes and removing overlapping samples.

The quality-controlled summary statistics were used as starting point for the PRS model that was built and applied with the PRSice framework (http://www.prsice.info). This software tool uses the available collection of effect sizes of candidate SNPs to form single-subject predictions of genetic predisposition for a phenotype of interest. In particular, this tool determined the optimal PRS model based on the UK Biobank participants (training data, n = 253,295) of European ancestry who did not provide any brain imaging data (at the time of study). This model training step involved automated adjustments, such as identifying ideal clumping and pruning choices, to select the thresholds that decide which SNPs are included in the PRS model. Subsequently, once optimized, the final PRS model was then used to predict the genetic predisposition for each of 23,423 UK Biobank participants of European ancestry with brain imaging data (test data). This PRS model consisted of the additive effects of weighted SNPs, whereby the weighted sum of the participants’ genotypes was computed as follows: where g denotes an individual’s genotype at SNP i (value 0, 1, or 2), is the obtained point estimates of the per-allele effect sizes at SNP i, and j is a particular individual (68).

Finally, Bayesian linear regression was used to regress the subject-specific predictions of genetic liability for the loneliness trait onto the participant expressions of HC-DN covariation modes. More specifically, the individuals in the top 5% predictions (i.e., highest PRS estimates) and the individuals in the bottom 5% predictions (i.e., lowest PRS estimates) were considered in a Bayesian logistic regression model with mode expressions serving as input variables (72–74). In this multiple regression setup, PRS was regressed against each of the 25 canonical variates (linearly uncorrelated by construction) for each individual on the hippocampus side. An analogous multiple regression model was estimated for the (uncorrelated) 25 canonical variates from the DN side. The fully Bayesian model specification for these regression analyses was as follows: where β denote the slopes for the subject-specific 25 mode expressions as canonical variates x; yprs denotes the PRS estimates of each participant. Potential confounding influences were acknowledged by the nuisance variables α, which accounted for variation that could be explained by sex and (z scored) age. Once the Bayesian model solution was approximated with Markov chain Monte Carlo sampling, it yielded fully probabilistically specified posterior parameter distributions for each β coefficient corresponding to one of the signatures of allocortical-neocortical covariation (cf. above). The association with trait heritability of a mode expression was then determined based on how robustly their corresponding model coefficients deviated from 0 (e.g., >95% of model coefficient posterior probability excluded a value of 0) (173, 174).

RESULTS

Previous brain imaging studies aimed at brain parcellation have typically investigated either the default network alone (e.g., Ref. 37) or the hippocampus alone (e.g., Ref. 38). Additionally, there is still insufficient work relating animal studies on anatomically defined subregions of the hippocampus to what can be reliably measured in the human medial temporal lobe with brain imaging techniques. Advances in the automatic segmentation of the hippocampus using ex vivo brain imaging (39, 40) now allow the reliable assessments of microanatomically defined hippocampus subregions in a way that scales to the ∼40,000 UK Biobank Imaging cohort. This enables deeper analyses of the principled interrelationships between the evolutionarily more conserved hippocampus of the allocortex and the default network of the recently expanded neocortex. By leveraging a framework for high-dimensional codecomposition at a fine-grained subregion resolution, we here investigate the structural deviations of the HC-DN covariation signatures that characterize loneliness.

Structural Variation Relationships of the Fornix

In a preliminary set of exploratory analyses, we sought to elucidate whether individual subregion volumes of the hippocampus could distinctly explain microstructural variation of the fornix white matter tract—an important fiber pathway known to carry axons from the allocortical hippocampus to the DN. To this end, we performed a multiple regression analysis to regress the structural integrity of the fornix (mean FA, diffusion MRI) onto 38 hippocampal subregion volumes. We found that the bilateral CA1 body dominates in explaining variation in fornix integrity, followed by the molecular layer (ML) head of the left hippocampus (Fig. 1). Smaller positive effects were found in the bilateral subiculum body, bilateral presubiculum body, bilateral fimbria, right ML head, and right CA2/3. We found negative effects in the bilateral hippocampal fissure, left parasubiculum, bilateral CA4 head, as well as bilateral CA4 body. These findings indicate that several hippocampal subregions have unique structural relationships with fornix integrity. In a variant of this analysis, rather than marginal or partial association strengths we computed pairwise Pearson correlation strengths between fornix integrity and each hippocampal subregion. The hippocampus subregions with the largest correlation coefficient rho were also the subregions with the greatest amount of explained variance in the multiple regression analysis (Supplemental Fig. S1; all Supplemental Material can be found at https://doi.org/10.6084/m9.figshare.15060684). Thus, from convergent evidence across two different analyses of how hippocampus subregions track fornix variation, the CA1 body showed the strongest positive association with the fornix, whereas the hippocampal fissure displayed the strongest negative association.

Figure 1.

Structural correlates of the fornix white matter tract. How can microstructural variation of the fornix be explained as a function of subregion volume variations in the hippocampus (left) or the default network (DN) (right)? To address these two questions, we have computed separate exploratory multiple regression analyses based on 38,701 UK Biobank participants. The regression parameter weights corresponding to each specific subregion (hot and cold colors) outline the structural associations with fornix integrity exclusively explained by that particular subregion. A shows the hippocampus-fornix model based on 38 subregions from the hippocampus segmentation (40). The parameter weights indicate the variation explained in the fornix specifically by each hippocampal subregion, mapped onto 8 consecutive coronal sections of the left and right hippocampus from anterior (top) to posterior (bottom) (hot/cold colors = positive/negative volume association). The subregions with the strongest positive volume effect are CA1 body of the left and right hippocampus as well as left molecular layer head. The strongest negative relations with fornix microstructure are found in left and right hippocampal fissures. B shows the DN-fornix model based on 91 subregions in the Schaefer–Yeo parcellation of the DN. The parameter weights indicate the fornix variation explained specifically by each DN subregion. The subregions with the strongest marginal volume effect among DN subregions are bilateral posterior cingulate cortex (PCC), bilateral retrosplenial cortex (RSC), and left dorsal-medial prefrontal cortex (dmPFC). Overall, these results suggest that fornix architecture has differentiated structural relationships with subregions of the hippocampus and midline DN. DG, granule cell layer of dentate gyrus (GC-DG-ML); ML, molecular layer; OFC, orbitofrontal cortex; Para, parasubiculum; PrS, presubiculum; Sub, subiculum, vmPFC, ventromedial prefrontal cortex.

We next sought to regress variation of fornix integrity onto 91 subregion volumes of the DN (Fig. 1). The subregions with the strongest positive effects included the left dorsal mPFC (dmPFC) and the right retrosplenial cortex (RSC). There was also a small effect in the left RSC, bilateral ventromedial prefrontal cortex (vmPFC), and left dorsolateral prefrontal cortex (dlPFC). Conversely, we located the strongest negative effects to the bilateral posterior cingulate cortex (PCC). Furthermore, in the companion analysis using Pearson correlation between each DN subregion and the fornix, mPFC subregions had in general strong positive effects (Supplemental Fig. S2). On the other hand, PCC/RSC subregions had weak or negative effects (Supplemental Fig. S2). These findings thus made apparent the subregion-specific relationships of how DN subregions track fornix integrity in our UKB cohort, with lateral DN subregions showing more subtle effects and midline DN subregions showing the most prominent effects.

Structural Variation Relationships of the Nucleus Accumbens

Using the identical analysis approach (cf. above), we next explored how the volume of the nucleus accumbens (NAc) can be explained as a function of hippocampal subregion volumes. These multiple regression analyses revealed that bilateral CA1 body and right molecular layer head subregions explained the most variation in NAc volume (Supplemental Fig. S2A). Small positive effects were also found in CA2/3 head of the right hippocampus, left molecular layer head, and bilateral hippocampus amygdala transition area (HATA). Conversely, there was a strong negative association of the NAc with bilateral hippocampal fissure as well as with left parasubiculum. Interrogating the HC-NAc relationships thus showed distinct explanatory effects among subregions. In particular, the hippocampal subregions that explained the most variation in NAc recapitulated those that also best explained fornix integrity (cf. above). To further investigate the associations of the NAc, we next computed gross pairwise Pearson correlation coefficients of each hippocampal subregion with the NAc. We observed a noticeably stronger correlation for CA1 head than any other subregion (Supplemental Fig. S3), although CA1 head showed a negligible effect in the multiple regression analysis of the NAc on hippocampus subregions. Thus, although CA1 head on its own has a strong correlation with the NAc, when considered in a joint model with all other hippocampal subregions in our atlas it explained little NAc variation in a unique way. This observation suggests that distinct hippocampal subregions, other than CA1, make separate contributions that together better explain volume variation in the NAc.

To complement these hippocampus analyses, the relationship between the NAc and DN subregions was also investigated with a dedicated multiple regression analysis. The subregions with the largest explained variance in NAc volume included right orbitofrontal cortex (OFC), left OFC, left posterior superior temporal sulcus (pSTS), left ventrolateral prefrontal cortex (vlPFC), and left vmPFC (Fig. 2). Conversely, strong negative associations with NAc volume were located to the left and right vmPFC, left and right PCC, as well as right temporal pole. In a Pearson correlation analysis of each DN subregion with the NAc, we also observed that mPFC subregions had generally strong correlations with the NAc and PCC/RSC subregions had generally weak correlations (Supplemental Fig. S4). Overall, in contrast to our HC-centric analyses for the fornix and NAc (cf. above), DN subregions were found to have individually varying relationships with the NAc that noticeably drew a different picture than our results obtained for the fornix-DN relationship.

Figure 2.

Structural correlates of the nucleus accumbens (NAc). Regression parameter weights of various subregions in separate exploratory multiple regression analyses outline the volume effects on the NAc, due exclusively to each subregion in either the hippocampus or default network. A shows the hippocampus (HC)-NAc model mapped onto 8 consecutive coronal slices of the left and right HC in the anterior (top)-to-posterior (bottom) direction. The parameter weights indicate the structural variation uniquely explained in the NAc specifically by each hippocampal subregion (hot/cold colors = positive/negative volume association). The subregions with the strongest positive volume effect are CA1 body of the left and right hippocampus as well as right molecular layer (ML) head, in line with our analogous analyses on the fornix (Fig. 1). The strongest negative volume effects are left and right hippocampal fissures. B shows the default network-NAc model. The parameter weights indicate the variation explained in the NAc specifically by each default network subregion (warm/cold colors = positive/negative volume association). The subregions with the strongest volume effect among default network subregions are bilateral orbitofrontal cortex (OFC), left posterior superior temporal sulcus (pSTS), and left ventrolateral prefrontal cortex (vlPFC). Overall, these results suggest there are diverse structural relationships between the NAc and hippocampal subregions, which are substantially similar to the structural relationships of the fornix (Fig. 1). Yet, the relationship of the NAc with default network subregions shows substantially different volume effect sizes for the most prominent cortical regions compared with the fornix. DG, granule cell layer of dentate gyrus (GC-DG-ML; HATA), hippocampal amygdala transition area; ML, molecular layer; Para, parasubiculum; PrS, presubiculum; vmPFC, ventromedial prefrontal cortex.

Structural Covariation between Hippocampus and Default Network at the Subregion Level

In our core analysis, we explored the principled signatures of structural covariation between the full set of 38 hippocampal subregions and the full set of 91 DN subregions. The concurrent patterns of subregion variation within the hippocampus and within the DN were computed using doubly multivariate pattern-learning analysis. In so doing, we achieved a codecomposition of hippocampal subregion volumes and DN subregion volumes. Each of the top 25 modes of covariation was characterized based on how much of joint variance a particular signature explained, with the most explanatory signature (mode 1) achieving a canonical correlation of rho = 0.51 (measured as Pearson’s correlation coefficient) (Supplemental Table S1). The second most explanatory signature (mode 2) achieved a canonical correlation of rho = 0.42, the third rho = 0.39, the fourth rho = 0.31, the fifth rho = 0.27, and the sixth rho = 0.23, through to the twenty-fifth signature, which had rho = 0.06 (see Supplemental Table S1 for full list). This analysis thus established the scaffold for all subsequent analyses that delineates how multiple complementary hippocampal patterns covary hand in hand with DN patterns.

Differences in the Hippocampus-Default Network Covariation in Loneliness

On the basis of the identified population signatures of HC-DN covariation, we investigated the neurobiological manifestations of loneliness in our UKB sample. This was accomplished by examining robust subregion-level divergences in how hippocampal patterns are coexpressed with DN patterns between groups of lonely versus nonlonely participants. To this end, we first analyzed loneliness by a rigorous group difference analysis between the structural patterns of covariation in the lonely and nonlonely groups. This revealed the precise subregions contributing to the structural HC-DN covariation that systematically diverged between the two groups, for each mode of the CCA.

We uncovered a multitude of modes with systematic group differences in either specific HC and/or DN subregions. We also found modes with no significant structural divergences in any hippocampal or DN subregion. From here on, a subregion that was observed to have a robustly different weighting within a mode’s canonical vector, between the lonely and nonlonely groups, is termed a “hit” (i.e., an observed structural divergence in lonely individuals). Across all 25 examined modes, the group contrasts amounted to 28 total hits for HC subregions and 40 total hits for DN regions. Most of these subregion hits occurred in earlier modes. For example, in the first mode we found 26 DN hits (60% of the DN total), and in the first three modes we found 22 HC subregion hits (78.6% of the HC total). We also noted specific subregions with repeated hits across modes: The largest number of hits was in molecular layer body (4 hits), molecular layer head (3), CA1 head (3), CA1 body (3), and presubiculum head (3). In the DN, the greatest number of hits was observed in midline subregions (77.5% of the DN total) such as the RSC and mPFC. We observed few hits in lateral temporal or parietal subregions, with a total number of hits in 5 temporal (12.5%), 13 PFC (32.5%), 2 parietal (5%), and 20 PCC/RSC (50%) subregions.

The constellation of structural divergences between lonely and nonlonely groups in mode 1 (Fig. 3) provided a rough portrait of the hits that typically emerged across the next 24 modes. In mode 1 we observed hits in bilateral CA1 head and left CA1 body, the adjacent left presubiculum head and body, as well as the bilateral subregions internal to the hippocampus (ML and fissure). Concurrently, the DN hits for mode 1 were clustered in adjacent subregions in the bilateral PCC/RSC region and mPFC (2 temporal, 8 PFC, 1 parietal, 15 PCC). Additionally, the subregions with hits that played an especially strong role in the dominant mode’s patterns were bilateral CA1 body, bilateral fissure, bilateral RSC, and left dmPFC. In contrast, for the second mode of the HC-DN covariation, we observed hits only in the left parasubiculum and left ML head, with no hits in DN regions (Supplemental Fig. S5). For mode 3, we observed eight HC hits and no DN hits (Fig. 4); for mode 4, one hit in the HC tail and no DN hits for mode 5, one HC and one DN hit (Fig. 5); for mode 6, no hits; for mode 7, no HC hits and one hit in the right temporoparietal junction for mode 8, HC hits in three subregions, the left CA1 head, left CA2/3 head, and left granule cell layer of dentate gyrus (GC-ML-DG) head, and no DN hits (Supplemental Fig. S6); and for mode 9, one HC hit and seven DN hits (Fig. 6). We observed no hits in modes 10–13. In mode 14 we found no HC hits and three DN hits: left superior temporal sulcus (STS), right anterior cingulate cortex (ACC), and left RSC (Supplemental Fig. S7). In mode 15, we found no HC hits and two hits in DN subregions: bilateral PCC (Supplemental Fig. S8). Of note, we observed no hits in any of the modes between mode 16 and mode 25. Thus the dominant modes, which account for more population covariance, are most strongly coupled with loneliness.

Figure 3.

Loneliness is associated with distinct divergences in hippocampus (HC)-default network (DN) covariation. We explored the structural covariation between the 38 subregions of the HC and 91 subregions of the DN, by means of a codecomposition based on a canonical correlation analysis (CCA). We subsequently determined how the ensuing subregion patterns diverged in individuals with loneliness. Shown here are the subregion divergences in mode 1 of HC-DN covariation. Mode 1 of the CCA solution achieves the most explanatory covariation, with a canonical correlation of rho = 0.51. A shows the HC subregion pattern (left, one canonical vector of mode 1) with parameter weights that robustly diverge between lonely and nonlonely groups in mode 1; mapped onto 8 consecutive coronal slices of the left and right HC in the anterior (top)-to-posterior (bottom) direction. B shows the DN subregion pattern (right, other canonical vector of mode 1) that robustly diverged between the lonely and nonlonely groups. Overall, within the dominant structural covariation pattern between the HC and the DN there are specific subregions whose volumes systematically diverge in lonely individuals. The most pronounced structural divergences are in bilateral CA1 body and hippocampal fissure, as well as subregions bilaterally in the posterior cingulate cortex (PCC), retrosplenial cortex (RSC), and dorsomedial prefrontal cortex (dmPFC). Thus, specific HC and DN anatomical subregions are preferentially linked loneliness. ML, molecular layer; PrS, presubiculum; STS, superior temporal sulcus; Sub, subiculum; TPJ, temporoparietal junction.

Figure 4.

Divergences in CA1, subiculum (Sub), presubiculum (PrS), and molecular layer (ML) subregion volumes are preferentially associated with loneliness. Shown are the subregion divergences in mode 3 of hippocampus (HC)-default network (DN) covariation. Mode 3 of the canonical correlation analysis (CCA) solution achieves the third most explanatory hippocampal-DN covariation, with a canonical correlation of rho = 0.39. A exhibits the (HC) subregion pattern (left, one canonical vector of mode 3) with paramater weights that robustly diverge between lonely and nonlonely groups in mode 3. B shows the DN subregion pattern (right, other canonical vector of mode 3) of the DN (right) that robustly diverge between lonely and nonlonely groups. The most pronounced divergences in the loneliness group are in right ML head and body, right Sub head and body, as well as right CA1 head. Thus, similar to mode 1, within the third most explanatory pattern of HC-DN covariation, loneliness is concomitant with hippocampal divergences in CA1, Sub, PrS, and ML.

Figure 5.

Loneliness is associated with divergences in left hemisphere CA1 and retrosplenial cortex (RSC). Shown are the subregion divergences in mode 5 of hippocampus (HC)-default network (DN) covariation. Mode 5 of the canonical correlation analysis (CCA) solution achieves the fifth most explanatory hippocampal-DN covariation, with a canonical correlation of rho = 0.27. A exhibits the HC subregion pattern (left, one canonical vector of mode 5) with parameter weights that robustly diverge between lonely and nonlonely groups in mode 5. B shows the DN subregion pattern (right, other canonical vector of mode 5) of the DN (right) that robustly diverge between lonely and nonlonely groups. The only subregion divergences associated with loneliness are in left CA1 body and the left RSC. These results accentuate the selectivity of the subregion divergences within a particular mode and highlight CA1 and RSC in loneliness.

Figure 6.

Loneliness is associated with subregion divergences in parasubiculum (Para) and a distributed set of bilateral default network (DN) subregions. Shown here are the subregion divergences in mode 9 of hippocampus (HC)-DN covariation. Mode 9 of the CCA solution achieves the ninth most explanatory HC-DN covariation, with a canonical correlation of rho = 0.18. A exhibits the HC subregion pattern (left, one canonical vector of mode 9) with parameter weights that robustly diverge between lonely and nonlonely groups in mode 9. B shows the DN subregion pattern (right, other canonical vector of mode 9) of the DN (right) that robustly diverged between lonely and nonlonely groups. Overall, there are divergences for lonely individuals in left Para, as well as bilateral ventromedial prefrontal cortex (vmPFC), left retrosplenial cortex (RSC), bilateral orbitofrontal cortex (OFC), and bilateral superior temporal sulcus (STS). These results emphasize the bilateral tendency of the DN divergences observed for loneliness. The results also highlight a coherent set of reward-related subregions, as most of the subregions identified are also found to have a strong volumetric relationship with the nucleus accumbens (Fig. 2).

These results suggest that specific subregions that play a role in the HC-DN correspondence exhibit systematic diverges in individuals with loneliness. Additionally, our findings specify within which context of structural dependence (mode of covariation) a subregion diverged between the two groups. Overall, across loneliness’s relationships to our discovered population modes of joint volume variation, CA1 head and body, molecular layer head and body, and presubiculum head showed the largest number of hits for hippocampal subregions. On the other hand, there was a preponderance of hits located to midline cortical structures, such as the mPFC and PCC/RSC, for DN subregions. Thus, our population-level findings made apparent that within structural covariation modes, loneliness was mostly concomitant with divergences in DN midline subregions that were grounded in parallel divergences in the CA1 and ML subregions of our hippocampus atlas.

Genetic Predisposition for Loneliness

We finally sought to interrogate whether the uncovered HC-DN covariation expressions featured systemic relationships with the participants’ liability for loneliness (cf. materials and methods). For this purpose, we first computed polygenic risk score predictions of loneliness risk for our UK Biobank participants. We observed that there was a statistically relevant relationship between loneliness PRS score and participant expressions (i.e., canonical variates) of modes 7, 8, 9, and 22 in the hippocampus and the participant expressions of mode 7 in the DN (Fig. 7). On the whole, there was genetic predisposition linked to the expressions of later modes compared with the earlier modes achieving greater explained variance. This analysis suggested that the expression of the identified HC-DN signatures tracks the purely heritable components of loneliness due to single-nucleotide polymorphism. As a matter of course, the signatures not identified likely have contributions that are not due to genetic factors.

Figure 7.

The genetic predisposition for loneliness is associated with the expression of specific hippocampus(HC)-default network (DN) covariation patterns. A polygenic risk score (PRS) analysis was conducted to estimate the subject-specific heritable tendency for loneliness based on genome-wide effects in single-nucleotide polymorphisms. The subject-specific PRS estimates were then regressed against the expressions of each of the modes of structural covariation between the HC and DN. The relevance of the heritability effects was judged based on the posterior parameter distributions inferred by the Bayesian logistic regression model (histograms). The x-axis of each plot represents the β coefficient of the model parameter, and the y-axis of each plot represents the plausibility of each coefficient value. A value of 0 indicates no association between PRS and interindividual mode expression. Overall, loneliness predisposition was related to the expression of modes 7 and 9 on the HC side and to mode 7 on the DN network side. These results demonstrate that specific mode expressions are robustly linked to the predisposition of becoming lonely. The three posterior distribution histograms show convergence across three different Markov chain Monte Carlo runs. Overall, this pinpoints the specific neurobiological signatures in the HC and DN that relate to the heritable components of loneliness.

DISCUSSION

We have tailored an analytical framework to examine how the structural substrates of the HC-DN correspondence systematically deviate in loneliness. We work toward this goal by directly estimating principled covariation signatures that delineate how hippocampus subregion volumes are coexpressed with DN subregion volumes in ∼40,000 participants from the UK Biobank cohort. In so doing, our study aimed to refine the understanding of the interrelationship of the DN with loneliness. By anchoring our analysis in HC-DN covariation, we aimed to facilitate the interpretation of potentially human-specific DN brain regions based on their hippocampal subregion affiliates. Indeed, the evolutionarily conserved hippocampus has been extensively studied in animal species at a single-cell resolution. This approach can detect patterns in the structural brain scans that are not discoverable by any analysis of the HC or the DN alone.

The DN has recently been found to be the major brain network that is most closely associated with loneliness (29). However, we found a high degree of heterogeneity in the relationship of individual DN subregions with loneliness in the context of population modes of HC-DN volume covariation. In particular, we found that especially midline structures of the DN tended to diverge in lonely individuals. Analogously, the hippocampal subregions with the most divergences across modes were CA1 and molecular layer. The allocortical-neocortical divergences implicating midline DN as well as the CA1 and molecular layer of the HC stood out, especially in modes 1, 3, and 5 of our analysis. The concomitant structural divergences of HC and DN subregions in lonely individuals may in part be explained by the varieties of internally directed cognition that have been consistently associated with each neural system—such as episodic memory processing, sense of self, and prospective cognition (21, 75–77). Especially midline structures of the DN such as the vmPFC (78), posteromedial cortex (79, 80), and the RSC (75) have been found to serve these types of neural processes. In line with this notion, the fornix—known to be a direct mediator between the HC and midline DN subregions—has been found to be the white matter structure that is most strongly linked to loneliness profiles (29). This observation is of relevance to our present considerations, as direct neuron-to-neuron connections between CA1 and subiculum of the hippocampus with mPFC, OFC, and RSC subregions of the DN have been demonstrated through invasive axonal tracing studies in monkeys (28, 81–84). Human studies have also found that the functional connectivity strength between the hippocampus and mPFC correlates with the imagination of mental scenes (78, 85, 86) and precommissural fornix integrity correlates with episodic richness (87, 88). Indeed, these earlier findings may explain various aspects of the psychological alterations that are characteristic for lonely individuals.

For example, lonely individuals have been described to more frequently reiterate social events from the past, imagine hypothetical or future encounters with others, and reminisce on nostalgic memories (89–93). The joint divergences of CA1 and midline DN subregions thus suggest an intimate link between these spatially distributed subregions, a link that 1) is supported by the fornix, 2) implicates internally directed cognition, and 3) is systematically altered in lonely individuals. The notion of CA1, molecular layer, OFC, and mPFC subregions having preferential structural relationships with the fornix receives additional support from our multiple regression and pairwise Pearson’s correlation analyses. Indeed, each of these subregions showed especially strong relationships with the fornix in both analyses. In accord with the growing knowledge of the fornix, it is possible that these concordant findings reflect a difference in the neurocognitive processes necessary for episodic memory retrieval (88, 94). For example, the CA1 and subiculum, which are the primary origin of precommissural fornix projection fibers to the PFC and NAc (28, 30, 95), have been found to be crucial for proto-episodic memory in rodents (27, 96) and detailed episodic autobiographical memory in humans (97, 98). Hence, the known cognitive upregulation in memory retrieval and social perception biases of lonely individuals may be rooted in how the fornix pathway assists HC-DN communication.

The cognitive biases of the loneliness trait further include perceiving one’s social world as a more threatening place (4, 99), such as remembering more negative information from past social encounters and expecting more negative social interactions (5). On the basis of these earlier psychological insights, we expected that the subregions with multiple divergences in lonely individuals would show a strong relationship to the NAc, a key node of the brain’s reward circuitry that has previously been associated with loneliness (31). In support of this idea, our multiple regression analysis revealed that CA1, molecular layer, vmPFC, and OFC each showed strong relationships with the NAc. The seemingly separate psychological and behavioral characteristics of loneliness may thus be a product of the remodeling of integrated neural systems. This notion receives further support from the similarity in the subregion associations revealed through our separate multiple regression models of the hippocampus. In these, we observed a strong overlap between the hippocampal subregions having strong relationships with either the NAc or fornix. Indeed, the NAc receives the bulk of its hippocampal inputs through the precommissural fornix pathway (83, 100). Hence, our results substantiate previous findings on cross-associations between the CA1 of the hippocampus, fornix, mPFC, and NAc.

In addition to the collective divergences in hippocampal and mPFC subregions in lonely individuals, our analysis pinpointed parallel findings in other DN subregions. In particular, the majority of DN subregion divergences we found were in the PCC, and especially RSC. In general, these subregions also showed strong structural relationships with the fornix in our regression analyses. In line with our results from brain imaging modeling, direct axon tracing studies in rodents have documented CA1 and subiculum as the main hippocampal projection sites to the RSC (101, 102). Monkey experiments have additionally reported that direct HC-RSC projections are chiefly mediated through the subiculum, a close interaction partner of hippocampal CA1 (84, 102–104). In conjunction with the CA1 structural divergences we observed in lonely individuals, we also observed four total subiculum divergences across all examined signatures. Hence, our population neuroscience evidence for interdependence between CA1 and subiculum in the allocortex and RSC subregions in the neocortex directly confirms established knowledge from invasive animal experiments (105). Thus, an attractive interpretation of these findings of structural divergences in RSC subregions for lonely individuals (in modes 1, 5, 9, 14, and 15) invokes a difference in the capacity to visualize internally centered thoughts. Indeed, among the limited knowledge regarding the function of the human RSC, this part of the posterior midline is believed to assist in spatial navigation, episodic memory, visualizing details of mental scene construction, and perhaps mentalizing capability (21, 26, 106–108).

Moreover, RSC-mediated mental scene construction could be one aspect of a more general role the DN has in spontaneous cognition, such as task-free “mind-wandering” (109, 110). In fact, the hippocampus has been found to support stimulus-independent neural activity (111), which would be a perfect partner for these DN functions, via so-called sharp-wave ripples (SWRs) (112, 113). These phenomena are manifestations of spontaneous neural activity synchronizations in the hippocampus circuit that occur independent of any environmental cues or externally structured tasks (114, 115). SWRs have been carefully studied through electrophysiological experiments in rodents and other animals, as opposed to most of the neural activity specific to regions of the human DN (27, 116). Of relevance for our present considerations, CA1 pyramidal neurons are purported to underlie the expansive “ripples” of SWR activity. The preceding “sharp wave,” on the other hand, is most pronounced in the apical dendrites of CA1 neurons (27). Along this line, the hippocampus atlas used for our study collapsed the apical dendrites (both the stratum radiatum and stratum lacunosum moleculare) of subiculum and CA neurons into one composite subregion: the molecular layer (40). Thus, the molecular layer and CA1 of our set of human hippocampus subregions offered a way to analyze systematic alterations of the neural substrates shown to be highly visible locations of SWR activity in the hippocampus of rodents and monkeys.

Indeed, the CA1 and molecular layer had the greatest number of total hits in lonely individuals. More specifically, there were seven total hits in the molecular layer and six total hits in CA1 across all examined modes. Since our analyses revealed co-occurring divergences in midline DN subregions, SWR activity suggests itself as a compelling interpretation that can accommodate these collective results. The sensory-distal processing regimes realized by the DN may be perfectly attuned to integrating spontaneous SWR activity. For example, in rodents, CA1 cell activity and coupled SWR have been shown to assist in prospective cognitive processes (117), such as simulating spatial trajectories that were never experienced before (118–120). A functional interaction between hippocampal SWR activity and association cortex regions, such as the RSC, has even been demonstrated in animal species (105, 121). For example, in rats it has been shown that auditory cortex activity predicts subsequent SWR activity (122). Conversely, hippocampal activity also predicted subsequent activity in the auditory cortex. On the basis of these results, the authors proposed that immediately before and after SWR occur there is rapid information flow in a cortical-hippocampal-cortical loop (122). Similarly, studies in primates (123, 124) have also identified the RSC and PCC to be the regions of the neocortical mantle that are most closely coupled with hippocampal SWR (124). Our findings of selective structural divergences in both the CA1 and molecular layer as well as RSC subregions across modes is thus very suggestive of a difference in the allocortical-neocortical information processing pathways that involve SWR activity.

Facilitation of mental scene construction is discussed as one of the main functional contributions of the hippocampus (21, 109). It thus appears plausible that a systematic alteration of SWR activity may affect the vividness of an individual’s episodic memories and imaginings (125). In fact, SWR activity has been proposed to underlie trajectory sequence replay in rodents as well as “offline” states of human cognition, such as daydreaming (32). Thus, the increased frequency and intensity of nostalgic reminiscence that has been associated with loneliness may potentially also be associated with altered SWR activity (91). A strong relationship between SWR and loneliness is additionally suggested by their shared association with reward processing. Indeed, loneliness alters the subjectively perceived valence of social and nonsocial reward in humans (7), while SWR activity has been demonstrated to be responsive to reward-contingent processing in animals (126–128). Similarly, the fundamental place cell representations of the CA1 pyramidal cells, which underlie SWR activity (cf. above), were reported to intimately relate to reward processing (129, 130). Certain populations of CA1 cells have even been found to act as a dedicated channel for reward-related processing across environments (131). These findings and conclusions from previous animal experiments therefore appear to dovetail with our analyses in people highlighting the NAc with CA1 and molecular layer. We hence wonder whether the attempt to fill the perceived social void through mental imagery in lonely individuals coincides with interlocked structural alterations in the allocortical-neocortical covariations.

However, the relationship between the structure of the HC and DN and the tendency to feel socially disconnected from others may only be partially determined by one’s environment and life experience. Indeed, our genome-wide analyses showed that participant-specific expressions of HC-DN covariation have distinct links to the heritable components of loneliness. This is in accord with previous research showing that loneliness has consistent, yet subtle genetic underpinnings (132). For example, twin studies have identified genetic contributions to individual differences in feelings of loneliness that are as high as ∼48% (133). Recent genome-wide association studies with large sample sizes [e.g., n = 452,302 (134)] have also reported the contribution of common genetic variants to loneliness to range from 4% to 27% (134–136). These earlier genome-wide analyses have additionally pointed out a small number of specific gene loci that are significantly associated with loneliness and regular participation in social groups (134).

An underlying genetic contribution to the experience of loneliness is further supported by recent genome-wide correlation analyses. For example, one study has demonstrated that our UK Biobank loneliness trait significantly shared underlying genetic factors with each of 264 different demographic, lifestyle, and disease phenotypes (29). Another study has demonstrated a strong genetic correlation between loneliness and neuroticism (136). Taken together, these recent population neuroscience studies indicate that the genetic determinants underlying loneliness are probably quite polygenic and involve complex gene-environment interactions. In our study, these previous insights were extended by relating PRS for loneliness with the interindividual expression of multiple spatially overlapping signatures of HC-DN covariation. Importantly, we have identified a select subset of HC-DN covariation modes that had a significant relationship with genetic predisposition for loneliness.

In particular, we found a significant heritability effect for the expression of the loneliness-specific signature in modes 8 and 9, both of which also showed multiple subregion divergences in lonely individuals. We also noted a particular concordance between our PRS analysis of loneliness with our uncovered subregion-specific divergences in mode 9. Intriguingly, in our multiple regression analyses of the NAc, each subregion that showed divergence in mode 9 for loneliness (i.e., OFC, vmPFC, pSTS, and left parasubiculum) was found to have a strong relationship with the NAc. These convergent results may suggest that a specific aspect of the genetic predisposition for loneliness consists in innate tendencies involving reward-related processing. This view would be in line with past findings associating loneliness with social reward valence (7), substance abuse or dependence (12, 137), executive function (138, 139), and impulsive behaviors (140). Recently, it has even been reported that social isolation is linked to altered neural activity responses in midbrain regions and induces social craving in a way that may be similar to how fasting causes hunger (141). These reports further motivate a suspicion that the unmet desire for social connection in people is related to the same reward mechanisms in nonsocial reward. More specifically, the β-endorphin system, which is involved in pain management and reward valence, has also been shown to be integral to social bonding through touch in multiple mammalian species (142–144).

However, our findings on how structural HC-DN covariation relates to the genetic predisposition for social isolation may be better framed in terms of susceptibility to the environment, such as the previously proposed notion that the driving forces behind genetic heritability are not direct but rather come into effect through an altered sensitivity to environmental conditions (145). For example, whereas oxytocin and social support have been shown to interact in the stress response (146), a single-nucleotide polymorphism in the oxytocin receptor gene (OXTR) has been found to differentially affect stress response depending on adequate social interactions with others (147). Indeed, the heightened levels of stress in humans with social isolation has been argued to be fundamentally different from a simple general and diffuse stress response (6, 148).

In fact, animal studies have afforded detailed accounts of the effects of social isolation stress on the hippocampus (149), focusing on the apical dendrites of CA neurons (150). When considered in conjunction with this animal literature, our analysis could therefore extend our incomplete understanding of the brain basis of social isolation-related stress by identifying the precise allocortical and neocortical anatomical subregions at play. A possible manifestation of the effects of social isolation-related stress in our study is the pronounced total number of divergences in the molecular layer. In particular, we found seven total hits in the molecular layer across all examined modes for lonely individuals. As mentioned above, this molecular layer subregion consisted of a combination of the apical dendrites of subiculum and CA neurons, with CA1 composing the largest portion (40). In light of this, the observation of preferential hits in the molecular layer of lonely individuals broadly aligns with past studies on chronic stress in rodents, tree shrews, and nonhuman primates. Indeed, experimental studies in these species have reported detailed cellular and physiological consequences of induced stress on the hippocampus and other related brain regions.

In the hippocampus specifically, it has been shown that chronic stress results in the selective atrophy of the apical dendrites of CA3 cells but not CA3 basal dendrites (151–155). In contrast, elevated chronic stress in the face of multiple stressors results in additional atrophy of the apical dendrites of CA1 cells yet not CA1 basal dendrites (156, 157). On the basis of these findings, researchers have proposed that CA3 apical dendrites are the most susceptible hippocampal structure to chronic stress. Yet, stress affects CA1 apical dendrites only with a more severe experimental stressor than that needed to affect CA3 apical dendrites (150). Our findings of preferential hits in the molecular layer of lonely individuals may thus be a marker of elevated stress and relate to the physiological characterization of loneliness as a trait coinciding with elevated stress and immune response levels (5, 8, 158–161).

Building upon this, our concerted analysis framework for HC and DN subregion variation thus allowed us to relate the potential effects of chronic stress on neocortical partners in the DN. For example, in mode 1 we observed four hits in the molecular layer with coincident hits in midline DN subregions such as the dmPFC, PCC, and RSC. The collection of parallel findings of divergences between the stress-susceptible molecular layer in the allocortex and specific midline cortical subregions might therefore reflect the effects of stress on the higher associative cortex.

Indeed, previous studies have shown that the effects of chronic stress on the hippocampus percolate into interaction partners with known axonal connections, such as target areas in the mPFC (28, 81, 82, 162, 163). For example, in one study where rats were submitted to 4 wk of chronic stress, they were found to have a deficit in hippocampal-PFC synaptic strength compared with control animals (162). It was also found that the stressed rats had lower PFC volume and poorer performance on working memory and behavioral flexibility tasks—activities that are thought to depend on mPFC activity (162). Another more recent study has demonstrated that among squirrel monkeys with bouts of prepubertal social isolation, there was an association between functional hyperconnectivity in PFC-subcortical circuits (e.g., dmPFC, vmPFC, OFC; amygdala, ventral striatum, HC) with reduced anxiety-like behavior in later adulthood. In line with these previous experiments, our analytical strategy in the UK Biobank population allowed us to identify concomitant stress-related divergences in the molecular layer of the hippocampus and divergences in the mPFC, PCC, and RSC of the human DN.

It is well known that chronic stress causes measurable neural consequences for the hippocampus in different animal species (164–166). Yet, the impact of social isolation and chronic stress on the hippocampus has repeatedly been argued to be reversible in a matter of a few weeks with social rehabilitation or reintroduction to an enriched environment (33, 34, 167, 168). In a similar vein, in the higher association cortex the effects of chronic stress on PFC dendritic morphology in rats and functional connectivity in humans have been shown to be reversible after stressful experiences (169–172). In light of this past literature, there is likely to be elasticity to the effects of social isolation and perhaps the stress-related volumetric divergences that we uncovered in the hippocampus and its neocortical interaction partners. Indeed, social rehabilitation by returning to the usual social environment may remold the examined allocortical-neocortical brain circuits.

SUPPLEMENTAL DATA

Supplemental Figs. S1–S9 and Supplemental Table S1: https://doi.org/10.6084/m9.figshare.15060684.

GRANTS

This project has been made possible by the Brain Canada Foundation, through the Canada Brain Research Fund, as well as by NIH Grant R01AG068563A and the Canadian Institutes of Health Research. D.B. was also supported by the Healthy Brains Healthy Lives initiative (Canada First Research Excellence fund) and by the CIFAR Artificial Intelligence Chairs program (Canada Institute for Advanced Research) as well as a Research Award and Teaching Award by Google.

DISCLOSURES

No conflicts of interest, financial or otherwise, are declared by the authors.

AUTHOR CONTRIBUTIONS

C.Z., R.N.S., and D.B. conceived and designed research; D.B. performed experiments; C.Z. and D.B. analyzed data; C.Z. and D.B. interpreted results of experiments; C.Z. and D.B. prepared figures; C.Z. and D.B. drafted manuscript; C.Z., R.N.S., and D.B. edited and revised manuscript; C.Z., R.N.S., and D.B. approved final version of manuscript.