Psychological resilience correlates with EEG source-space brain network flexibility

We aimed at identifying the potential relationship between the dynamical properties of the human functional network at rest and one of the most prominent traits of personality, namely resilience. To tackle this issue, we used resting-state EEG data recorded from 45 healthy subjects. Resilience was quantified using the 10-item Connor-Davidson Resilience Scale (CD-RISC). By using a sliding windows approach, brain networks in each EEG frequency band (delta, theta, alpha, and beta) were constructed using the EEG source-space connectivity method. Brain networks dynamics were evaluated using the network flexibility, linked with the tendency of a given node to change its modular affiliation over time. The results revealed a negative correlation between the psychological resilience and the brain network flexibility for a limited number of brain regions within the delta, alpha, and beta bands. This study provides evidence that network flexibility, a metric of dynamic functional networks, is strongly correlated with psychological resilience as assessed from personality testing. Beyond this proof-of-principle that reliable EEG-based quantities representative of personality traits can be identified, this motivates further investigation regarding the full spectrum of personality aspects and their relationship with functional networks.


INTRODUCTION
An evolving field of neuroscience aims to reveal the neural substrates of personality, referring to the relatively stable character of an individual that influences her or his long-term behavioral style (Dubois et al., 2018). A key personality character is resilience, defined as the ability to adapt to stress, adversity, and negative events and cope actively with life challenges (Fletcher & Sarkar, 2013;Luthar, 2003;Rutter, 2006). Recently, a multisystem model of resilience has been proposed by Liu et al. (2017), in which resilience comprises three structures: (1) the innermost layer, which comprises the physiological, biological, and demographic profiles of an individual; (2) the intermediate layer, which includes internal factors such as family, friends, and personal experiences; and (3) the outermost layer, which corresponds to external resilience such as access to healthcare and social services. However, except a few efforts (Kong et al., The last decade has witnessed an increase of studies that consider the human brain as a large-scale network. It is thus unsurprising that network neuroscience, which uses tools from graph theory to better understand neural systems (Bassett & Sporns, 2017), has become one of the most promising approaches to link behavior to brain function, including personality traits (Markett et al., 2018). In the network neuroscience model, the human brain is summarized by a set of nodes representing brain regions and a set of edges representing the connections between these brain regions. The magneto/electro-encephalography (M/EEG) source-space EEG: Electroencephalography is a neuroimaging technique consisting of recording with scalp electrodes the small currents induced by neuronal activity at the cortex level.
networks provide a unique direct and noninvasive access to those electrophysiological brain Source-space networks: Networks formed by interacting anatomical brain regions, as opposed to sensor-level networks.
networks, at the milliseconds temporal scale Mheich et al., 2018). The excellent time precision of this method allows the tracking of the brain networks dynamics at an unprecedented timescale. In this paper, we aimed to test the hypothesis that metrics derived from brain networks dynamics can be correlated to resilience scores.
To test this hypothesis, we used resting-state EEG data recorded from N = 45 healthy subjects. Resilience was quantified using the Connor-Davidson Resilience Scale (CD-RISC) (Campbell-Sills & Stein, 2007), where higher scores correspond to greater resilience. Brain networks in each frequency band (delta, theta, alpha, and beta) were constructed using the EEG source-space connectivity method. By applying multislice modularity algorithms on the dynamic networks, the reconfiguration of EEG source-space networks was quantified using flexibility, defined as how often a given node changes its modular affiliation over time, com-Flexibility: Metric from graph theory quantifying how often one node of a network changes its affiliation with a module over time.
puted at the level of each brain region. The choice of flexibility as an appropriate metric was made following a recent study discussing the relationship between functional networks and personality traits (Tompson et al., 2018). Indeed, flexibility captures, in an intuitive manner, the dynamics of functional networks. Results revealed essentially a negative correlation between psychological resilience and EEG-based brain functional network flexibility.

RESULTS
The study is summarized in Figure 1. First, EEG data were recorded and preprocessed. Second, the dynamic networks were estimated using the EEG source-space approach, giving a set of brain network at the given time period for each frequency band. Third, the flexibility of each brain region in each frequency band were computed for each subject. Finally, the correlation between the brain regions flexibility and the resilience score (CD-RISC) was computed.
Then, we computed the Pearson correlation between the resilience score and the global flexibility (averaged over all brain regions for each subject in each frequency band). The corresponding results are presented in Figure 2. For all frequency bands, a negative correlation was observed. This negative correlation was significant for delta (R = −0.51, p = 0.0003), alpha (R = −0.41, p = 0.004), and beta (R = −0.43, p = 0.002) bands, but nonsignificant for the theta band (R = −0.19, p = 0.2).
We then focused on the correlation between flexibility and the resilience score at the level of each brain region. In the delta band, brain regions that were significant and correlated with the CD-SCORE are illustrated in      (e.g., the brain regions that had the highest 90% correlation between the resilience score and the network-based score).
The results for the alpha band are presented in

DISCUSSION
Our results reveal a robust and direct relationship between the network flexibility (as quantified through graph theory) of human functional brain networks and resilience (as measured with neuropsychological testing). This result has been obtained using noninvasive brain recordings (EEG), in the absence of any specific task or stimulus (resting-state recordings), for a relatively large sample size (N = 45). Specifically, the network flexibility of specific brain regions is shown to be significantly decreased with a higher resilience score, which can appear nonintuitive: if resilience refers to a capacity to cope with external stress and adapt, it could imply that brain networks can reconfigure to adapt to these changes. Our results point at the opposite: most brain regions negatively correlated with the resilience score are part of the "rich club" ("core" functional network of brain regions [Van Den Heuvel & Sporns, 2011]) and illus-Rich club: Subset of cortical regions that are extensively connected with other brain regions.
trate that the brain core network is less flexible in resilient subjects. This suggests that a stable core network would allow preventing excessive reconfiguration of functional brain networks following external stress factors, for example. From a fundamental perspective, these results also suggest that there exists a brain network of critical importance, which underlies one of the most prominent personality traits (resilience). This raises the possibility to probe further potential brain network correlates of other personality traits in the future. Furthermore, this link between core network structures and resilience is fundamental, since it can be detected noninvasively even in resting-state recordings. Another insight from our analysis is that dynamic but not static network analysis was able to reveal this association between resilience and a Static network analysis: Analysis of the brain network involving some form of averaging over time.
graph theory metric. Static analyses were realized using clustering coefficient and participation coefficient, and did not reveal any significant correlation with the resilience score (see Supporting Information Figure S1 and Figure S2; Paban, Modolo, Mheich, & Hassan, 2019). Therefore, this highlights the importance of the dynamic reconfiguration of brain networks in shaping behavior. Recent studies have emphasized the interest to characterize the architecture underlying these dynamics, for example, by identifying the "life duration" and transition rates of these frequency-specific networks (Vidaurre et al., 2018). Overall, as emphasized in a recent perspective paper by Tompson et al. (2018), there is a significant potential regarding the use of network neuroscience tools to unveil the brain networks underlying personality traits, with flexibility and integration being especially promising.

Neural Substrates of Resilience
Functional neuroimaging studies examining brain networks associated with resilience are limited. Most of them mainly focused on patient populations, such as depressive, traumatized, and posttraumatic stress disorders, in which alterations where described in regions of the brain involved in emotion and stress regulation circuitry (see van der Werff et al., 2013, for review). In healthy adults, to our knowledge, only one resting-state fMRI (rs-fMRI) experiment has been reported to date. Interestingly, Kong et al. (2015) reported that psychological resilience had significant negative correlations with the rs-fMRI signals, but in regions such as the bilateral insula, right dorsal, and rostral anterior cingulate cortex. Our study did not highlight the same brain regions. Our dynamic network analysis identified regions belonging to the "core" func-Dynamic network analysis: Analysis of the brain network made using time-dependent estimates of the network, as opposed to static network analysis. tional network described by Van Den Heuvel and Sporns (2011). Those regions are supposed to play a central role in establishing and maintaining efficient global brain communication. They have been shown to be involved in cognitive processes related to top-down attentional control (superior parietal cortex [Sestieri et al., 2017]), decision-making (subdivision of the orbitofrontal cortex [Schuck et al., 2018]), cognitive regulation of behavior (caudal anterior cingulate [Bush et al., 2000]), and spatial location (entorhinal cortex [Kim & Maguire, 2018]). Our data also highlighted the regions participating in the visual network (e.g., cuneus, the lingual gyrus, and the pericalcarine cortex). This network oscillated in fast frequency band (beta), which is thought to reflect different aspects of sensory information processing (Hong et al., 2008). Another region negatively correlated with the resilience score comprised the inferior parietal cortex. It has been shown to be implicated in major cognitive functions, including visuo-spatial attention and visual memory (Egner et al., 2008).
The present study provided evidence of neuronal substrates of resilience. In healthy adults, it seems to make sense that the ability of individuals to cope actively with life's challenges (Fletcher & Sarkar, 2013) requests brain regions involved in high-level cognitive processes. To be efficient, these regions should have a modular organization stable over time. As suggested by Betzel et al. (2016), such network architecture limits information exchange across modules, which allows a specialized information processing. Such modular organization seemed to be an important factor in maintaining a stable equilibrium of psychological and physical functioning when facing adversity, ranging from daily problems to major life events (Bonanno, 2004). When this is not the case, that is, when few brain regions showed high modular affiliation exchanges, one may assume that they multiply their contributions and so participate in a multitude of brain functioning. In this context, subjects had low-resilient scores, meaning that they exhibit passive coping strategies, high emotional load, ruminative and depressogenic thinking, and low life satisfaction (Davydov et al., 2010).

Methodological Considerations
First, in the present study we used a template source-space, instead of a subject-specific one, that is, the same structural MRIs of healthy subjects were used in our EEG functional connectivity analysis. In the case of healthy subjects, it was reported that coregistration with the template yielded largely consistent connectivity and network estimates compared with native MRI (Douw et al., 2018). However, in the case of patient-related studies, it is more preferred to use patient-specific MRIs.
Second, it is important to keep in mind that computing the functional connectivity at the scalp level is severely corrupted by the volume conduction problem (Brookes et al., 2014;Lai et al., 2018;Schoffelen & Gross, 2009). M/EEG connectivity analysis at the source-space was shown to reduce this effect as connectivity is applied to "local" time series generated by cortical neuronal assemblies modeled, for instance, as current dipole sources. Yet, the "mixing effects" can also occur in the source-space. To address this issue, a number of methods were developed based mostly on the rejection of zero-lag correlation. In particular, "unmixing" methods, called "leakage correction" (such as the orthogonalization approach [Colclough et al., 2015]), have been reported, which force the reconstructed signals to have zero crosscorrelation at lag zero. Here we preferred to use the phase locking value (PLV) metric and keep these zero-lag correlations, as several experimental studies reported the importance (in many conditions) of the zero-lag synchronizations in the human brain. Nevertheless, we believe that there is no ideal solution yet for this methdological issue and that further efforts are needed to completely solve the spatial leakage problem.
Third, we have performed the whole analysis in all EEG frequency bands in order to avoid making any a priori judgments about the involved rhythms. The reported effects were the same when other thresholds were used. However, it should be kept in mind that standard EEG frequency bands are increasingly associated with specific functional roles in brain-scale information processing. As an exemple, alpha activity appears to mediate the so-called "gating by inhibition" mechanism involved in information routing (Bonnefond et al., 2017) through specific nodes of the cortical network while inhibiting irrelevant regions. On the other hand, theta oscillations appear to be involved in locking distant brain regions to enable further processing (Lisman and Jensen, 2013). Since different frequency bands appear linked with specific functional roles, there is no reason a priori that psychological resilience results from an interplay between all of those frequency bands. Therefore, we believe that our results point at those oscillations that are functionally linked with resilience, and that theta oscillations do not appear to have a specific role in this case.

Psychological resilience and functional networks flexibility
Finally, although EEG source connectivity has improved the spatial resolution of EEG, the networks identified are, however, still limited to cortical gray matter. In fact, the localization of subcortical structures remains extremely challenging using EEG technology, namely because of anatomical and analytical reasons. For example, unlike the layered cortex, a subcortical region does not have a sufficient organization of pyramidal cells to give rise to localizable scalp-recorded EEG.

Participants
A total of 45 healthy subjects were recruited (22 women). The mean age was 34.7 years old (SD = 9.1 years, range = 18 − 55). Education ranged from 10 years of schooling to a PhD degree. None of the volunteers reported taking any medication or drugs, nor suffering from any past or present neurological or psychiatric disease. After receiving approval from the Aix-Marseille university ethics committee according to the Declaration of Helsinki, participants filled out the CD-RISC questionnaire at home approximately 1 week before the EEG experiment. Written informed consent was obtained from all participants prior to study onset.

Measuring Psychological Resilience
The Connor-Davidson Resilience Scale (CD-RISC; Connor & Davidson, 2003) is a 25-item scale that measures the ability to cope with adversity. A 10-items version (CD-RISC 10) of this scale has been developed by Campbell-Sills and Stein (Campbell-Sills & Stein, 2007). A 10-item version validated for French speaking populations was used in the present study (Guihard et al., 2018;Scali et al., 2012). The 10 items are rated on a five-point Likert scale that ranges from 0 (not at all) to 4 (true nearly all of the time). Higher scores correspond to greater resilience. This scale demonstrated good internal consistency and construct validity (Campbell-Sills & Stein, 2007). In our sample, the CD-RISC 10 exhibited a reliability of α = 0.90.

Data Acquisition and Preprocessing
Each EEG session consisted of a 10-min resting period with the participant's eyes closed (Paban et al., 2018). Participants were seated in a dimly lit room, were instructed to close their eyes, and then to simply relax until they were informed that they could open their eyes. Participants were instructed that the resting period would last approximately 10 min. The eyes-closed resting EEG recordings protocol was chosen to minimize movement and sensory input effects on electrical brain activity. EEG data were collected using a 64-channel Biosemi ActiveTwo system (Biosemi Instruments, Amsterdam, The Netherlands) positioned according to the standard 10-20 system montage, one electrocardiogram, and two bilateral electro-oculogram electrodes (EOG) for horizontal movements. Nasion-inion and preauricular anatomical measurements were made to locate each individual's vertex site. Electrode impedances were kept below 20 kOhm. EEG signals are frequently contaminated by several sources of artifacts, which were addressed using the same preprocessing steps as described in several previous studies dealing with EEG resting-state data Rizkallah et al., 2018). Briefly, bad channels (signals that are either completely flat or contaminated by movement artifacts) were first identified by visual inspection, complemented by the power spectral density. These bad channels were then recovered using an interpolation procedure implemented in Brainstorm (Tadel et al., 2011) by using neighboring electrodes within a 5-cm radius. Epochs with voltage fluctuations more than +80 μV and less than −80 μV were removed. Five artifact-free epochs of 40-s lengths were selected for each participant. This epoch length was used in a previous study, and was considered as a good compromise between the needed temporal resolution and the results reproducibility . By using a sliding windows approach to compute the functional connectivity, a large number of networks (depend on the analyzed frequency band) were obtained for each 40-s epoch.

Brain Networks Construction
First, brain networks were reconstructed using the "EEG source-space connectivity" method (Hassan et al., 2014;, which includes two main steps: (1) reconstruct the dynamics of the cortical sources by solving the inverse problem, and (2) measure the statistical couplings (functional connectivity) between the reconstructed time series. In summary, EEG source connectivity links the recorded EEG signals with the functional relationship between anatomical brain regions (e.g., networks), through the EEG inverse problem that provides the localization of the cortical sources originating these EEG signals. Let us mention that numerous methods exist for both of these two steps (EEG inverse problem and functional connectivity measures). EEGs and MRI template (ICBM152) were coregistered through the identification of anatomical landmarks by using Brainstorm (Tadel et al., 2011). A Desikan-Killiany atlas-based segmentation approach was used, consisting of 68 cortical regions (Desikan et al., 2006). The OpenMEEG (Gramfort et al., 2010) software was used to compute the head model. Here, we used the weighted minimum norm estimate (wMNE) algorithm as an inverse solution. The reconstructed regional time series were filtered in different frequency bands (delta,(1)(2)(3)(4)theta,(4)(5)(6)(7)(8)alpha,(8)(9)(10)(11)(12)(13)and beta,(13)(14)(15)(16)(17)(18)(19)(20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30). For each frequency band, functional connectivity was computed between the regional time series using the PLV measure (Lachaux et al., 1999). This combination wMNE/PLV was chosen according to a recent model-based comparative study of different inverse/connectivity combinations  Using PLV, dynamic functional connectivity matrices were computed for each epoch using a sliding windows technique. It consists of moving a time window of certain duration δ, and PLV is calculated within each window. As recommended in Lachaux et al. (2000), we selected the smallest window length equal to 6 central frequency , where 6 is the number of "cycles" at the given frequency band. For instance, in the theta band, since the central frequency (Cf) equals to 6 Hz, δ equals 1 s, and δ = 279 ms (Cf = 21.5 Hz) in beta band. Thus, for each epoch, 33 networks were obtained for the theta band, and 130 networks in the beta band. The same calculation was adopted for other frequency bands. Finally, we kept only the strongest 10% of connections.

Network Modularity and Flexibility
The obtained dynamic matrices were divided into time-dependent modules by using the multislice community detection approach described in Mucha et al. (2010). It consists of introducing a parameter that associates nodes across time, before applying the modularity procedure, and is defined as (one modularity value is computed for each frequency band of interest; : where nodes i and j are assigned to modules σ is and σ js in window s, respectively. A ijs represents the weight of the edge between these two nodes in window s; γ is the structural resolution parameter; ω is called the "inter-layer coupling strength," and links the same node j within the network evaluated at two different times ("layers"): therefore, ω jrs links the node j at times r and s, that is, in layer r and in layer s. Here, we chose γ = ω jrs = 1; p ijs represents the expected number of links according to a null model. As a reminder, the Kronecker δ-function is such that δ(x, y) is 1 if x = y, and 0 otherwise.
The multilayer network modularity was computed 100 times since Q may vary from run to run (degeneracy). This step is mandatory precisely because the degeneracy intrinsic to the community detection algorithm results in different networks structures depending on the run. Classically, this issue is dealt with using the so-called "co-classification matrix" or "consensus matrix," whose elements indicate the ratio of each region to be in the same module with the other regions among these 100 partitions. Only elements in the consensus matrix higher than an appropriate random null model were considered. The randomized networks (generated with the null model) have the same strength distribution as the true brain networks. To quantify the dynamics of brain networks, we used the flexibility metric proposed in Bassett et al. (2011). The flexibility of a brain region is defined as the number of times that a brain region changed modular assignment throughout the session, normalized by the total number of changes that were possible.

Statistical Analysis
We assessed the correlation between the two conditions (resilience score vs. network were assessed using the Pearson correlation). We applied a FDR correction for multiple comparisons across regions (Benjamini & Hochberg, 1995).

Software
The functional connectivity, network measures, and network visualization were performed using BCT (Rubinov & Sporns, 2010), EEGNET (Hassan et al., 2015), and BrainNet viewer (Xia et al., 2013), respectively. The Network Community Toolbox (http://commdetect.weebly.com/) was used to compute the consensus matrices as well as the values provided by the flexibility metrics.