Frequency-Specific Resting Connectome in Bipolar Disorder: An MEG Study

Background Bipolar disorder (BD) is a serious psychiatric disorder that is associated with a high suicide rate, and for which no clinical biomarker has yet been identified. To address this issue, we investigated the use of magnetoencephalography (MEG) as a new prospective tool. MEG has been used to evaluate frequency-specific connectivity between brain regions; however, no previous study has investigated the frequency-specific resting-state connectome in patients with BD. This resting-state MEG study explored the oscillatory representations of clinical symptoms of BD via graph analysis. Methods In this prospective case-control study, 17 patients with BD and 22 healthy controls (HCs) underwent resting-state MEG and evaluations for depressive and manic symptoms. After estimating the source current distribution, orthogonalized envelope correlations between multiple brain regions were evaluated for each frequency band. We separated regions-of-interest into seven left and right network modules, including the frontoparietal network (FPN), limbic network (LM), salience network (SAL), and default mode network (DMN), to compare the intra- and inter-community edges between the two groups. Results In the BD group, we found significantly increased inter-community edges of the right LM–right DMN at the gamma band, and decreased inter-community edges of the right SAL–right FPN at the delta band and the left SAL–right SAL at the theta band. Intra-community edges in the left LM at the high beta band were significantly higher in the BD group than in the HC group. The number of connections in the left LM at the high beta band showed positive correlations with the subjective and objective depressive symptoms in the BD group. Conclusion We introduced graph theory into resting-state MEG studies to investigate the functional connectivity in patients with BD. To the best of our knowledge, this is a novel approach that may be beneficial in the diagnosis of BD. This study describes the spontaneous oscillatory brain networks that compensate for the time-domain issues associated with functional magnetic resonance imaging. These findings suggest that the connectivity of the LM at the beta band may be a good objective biological biomarker of the depressive symptoms associated with BD.


INTRODUCTION
Bipolar disorder (BD) affects up to 60 million individuals globally (1). BD can lead to cognitive and functional impairments (2) and is associated with a high mortality rate, particularly as a result of suicide (3). Although suicide attempts can be triggered by stressful life events such as bereavement (4) and biological abnormalities (5), the presence of a psychiatric disorder is the most consistent risk factor for suicidal behavior (6). Therefore, it is important to diagnose BD in the early clinical stage. However, it is difficult to diagnose BD quickly and accurately in clinical practice because its onset is most commonly characterized by a depressive episode, which often resembles unipolar depression (7). The efficacy of pharmacological interventions in the treatment of BD depends on an accurate and early diagnosis (8); however, a diagnostic biomarker has yet to be identified (9). Accordingly, for the accurate and timely diagnosis of BD, there is an urgent need to establish objective clinical biomarkers that are based on a better understanding of BD neuropathology. Studies investigating the resting-state functional connectivity have revealed intrinsic properties of brain disorders (10)(11)(12)(13) and are expected to help identify such biomarkers.
Until recently, most studies investigating the functional connectivity in psychiatric patients have used functional magnetic resonance imaging (fMRI) (14). Although fMRI has a high spatial resolution, it measures an indirect neural activity and can only investigate slow oscillations in brain activity (< 0.1 Hz) (15). Compared to fMRI, a major advantage of magnetoencephalography (MEG) is its high temporal and spatial resolution, which enables synchronizations of the local and long-range neural activities to be assessed (16). In addition, the source reconstruction approach using MEG provides valuable spatial information that characterizes neural network modulations (14). Although MEG has many advantages, to our knowledge, there have been no studies of resting-state MEG for evaluating functional connectivity in patients with BD.
Unfortunately, a lack of established methodologies has prevented its widespread use for brain resting-state connectivity studies. In recent years, however, two British groups have developed methodological approaches to evaluate the frequencyspecific resting-state connectivity using MEG (17,18). These advances have enabled the investigation of frequency-specific alterations of resting-state networks via graph theory in autism (19) and age-related changes in these networks (20).
Simultaneously, graph theory analysis has been applied to explore brain network connectivity in psychiatric patients (21)(22)(23)(24)(25). Graphs are defined as sets of nodes and edges that represent system components and their mutual relationships (26). Using graph theory, brain networks can be represented as graphs that consist of nodes (which correspond to brain regions) and edges (which represent their connections) (27,28). Therefore, graph theory is ideally suited for investigating complex networks such as human brain networks (26) and may be a potential tool for investigating the network characteristics of psychiatric disorders that alter connectivity within brain networks (29).
The main objective of this study was to investigate the restingstate frequency-specific functional connectivity to glean deeper insights into the pathophysiology of BD and help identify objective clinical biomarkers for the disorder. Here, we describe a novel MEG analysis that explores the resting-state functional connectivity using graph analysis in patients with BD. We hypothesized that resting-state MEG studies would provide important information about the spontaneous oscillatory brain networks, thereby representing an advantage over the sole use of fMRI-based resting-state investigations. We believed that this approach could also help identify clinical objective biomarkers for psychiatric disorders. To achieve this goal, we explored oscillatory representations of subjective and objective depressive symptoms in patients with BD using graph analysis.

Participants
Seventeen patients with BD and 22 healthy controls (HCs) participated in this prospective case-control study. All patients with BD were recruited from the Gunma University Hospital in Gunma Prefecture, Japan ( Table 1). All enrolled participants were right-handed, as assessed using the Edinburgh Handedness Inventory (30). To diagnose BD and to exclude participants with a history of psychiatric disorders from the HCs, we administered the Structured Clinical Interview of the Diagnostic and Statistical Manual of Mental Disorders, Fourth Edition for Axis I (31) and Axis II Disorders (32). Subjective and objective depressive and manic symptoms were assessed using the Beck Depression Inventory-Second Edition (BDI-II), Altman Self-Rating Mania (ASRM) scale, 17-item Hamilton Depression Rating (HAMD-17) scale, and Young Mania Rating Scale (YMRS). Patients older than 60 years were excluded to rule out the possible interference of additional pathophysiological factors, such as aging and cerebrovascular changes. Most of the patients with BD (14/17) were using medications, including antipsychotics, mood stabilizers, antidepressants, anxiolytics, hypnotics, and/or antiparkinsonian drugs. The equivalent doses of antidepressants (imipramine equivalent dose in mg/day), antipsychotics (chlorpromazine equivalent dose in mg/day), anxiolytics (diazepam equivalent dose in mg/day), and hypnotics (flunitrazepam equivalent dose in mg/ day) were calculated for each patient (33). HCs had no history of major psychiatric or physical illnesses and medication use. The exclusion criteria for both groups included a clear abnormality in brain magnetic resonance imaging (MRI) findings, neurological illness, traumatic brain injury with a known cognitive consequence or loss of consciousness for more than 5 min, substance use or addiction, and the presence of hearing or vision impairment. The study protocol was approved by the Ethics Committee of the Gunma University Graduate School of Medicine, Japan. All participants provided written informed consent before the study. This study was carried out in accordance with the latest version of the Declaration of Helsinki.

MEG Data Acquisition and Preprocessing
For each participant, we acquired the MEG data for 7 min. A crosshair was presented as a fixation point in front of the participants. The only instruction provided to participants was to keep their eyes open, to relax, and remain awake. To confirm the awareness level and evaluate the extent of sleepiness during MEG acquisition, we used the Stanford Sleepiness Scale (34).
During the experiment, participants sat in an upright position. The neuromagnetic signals were acquired inside a magnetically shielded room (JFE Mechanical Co., Tokyo, Japan) using a wholehead Elekta Neuromag Vector View system (Oy, Helsinki, Finland), which was composed of 306 sensors arranged as 102 triplets (each unit consisting of two orthogonal planar gradiometers and one magnetometer). First, the locations of three fiduciary points (the nasion, and left and right auricular points) defining the head frame coordinate system and a set of head surface points (approximately 150-200 points) (35), as well as those of four head position indicator coils, were digitized using an Isotrak three-dimensional digitizer (Polhemus ™ , Colchester, VT, USA). During acquisition, the sampling rate was 1000 Hz, and the signals were bandpass filtered between 0.1 and 200 Hz. All off-line analyses were conducted using continuous raw data.
The current source estimate was carried out as follows. For noise suppression, we applied the "oversampled temporal projection" method to improve the sensor noise levels (36). Subsequently, the data were spatially filtered using the signal space separation method (37,38) with Maxfilter (Elekta, Stockholm, Sweden) to suppress the noise generated by sources outside the brain. A notch filter was used to eliminate noise from the power line (50 Hz) and its harmonics (100 and 150 Hz). We eliminated heartbeat-and eye movement-related artifacts from the raw data using independent component analysis (ICA). The ICA was applied to the MEG sensor signals, and the eye movement signals were isolated based on the hybrid approach of component discrimination using the time (signal waveform) and space (topology) (39).

MEG Source Imaging
For the MEG source localization analysis, we used the minimum norm estimates (MNE) suite and MNE python (43). We employed the boundary-element method, with a single compartment bound by the inner surface of the skull for the forward head modeling (44). The source current distribution was estimated using the noise-normalized dynamic statistical parameter mapping (dSPM) method (45,46). The noise covariance matrix was computed from the empty room acquisition data. To reduce the bias of the estimated source locations towards the superficial currents, we incorporated the depth-weighting parameter (depth = 0.8) by adjusting the source covariance matrix (47).

Data Analysis Procedure for Resting-State Networks
Schematic images of the data analysis procedure used in this study are illustrated in Figure 1. We started with MEG sensor detection level analysis and then evaluated the source distribution as described in subsection MEG Data Acquisition and Preprocessing ( Figure 1A). Subsequently, we imported the commonly employed cortex segmentation (Destrieux Atlas) from the FreeSurfer software, which automatically divided the cortex into 150 regions (48). After discarding the "Unknown" regions, the regions were further subdivided into a total of 446 cortical labels (20) (Figure 1B), such that each label covered almost identical areas through all the individual and standard brains. We employed this high-resolution parcellation because the cortical surface was convoluted, and the estimated source data were averaged across a large label, which crosses multiple sulci and gyri and possibly results in signal cancellation across the parcel (20). We first applied the above process to the standard brain (Montreal National Institute 305; provided by FreeSurfer software) (49), and then projected this onto each of the individual brains. Source leakage correction was essential for the evaluation of long-range correlations with MEG (50) and was determined in this study as described below. As shown in Figure 1C, we applied the orthogonal correlations to compute the connectivity metrics involved in handling source leakages using the "envelope_correlation" function in MNE (17,20) at delta (1-4 Hz), theta (4-8 Hz), alpha (8)(9)(10)(11)(12), low beta (13-20 Hz), high beta (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30), and gamma (30-80 Hz) bands (17). The "envelope_correlation" obtains the correlation values according to the following calculation flow. We started with x(t) and y(t) as the extracted source waveforms from specific ROI, where t denotes time stamps. First, the extracted source waveforms x(t) and y(t) were bandpass filtered at each frequency band and acquired x (t,f) and y(t,f). For each step-frequency band, the complex signal, x(t,f) and y(t,f), was generated with the Hilbert transform from x (t,f) and y(t,f). The connectivity between the two complex signals x(t,f) and y(t,f) was calculated by orthogonalizing one signal with respect to the other Y ⊥X (t,f).
The superscript * denotes the conjecture of the complex value, and "Im" denotes the imaginary part of the complex value. The power envelopes Env X (t,f) were subsequently determined by computation of the absolute value of X(t,f), as follows: Similarly, Env Y ⊥ X (t, f ) was obtained from Y ⊥X (t,f). We calculated the linear correlation between Env Y ⊥ X (t, f ) and Env Y ⊥ X (t, f ). Given the slow time course of these envelopes, and to ensure that enough independent samples were available in the correlation window (17), we calculated the orthogonal correlation using an overlapping sliding window of 30 s with a stride of 1/6 of the window size. Finally, we calculated the median values of sliding correlation values and applied the graph-theoretical analysis to these median correlation values ( Figure 1C) as described below.

Graph Analysis
We used the graph-theoretical measures (intra-community edges and inter-community edges) to quantify the basic property of resting-state networks (RSNs). We used the connectivity matrix that contained the orthogonal correlations between all N × N node pairs for each frequency band and acquired a connectivity array with the dimensions of N Nodes × N Nodes × N Bands for each subject. Thresholding and binarizing the connectivity matrices resulted in adjacency matrix A ( Figure 1D). We used a threshold proportional scheme to retain a given proportion of the strongest connectivity matrix entries in the connectivity matrix. As there is no rationale for using a cost threshold, we compared the graph network properties for a wide range of costs and used a thresholding range of 5% to 30% at increments of 5%. There was no significant difference among thresholding values ( Figure S1), therefore, we set a threshold at the top 15% to ensure that the number of network connections (14,851 connections among all regions of interest) was the same across all individuals.
We investigated the changes of RSNs in BD from the existing network understanding framework provided by previous fMRI studies. According to the network template, we separated network modules into seven left and right network modules, including the frontoparietal network (FPN), limbic network (LM), salience network (SAL), dorsal attention network (DAN), default mode network (DMN), somatomotor network (SMN), and visual network (VN) (51) (Figures 1E and 2). In this way, we investigated differences in the connection patterns of intra-and inter-networks community between BD and HC groups. These network templates do not contain sub-cortical regions. We subsequently subdivided all edges into intra-community edges, which connect regions within each network, and inter-community edges, which connect inter-network regions, and measured the number of intra-community and inter-community edges for each frequency band ( Figure 1F).
We compared the numbers of inter-community and intracommunity edges for each frequency band between the HC and BD groups using an independent t-test in SciPy (52). Group differences were considered to be statistically significant at a pvalue of less than 5%. We examined the relationship between significantly different edges (p-value less than 1% significance level) in patients with BD and clinical measures (the HAMD-17, BDI-II, ASRMS, and YMRS) using Spearman's rank correlation analysis. The correlation values were determined at a significance level less than 5%. Table 1 summarizes the demographic characteristics of the sample. Age and sex ratios were not significantly different between the two groups (t = -1.4643, p = 0.152; chi-square = 0.0198, p = 0.888). At the time of the study, the patients with BD were euthymic or depressive as indicated by their HAMD-17 and YMRS scores. Figure 3 shows the results of the graph theory analysis, in which the number of inter-and intra-community edges were measured. Figure 3A depicts the summarized degrees of inter-and intracommunity edges in six frequency band ranges (delta, theta, alpha, low beta, high beta, and gamma), including the plot box bars at a significance level < 1% ( Figure 3B). We compared the number of edges between HC and BD using an independent ttest ( Figure 3A, left).

Inter-and Intra-Community Edges
In the delta band range, we found lower inter-community edges among the FPN and other network modules in the right To summarize, among the six frequency band ranges, the LM exhibited higher community edges irrespective of the A B FIGURE 3 | The results of the graph theory analysis for inter-and intra-community edges. (A) Graph theory was applied to six frequency band ranges. The left-most bar indicates the color scale of the t-value. All figures in (A) indicate t-values of results that compared inter-community edges (line thickness and color) and intracommunity edges (circle size and color) between healthy controls and patients with bipolar disorder (BD) at each frequency band. Thicker and deeper lines indicate higher inter-community edges, and larger circles indicate higher intra-community edges in patients with BD. If the p-value of the inter-or the intra-community edges was below the threshold for significance, the circles or lines are bordered by a thick line (p < 0.01) or thin line (p < 0.05). (B) The results of further statistical analysis for stage (A). The green boxplot shows the results from the HC group, and the purple boxplot shows the results from the BD group. Each vertical boxplot bar indicates the standard error. hemispheric differences, while the SAL and FPN exhibited a lower community edge among lower frequency band ranges. In particular, the number of inter-and intra-community edges were significantly different between the following groups at each frequency band (all, p < 0.01): the right FPN-right SAL at the delta, left SAL-right SAL at the theta (HC > BD), left LM at high beta, and right LM-right DMN at gamma (HC < BD) ( Figure 3B).

Relationship Between Characteristics of Graph Indices and Clinical Symptoms
We further evaluated the clinical interpretations and relationships via commonly used approaches in psychiatric diagnosis. A Spearman's rank correlation analysis was used to assess the relationship between significantly different graph indices (interand intra-community edges) and clinical measures (BDI-II, HAMD-17, ASRM, and YMRS scores). Figure 4 shows the relationship between graph indices and depressive symptoms. The Spearman's rank correlation analysis showed that the number of intra-community edges in the left LM at the high beta band was positively correlated with the HAMD-17 (r = 0.519, p = 0.032; Figure 4, left) and BDI-II (r = 0.522, p = 0.031) scores within the BD group (Figure 4, right). In contrast, age, drug dosages, and the JART, GAF, ASRM, and YMRS scores had no significant relationship with the graph indices.

DISCUSSION
In this study, we performed an analysis using graph theory and observed frequency-specific alterations of RSNs in patients with BD compared to HCs. We observed higher LM community edges at higher frequency band ranges, especially inter-community edges of the right LM-right DMN at the gamma band, and intracommunity edges of the left LM at the high beta band, in the BD group. We also observed decreased inter-community edges at lower frequency band ranges, especially in those of the right SAL-right FPN at the delta band and the left SAL-right SAL at the theta band, in the BD group (Figure 3). The increase in the number of connections in the left LM at the high beta band was positively correlated with subjective and objective depressive symptoms in the BD group (Figure 4).
Several fMRI studies investigating the resting-state functional connectivity in the BD population have shown altered connectivity between parts of the DMN and the limbic system (14,22,(53)(54)(55)(56)(57)(58). The hyper-connectivity between the right LM and the right DMN of the BD group in the present study supports these previous findings and may indicate that alterations in resting-state functional connectivity between the LM and DMN are the essential elements associated with BD. The LM deals with emotional regulation/processing and memory, and previous studies have reported that patients with BD showed increased functional connectivity in the mesoparalimbic network (59,60). The DMN is thought to be involved in self-referential activities such as inner speech, mental images, remembering past events, and planning future events (61,62). The higher intercommunity edges of LM-DMN suggest that self-referential activity tends to be affected by emotional activity, and we speculate that this tendency induces ruminative thoughts in patients with BD.
We also observed attenuation of the number of intercommunity edges of the SAL-FPN at the delta band ( Figure 3A, top left) and left SAL-right SAL at the theta band ( Figure 3A, top middle) in the BD group. A previous study also reported a significantly altered functional connectivity between the bilateral dorsal anterior insula and the left inferior parietal lobe in patients with BD (63). The SAL has been implicated in playing an important role in switching between the central executive network and DMN (64). Furthermore, the SAL activates the central executive network to manipulate information stored in working memory, while suppressing the DMN to keep attention focused on task-relevant goals. The hypoconnectivity between the SAL-FPN and left SAL-right SAL may contribute to the difficulty in maintaining a cognitive set and to the focusing of attention on task-relevant goals in patients with BD.
Remarkably, we found that the number of connections in the left LM at the high beta band ( Figure 3A, lower middle) had positive correlations with the subjective and objective depressive symptoms in the BD group. This finding is consistent with previous reports that alterations in connectivity in the beta band are associated with depressive symptoms in major depressive disorder (65,66), which demonstrates that excessive connectivity of the LM in the beta band is associated with depressive symptoms. Salvadore et al. investigated the patterns of the anterior cingulate cortex (ACC) and amygdala activity in patients undergoing antidepressant treatment with ketamine using the task of fearful faces (65) and found that pretreatment amygdala activity was negatively correlated with changes in depressive symptoms. These results suggest that LM activity increases with an increase in the severity of depressive symptoms and that the connectivity of the LM in the beta band may be a good objective biological biomarker of the depressive symptoms.
A resting-state EEG study in a BD population found frequency-specific alterations of the network topologies (67). Although no resting-state MEG study investigating functional connectivity has been reported so far, frequency-specific intercommunity edges, including the SAL change in the delta and theta bands in patients with BD, were clearly demonstrated in our study. These differences in inter-community edges in the SAL were observed as changes in the low-frequency bands such as delta and theta. A recent study reported that greater theta band connectivity between the right ACC and right anterior insula, a key region within the SAL, predicted a greater reduction in the severity of depression (68). Altered networks at the theta and delta bands have been reported to be associated with pain control and auditory hyperacusis (69,70). Parker reported that patients with BD exhibited various suprasensory changes experienced during manic/hypomanic "high" states, which generally attenuate or disappear during depressive and euthymic periods (71). Network alternations at the theta band in the SAL may be related to the altered sensory perception of BD. Furthermore, patients with BD experience significantly increased levels of pain, particularly chronic pain and migraine. Difficulty in controlling pain in patients with BD may be related to a slow change in the oscillatory network in the SAL. In future studies, we aim to investigate the relationship between symptom formation and frequency-specific networks in BD.

LIMITATIONS
Unfortunately, statistical analysis with a multiple comparison correction was not possible due to the insufficient number of subjects. Patients with BD in a hypomanic or manic state were not included; therefore, an investigation of the relationship between the manic symptoms and the connectome was not performed. Accordingly, large-scale prospective studies using MEG measurements in patients with BD are needed.

CONCLUSION
In conclusion, we observed frequency-specific alterations of RSNs in patients with BD via MEG analysis using graph theory. To the best of our knowledge, this is the first resting-state MEG study investigating functional connectivity using graph theory in a BD population. We show that frequency-specific alterations in the inter-community and intra-community edges in the BD group, and the intra-community edges in the left LM in the high beta band were significantly positively correlated with subjective and objective depressive symptoms in the BD group. We suggest that the connectivity of the LM in the beta band may be a good objective biological biomarker of depressive symptoms in patients with BD. Also, the higher inter-community edges of LM-DMN suggest that self-referential activity tends to be affected by emotional activity in patients with BD. Moreover, the hypoconnectivity between the SAL-FPN and left SAL-right SAL may contribute to the difficulty in maintaining a cognitive set and to the focusing of attention on task-relevant goals in the BD group, and may be related to the altered sensory perception of BD. Further studies based on MEG analysis using graph theory would help elucidate the frequency-specific alterations of RSNs in psychiatric patients and identify biological biomarkers for psychiatric disorders.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by Ethics Committee of the Gunma University Graduate School of Medicine, Japan. The patients/participants provided their written informed consent to participate in this study.

AUTHOR CONTRIBUTIONS
MT, MF, and YuiT designed the study. MS wrote the initial draft of the manuscript. YK, MT, and NH contributed to the analysis and interpretation of data and assisted in the preparation of the manuscript. All authors contributed to the article and approved the submitted version.