Functional network characteristics in anxiety-and mania-based subgroups of bipolar I disorder

.


Heterogeneity of bipolar disorder
Bipolar disorder I (BD-I) is a psychiatric condition characterised by manic and depressive episodes, and at least one manic episode is required to meet the diagnostic criteria (Phillips and Kupfer, 2013).It has long been suggested that BD-I is a heterogeneous disorder, and previous research has proposed that a subset of BD-I patients has a more manic presentation of the disorder (Angst et al., 2004).Furthermore, different combinations of manic symptoms can lead to the diagnosis of BD-I.However, research into subgroups based on mania symptoms in BD-I is limited.Most subgroup analyses in bipolar disorder conducted so far have been based on cognitive performance (Burdick et al., 2014;Lima et al., 2019), psychotic symptoms (Bora, 2018;Markota et al., 2018) and age of onset (Bolton et al., 2021).
Moreover, anxiety disorders are notably more prevalent in BD-I patients than in the general population, with close to half of BD-I patients experiencing an anxiety syndrome in their lifetime (Vázquez et al., 2014).Anxiety is associated with symptoms such as worry about potential future negative events, increased muscle tone, and avoidance (Craske et al., 2009).Comorbid anxiety is suggested to have a significant effect on the course of illness and treatment outcomes in BD-I.Vázquez et al. (2014) reported that for BD-I patients, the presence of anxiety symptoms or syndromes increases the likelihood of hospitalisation, poor treatment response, substance abuse, as well as suicidal ideation and suicide attempts.The researchers noted that our understanding of anxiety morbidity in BD is currently very limited and argued that an improved understanding of anxiety symptoms in BD could aid the selection of more targeted treatment.

Precision psychiatry and symptom-based subgroups
Psychiatric disorders are often considered to be homogeneous and distinct clinical conditions.According to Feczko and Fair (2020), this 'homogeneity assumption' during the diagnostic classification of psychiatric illnesses leads to two major issues: the heterogeneity problem and the comorbidity problem.On the one hand, distinct diagnostic entities can show considerable overlap in symptoms as well as neurobiological signatures, which alludes to a possible shared underlying mechanism or a common psychopathology factor (p) (Caspi et al., 2014), referred to as the comorbidity problem.On the other hand, different underlying factors may drive the same clinical manifestations, which is described as the heterogeneity problem.As mentioned by Marquand et al. (2016), heterogeneity in psychiatry is especially challenging since no neurobiological assessments are currently available to make meaningful individual-based inferences or identify meaningful subgroups within clinical conditions.A single diagnostic category can thus become an umbrella term for many different clinical constructs, which can be driven by divergent biological mechanisms.These challenges illustrate the limitations of case-control studies investigating biological substrates of psychiatric illness and the limited applicability of their research findings in clinical settings at the individual level.Possible ways to tackle the heterogeneity of clinical conditions include identifying subgroups in psychiatric populations or employing more dimensional, symptom-based approaches (Feczko et al., 2019).As proposed by Feczko and Fair (2020), identifying subgroups in diagnostic labels could benefit treatment selection in clinical settings as well as our biological understanding of pathophysiological mechanisms involved in psychiatric disease.Furthermore, subtyping may improve performance of prediction methods for symptom severity, such as in Autism Spectrum Disorder (ASD) as reported by Hong et al. (2018).Altogether, clinical and neurophysiological subgroups may shed light on dimensions of heterogeneity in psychiatric populations that are otherwise obscured.

Brain network disturbances in BD-I
Previous research has reported disturbances in functional brain network properties in bipolar disorder as assessed using graph theory approaches.Graph theory is a mathematical framework which is extensively used to quantify and analyse complex network topology (Bullmore and Sporns, 2009), including brain network interactions measured by resting-state functional magnetic resonance imaging (rsfMRI).Disruptions of large-scale neural networks could potentially reveal some underlying neural pathophysiology (Hulshoff Pol and Bullmore, 2013).A higher clustering coefficient and lower global efficiency in the right amygdala have, for instance, been reported in bipolar disorder compared to controls (Spielberg et al., 2016).In their study, Spielberg et al. (2016) reported that depressed mood was associated with lower assortativity, a measure of the correlation between the number of connections (degree) of connected nodes (Bullmore and Sporns, 2009).Elevated mood, on the other hand, was related to hyperconnectivity in the amygdala and a midbrain region encompassing the substantia nigra, red nucleus, and ventral tegmental area.
These findings hint at the possibility that different network disturbances are at play for the variety of clinical symptoms of bipolar disorder.Literature investigating the interactions between mania symptoms and topological properties of brain networks is limited.
Research by Spielberg et al. (2016) showed that hypomanic symptoms are associated with alterations of network parameters including the small-world network organisation (namely lower global efficiency), increased clustering around the amygdala, and heightened connectivity of the left posterior superior frontal gyrus (pSFG), especially with the association cortex.
The relationship between anxiety symptoms and network properties in BD has been investigated in a study by Lin et al. (2020).In their study, they subdivided a group of BD patients into a subgroup with and a subgroup without anxiety, based on a commonly used cut-off score of points on the Hamilton Rating Scale for Anxiety (HAM-A).Compared to the subgroup without anxiety, the BD group with anxiety symptoms showed significantly lower assortativity and characteristic path length, which the researchers suggested could be indicative of a more vulnerable and more randomised brain network (Lin et al., 2020).Given the limited body of research investigating anxiety in BD-I on a network level, the biological correlates of the clinical manifestation of different anxiety symptoms in BD-I remain unclear.

Aim of the study
To the best of our knowledge, no study to date has investigated the relationships between anxiety and mania symptoms and graph network characteristics in BD-I using data-driven clustering methods.More symptom-based approaches could improve our understanding of the biological mechanisms of mania and anxiety symptoms in BD-I on the brain network level.Therefore, this study investigated whether brain network characteristics differ between data-driven mania and anxiety subgroups in BD-I.

Participants
The data for the current study were obtained from the Psychosis and Affective Research Domains and Intermediate Phenotypes (PARDIP) study 2 through the NIH Data Archive Repository.Recruitment took place with the help of newspaper advertising at three study sites across the United States (Beth Israel Deaconess Medical Centre, University of Texas Southwestern Medical Centre, and Hartford Hospital).To meet the inclusion criteria, participants had to be aged between 18 and years old, and participants of the bipolar group had to meet the DSM-IV criteria for bipolar I disorder.Healthy control subjects were included if they had no history of a mood or psychotic illness and no first-degree relatives with a history of any psychotic or recurrent mood disorder.Exclusion criteria included the presence of a medical disorder that could affect CNS functioning, any major cognitive or neurological disorder, pregnancy, substance abuse or extensive history of substance use, and an estimated IQ score of lower than 70.Participants in the bipolar group were receiving stable drug treatment for at least a month and were considered clinically stable.Psychotropic medication was recorded and included in the analysis as a covariate.
Selecting volunteers with a resting-state fMRI and T1-weighted scan as well as clinical data available from both the Clinical Anxiety Scale (CAS) and the Young Mania Rating Scale (YMRS), a total of 104 participants with BD-I and 70 healthy control participants could be included from the PARDIP study.The bipolar group consisted of patients with a history of psychosis (BDP) and 57 patients without a history of psychosis (BDNP), determined on the basis of the Structured Clinical Interview for DSM-IV (SCID-IV).Given that this study did not investigate the role of psychosis in bipolar disorder, but rather investigated differences between data-driven anxiety and mania clusters, the group of bipolar patients is as of this point referred to as the BD-I group.For further demographic information see Table 1.

Ethics
The research was reviewed by the Cambridge Psychology Research Ethics Committee (PRE.2020.104)and deemed not to require PRE oversight as it consists of secondary analysis of de-identified primary datasets.Informed consent of participants (or their guardians) in primary studies was provided for the primary study from which this study has drawn data.

MRI data acquisition
At all three study sites, the resting-state and structural T1-weighted scans were acquired using 3T scanners, namely scanners by Siemens Allegra and Skyra (Hartford), Philips Achieva (Dallas), and GE Signa HDxt and Discovery MR750 (Boston).Slice thickness was 4 mm, and repetition time (TR) differed across subjects, there was a range of 1.5 s, 2 s, and 3 s.The echo time (TE) was either 27 ms or 30 ms.The flip angle was 60 • for the Siemens and Philips scanners, it was 90 • for the GE scanners.

Clinical clusters
The two clinical assessments used for the creation of clusters were the Clinical Anxiety Scale (CAS) and the Young Mania Rating Scale (YMRS).The CAS consists of six point-based subscores, the YMRS consists of eleven.Rather than taking only the total or average clinical score of these two scales for data-driven clustering, all subscores were included in our clustering analyses.Uniform Manifold Approximation and Projection (UMAP) was used as a dimension reduction method to map subjects in 2-dimensional space (McInnes et al., 2020).UMAP is an algorithm that can represent higher dimensional data on a lower dimensional scale.It is built on complex mathematical theory that is described elsewhere (McInnes et al., 2020).In short, data is first represented in what is called a fuzzy topological structure, which we can understand in simple terms as representing data in a way that extends beyond a binary value such as 0 or 1 and allows data to take a value in between.This fuzzy representation can then be reduced to a low-dimensional structure using mathematical techniques from Riemannian geometry and algebraic topology.UMAP is favoured over dimension reduction techniques such as t-SNE and PCA, because it maintains the local and global structure of the data provided (Dalmia and Sia, 2021).
After mapping subjects on this 2-dimensional scale, clusters were identified using partitioning around medoids (PAM) (Kaufman and Rousseeuw, 1990).This technique starts by randomly identifying data points as centers in your data, called 'medoids'.The algorithm then iteratively evaluates the distance from a data point to a medoid and computes the cost change of swapping the current medoid with this data point.To this end, the data point at minimal distance to other data  Numeric data are reported as mean (standard deviation).BDNP = bipolar disorder without psychosis, BDP = bipolar disorder with psychosis, HC = healthy controls.
points in the cluster, the 'ideal center', will be identified as the final medoid.The number of clusters k one wants to identify needs to be specified a priori.In our approach, we ran the algorithm for a k between 2 and 10, after which the ideal number of clusters was selected based on the highest Silhouette scores (Rousseeuw, 1987).Clustering results were assessed for consistency using a bootstrapping method.Specifically, a sample of 80 % of the BD-I population was taken five times and the same clustering approach was run on each subset.We then calculated an average percentage overlap score.

fMRI pre-processing
First, the raw fMRI data were converted to BIDS format, after which they were pre-processed using the fMRIPrep (Esteban et al., 2019(Esteban et al., , 2020;;version 20.2.0).Pre-processing steps included tissue segmentation, slice time correction, co-registration to the T1 scan, head-motion estimation, Susceptibility Distortion Correction (SDC), and confounds estimation.Component-based noise correction (CompCor) was performed after high-pass filtering of the time-series for the temporal (tCompCor) and anatomical (aCompCor) variants.For details please refer to the methods summary provided by fMRIPrep in the Supplementary Material.

Network construction
Denoised fMRI data were then parcellated by warping the Human Connectome Parcellation template (HCP-MPM-2.0;Glasser et al. (2016)) to each individual's subject space (i.e.mapped to volume space) and average time-series were estimated within each region identified in the parcellation.This resulted in a time-series matrix of 180 regions (or 'nodes') per hemisphere.Edges were calculated by the estimation of pairwise associations between nodes (i.e., Pearson correlations between average nodal time-series) generating correlation matrices.These weighted matrices were thresholded using a combination of computing the minimum spanning tree to ensure consistent matrix size and thresholding to include the top 10% strongest connections (see the Supplementary Material for analyses across multiple thresholds).Graph metrics were then calculated using the Brain Connectivity Toolbox (Rubinov and Sporns, 2010, https://sites.google.com/site/bctnet/).Characteristic path length, clustering coefficient, participation coefficient, assortativity and transitivity were computed as the global functional network measures and strength, degree, and participation coefficient were computed as regional functional network measures to be compared between subgroups and with healthy controls.These network measures have been described in detail elsewhere (Rubinov and Sporns, 2010).Characteristic path length, a measure of integration, refers to the number of average edges from one node to all other nodes.Clustering coefficient is a measure of functional segregation that describes the interconnectedness of neighbours of a node (i.e. to what extent neighbours of a node are also neighbours of one another), as does transitivity, though this measure is normalised collectively rather than individually.Participation coefficient refers to the diversity of a node's connections to different modules and is a measure of centrality.Assortativity is a measure of resilience that describes the correlation between the number of connections (degree) of connected nodes.Regional strength is the sum of the weights of all the node's edges and regional degree is the total number of edges of a node.Network metrics were standardised by dividing them by metrics obtained from random networks.

Statistical analysis
Firstly, the effects of diagnosis and sex on clustering results and the correspondence between anxiety and mania clustering were determined using Fisher's exact test.Two-Sample Kolmogorov-Smirnov Tests were used to assess if there were significant differences between the mean connectivity distributions of the identified subgroups.After this, linear regression was used to investigate functional network differences between subgroups while controlling for a number of covariates in four different models.In model 0, only sex, age, and study site were controlled for, in model 1, sex, age, study site, and framewise displacement were controlled for (technical model), in model 2, sex, age, study site, and medication use were controlled for (clinical model), and in model 3, sex, age, study site, framewise displacement, and medication use were controlled for (full model).Medication was controlled for using a binary variable, in which anti-epileptic medication, antipsychotic medication, antidepressants, benzodiazepines, lithium, anxiolytic medication, and stimulants, were considered as medication.Post-hoc ttests with Benjamini-Hochberg False Discovery Rate (FDR) correction were used to test pairwise differences between specific subgroup pairs.FDR-correction was also applied for the regional analyses, across all regions per network measure.P FDR values represent the q values.See Fig. 1 for an overview of the methods pipeline.

Cluster formation using UMAP and PAM
UMAP's projection of subjects in 2-dimensional space based on their CAS scores revealed two subgroups, see Fig. 2A.Data-driven clustering of these 2-dimensional projections using PAM also showed that two subgroups fit the data best, with a Silhouette score of 0.91.There was no significant association between subgroup label and sex (P = .39,Fisher's exact test), nor was there a significant association between subgroup label and history of psychosis (P = 1, Fisher's exact test), see Supplementary Figure A1.For the YMRS scores, UMAP's projection of subjects in 2-dimensional space revealed four subgroups, see Fig. 2B.Using PAM, the highest Silhouette score (0.81) was found for four subgroups, demonstrating that this number of clusters fit the data most appropriately.Again, there was no significant association between cluster label and sex (P = .46,Fisher's exact test).There was however a significant association between cluster label and history of psychosis (P < .001,Fisher's exact test), with cluster 1 only containing subjects with a history of psychosis.For details, see Supplementary Figure A2.See Supplementary Figures A7-A8 for plots of the first two principal components for the anxiety-and mania-based subgroups.Both the anxiety and mania subgroups showed high consistency when validated using five rounds of bootstrapping, for details see the Supplementary Material.

Differences in symptoms between subgroups
Overall, anxiety subgroup 0 generally scored higher on CAS symptoms than anxiety subgroup 1, see Fig. 2C.The most notable differences between the anxiety subgroups were found for psychic tension, muscular tension, and worrying.As for the mania subgroups, symptom scores show they were not defined by increasing YMRS symptom scores, but rather qualitatively differed from one another.The most notable differences can be seen on YMRS symptom scales for (speech) content, irritability, and speech (rate and amount), see Fig. 2D.Subjects in mania subgroup clearly differed from the others based on their high (speech) content scores, whereas mania subgroup 2 was most clearly distinguished from the rest by high scores for speech (rate and amount).Mania subgroup 0 scored relatively low on all three of these YMRS symptom scales, whereas mania subgroup 3 generally had low scores for (speech) content and speech (rate and amount), yet scored higher for irritability.Furthermore, there was a significant association between the anxiety-based subgroups and the mania-based subgroups (P < .001,Fisher's exact test).Mania subgroups 1, 2, and 3 all showed considerable overlap with anxiety subgroup 0, and most subjects in anxiety subgroup 1 were assigned to mania subgroup 0. For details see Supplementary Figure A9.

Network analysis findings
Based on our clustering analysis, there were 71 subjects in anxiety cluster 0 and 33 in anxiety cluster 1.There were 41 subjects in mania cluster 0, 14 in mania cluster 1, 14 in mania cluster 2, and 35 in mania cluster 3.There were 70 healthy control participants.Supplementary Figures A15 and A16 show the unthresholded matrices for all anxiety and mania subgroups, Figure A17 and A18 show the thresholded matrices after computing the minimum spanning tree and applying a cut-off point of the top 10 % connections.The mean connectivity distributions for each subgroup significantly differed from each other as measured by a Two-Sample Kolmogorov-Smirnov Test, for details see the Supplementary Material.

Global functional network characteristics
Multiple linear regression was used to test if anxiety subgroups were significantly associated with the functional network measures characteristic path length, transitivity, assortativity, clustering coefficient, and participation coefficient.The model was scaled by subtracting the mean network parameter score in the formula, by which the intercept could more clearly be interpreted.Age, sex, framewise displacement, study site, and medication were used as covariates in the full model (model 3).It was found that anxiety subgroup labels were significantly associated with assortativity (Tables 2 and 3).No significant main effect of anxiety subgroups on characteristic path length, transitivity, clustering coefficient, and participation coefficient was found.When comparing these findings to model 0 (controlled for age, sex, and study site), model 1 (technical model), and model 2 (clinical model), assortativity remained the only measure for which a significant main effect was found.Anxiety subgroup 0 was associated with participation coefficient scores (β = − 0.06, t(164) = − 2.71, P < .01).Furthermore, we found that age was significantly associated with participation coefficient scores, β < 0.01, t (164) = 2.03, P = .044.Despite harmonised processing, there were some lingering effects of study site on assortativity (β = − 0.03, t(164) = − 2.16, P = .032).It was also found that motion, quantified using framewise displacement scores, was associated with characteristic path length (β = 0.01, t(164) = 2.65, P < .01),transitivity (β = 0.05, t(164) = 2.09, P = .038),and clustering coefficient (β = 0.03, t(164) = 2.06, P = .041).
For the mania-based subgroups, subgroup assignment was not associated with any network metrics when controlling for age, sex, framewise displacement, study site, and medication (full model, model 3).When controlling for model 0 (controlled for age, sex, and study site), model 1 (technical model) and model 2 (clinical model), mania Fig. 2. Fig. 2A shows the clustering results based on the CAS scores (Silhouette score = 0.91).Fig. 2B shows the clustering results based on the YMRS scores (Silhouette score = 0.81).The axes in Fig. 2A and 2B represent the two dimensions to which the multidimensional data was reduced using the UMAP technique (rather than representing specific clinical constructs or scales) and the axis units are arbitrary.Fig. 2C shows the distribution of three anxiety symptoms from the CAS for the two anxiety-based subgroups shown using a box plot and violin plot.These three rating scores showed the biggest differences between the subgroups, the plots for the other three scores can be found in Supplementary Figure A5.Fig. 2D shows the distribution of three mania symptoms from the YMRS for the four mania-based subgroups shown using a box plot and violin plot.These three rating scores showed the biggest differences between the subgroups, the plots for the other eight scores can be found in Supplementary Figure A6.subgroups were significantly associated with assortativity.Multiple linear regression also showed that mania group 2 significantly predicted participation coefficient (β = − 0.04, t(162) = 2.02, P = .045).Again, study site had a significant impact on assortativity (β = − 0.03, t(162) = − 2.09, P = .038).Framewise displacement significantly affected characteristic path length (β = 0.01, t(162) = 2.62, P = .010),transitivity (β = 0.05, t(162) = 2.13, P = .035),and clustering coefficient (β = 0.03, t (162) = 1.98,P = .049).Table 2 summarizes linear regression findings for the full model (model 3), linear regression findings from the other three models can be found in the Supplementary material.
Based on these results, we conducted an exploratory post-hoc analysis for differences in assortativity between anxiety subgroups.Post-hoc comparisons using linear regression with the same covariates as listed above and corrected with Benjamini-Hochberg False Discovery Rate (FDR) correction revealed no significant differences between the anxiety subgroups (Fig. 3).

Regional functional network characteristics
Regional analysis using univariate linear regression for each cortical region showed that anxiety subgroups differed in regional strength in 24 regions, for regional degree in 15 regions, and for regional participation Significance codes: 0.01 < '*' < 0.05. in 5 regions (Fig. 4A-C, see Supplementary Table 1 for detailed linear regression results).Mania subgroups differed in strength in 24 regions, in degree in 10 regions, and in no regions for participation coefficient (Fig. 4D-F, see Supplementary Table 2 for detailed linear regression results).For more details on pairwise comparison between subgroups within the significant regions, see Supplementary Figures A27-A31.

Discussion
Taking a data-driven approach, this study identified subgroups based on anxiety and mania symptoms in a population of BD-I patients (i) and showed that anxiety subgroups were significantly associated with global assortativity but not global degree, transitivity, characteristic path length, and clustering coefficient, and mania subgroups were not associated with any global network metrics (ii).On a regional level, anxiety subgroups were associated with differences in strength, degree, and participation coefficient.Mania subgroups were associated with regional strength and degree, but not regional participation coefficient.

Symptom-based subgroups in BD-I based on anxiety and mania scores
Clustering analysis revealed two distinct subgroups in the BD-I group based on anxiety symptoms.The largest subgroup showed a general trend of higher anxiety symptom scores, which is in line with the high prevalence of anxiety comorbidity in bipolar populations (Vázquez et al., 2014).Approximately two-thirds of the study population was divided into the higher anxiety subgroup, which is more than the 48.2 % prevalence of anxiety disorders reported by Vázquez et al. (2014).This difference could reflect a bias in the initial recruitment process, or it could be that not all subjects in this subgroup would meet all diagnostic criteria for an anxiety disorder despite experiencing some anxiety symptoms.The finding of an anxiety-based subgroup in BD-I is in line with a prior study by Qiu et al. (2017), which also reported subgroup identification using anxiety-related factors.This supports the notion that there may be a distinct subgroup of individuals with BD-I that suffers from elevated anxiety symptomatology.In their work, Qiu et al. (2017) suggested that a low-anxiety subgroup may be associated with a more stable disease course and lower phenotype expression of BD-I compared to a high-anxiety subgroup.
The four mania-based subgroups identified in this study showed distinct mania symptom profiles, especially regarding speech and irritability.Compared to anxiety-based subgroups, literature on stratification of BD-I populations based on mania symptoms is limited.Work by Haro et al. (2006) identified three subclasses in a large population of more than three thousand patients treated for an acute manic episode using data-driven methods, which they described as 'typical mania', 'dual mania', and 'psychotic mania'.The psychotic mania group was found to have a more severe presentation of manic symptoms and more impairments in social functioning.Our findings for subgroup 1 are in line with their findings, given that subgroup 1, which consisted only of subjects with a history of psychosis, scored higher on mania symptoms such as elevated mood, motor activity, sleep disturbance, speech content, and language symptoms compared to the other subgroups.It would be of interest to investigate whether differences in social impairment could be seen between the four subgroups described in this study.Altogether, more research is needed to confirm or reject the current classification of a BD-I population into four subgroups on the basis of mania symptoms, but these findings suggest mania and psychotic symptoms could serve as indicators to divide BD-I patients into more homogeneous subpopulations that may be clinically relevant.4A-C show regional differences in strength, degree, and participation coefficient comparing anxiety subgroup 0, subgroup 1, and healthy controls to each other.Fig. 4D-F show regional differences in strength, degree, and participation coefficient comparing mania subgroup 0, subgroup 1, subgroup 2, subgroup 3, and healthy controls to each other.The blue colour represents the F statistic, orange indicates whether a region showed a significant difference in a network parameter across subgroups after multiple testing correction (FDR) with significance level P FDR < 0.05.

Network connectivity and graph theoretical measures
Linear regression revealed that only assortativity was significantly associated with anxiety subgroup label.Assortativity is a measure of the extent to which similar nodes are connected in a graph and has been proposed as an indicator of the resilience of a network (Spielberg et al., 2016).A prior study by Lin et al. (2020) showed that an anxiety subgroup of bipolar patients had significantly lower assortativity and characteristic path length than the subgroup without anxiety.However, in our study, linear regression revealed no association between anxiety subgroups and characteristic path length and post-hoc pairwise comparison did not reveal any significant differences between the anxiety-based subgroups for assortativity.Given that our subgroup sample sizes were small, larger sample sizes could have aided detection of differences in assortativity during pairwise comparison.Furthermore, differences between anxiety subgroups were much more apparent on the regional cortical level; strength and degree both showed differences in frontal and temporal regions, such as clearly seen in the right and left orbitofrontal cortex (OFC).Participation coefficients differed between anxiety subgroups in five regions, mainly located in the posterior cingulate cortex (PCC).Given that, to our knowledge, no other study has investigated regional differences in strength, participation coefficients, and degree in anxiety subgroups of bipolar disorder, these findings need replication in future studies and need to be interpreted with caution.
For the mania subgroups, the analysis yielded different results than the findings reported by Spielberg et al. (2016).In their study, they found that YMRS scores, and irritability in particular, were associated with lower global efficiency.This measure is inversely related to characteristic path length, but our linear regression did not reveal any significant differences in this global network property between mania subgroups.However, our sample sizes for the mania subgroups were small, especially for subgroups 1 and 2, which were most notably characterised by irritability, limiting direct comparisons.It should be noted that Spielberg et al. (2016) did not define YMRS-based subgroups, but only explored associations between network properties and symptoms.Again, differences between mania subgroups were more apparent on the regional level, with strength and degree showing significant differences in a large number of regions mainly located in temporal and frontal regions.No regional differences in participation coefficient were found between mania subgroups.Further research that stratifies BD-I populations into subgroups based on mania symptoms is necessary to replicate global and regional network findings of this study.Nonetheless, these findings give a first insight into understanding the neural correlates of symptom profiles within BD-I.

Limitations and future directions
There are a number of limitations to be taken into account when interpreting the findings of this study.Firstly, there was an issue of power, due to the small number of subjects in some of the subgroups.Clustering results need to be validated in a larger sample which could also aid the detection of network differences between subgroups.Secondly, the symptom-based subgroups identified in this study provide a categorical representation of mania and anxiety symptomatology in BD-I, whereas these symptoms might be understood better as existing on a continuous spectrum.One way to incorporate a more dimensional approach is to employ normative modelling.Mapping anxiety and mania symptomatology in BD-I onto a distribution that includes healthy participants may give more insight into the biological mechanisms involved, given that disease variation may overlap with population variation.Lastly, findings from this work represent a clinically stable population of BD-I patients, limiting the interpretability of findings and the extent to which they can be translated to a larger population of BD-I patients.Findings would need to be replicated in future work that would include BD-I populations with greater variance in clinical stability.
Future work could explore whether anxiety-and mania-based subtypes of BD-I respond differently to treatment, which could be informative for precision medicine in psychiatry.This could allow clinicians to choose which therapeutic intervention to start with based on the combination of symptoms a patient is experiencing, thereby avoiding long waiting times for patients due to ineffective treatments.Furthermore, future studies could investigate these subtypes in transdiagnostic cohorts, to see whether findings extend to other psychiatric conditions.For instance, future studies could explore whether datadriven clustering analysis reveals similar subgroups in MDD and BD-II.Lastly, it could be investigated whether mania-and anxiety-based subgroups correlate with other clinical features such as age of onset, cognitive functioning, and suicidality and other symptom dimensions could be considered in subgroup identification, such as a clinical assessment of depression.

Conclusion
This was the first study to investigate data-driven subgroups of BD-I based on mania and anxiety symptoms.The findings show that two robust anxiety-and four robust mania-based subgroups could be identified in this study population, which showed functional network differences that were most apparent on a regional level.The present study and its findings add to the body of research devoted to exploring more symptom-based approaches to psychiatric conditions and improving our understanding of different clinical profiles within these conditions.

Fig. 1 .
Fig. 1.Overview of the pipeline of the project.

Fig. 3 .
Fig. 3. Global network properties for the two anxiety-based subgroups (I) and the four mania-based subgroups (II).Post-hoc pairwise comparison showed no significant differences between anxiety subgroups.Abbreviations: cpl = characteristic path length, cc = clustering coefficient, pc = participation coefficient.

Fig. 4 .
Fig.4.Fig.4A-C show regional differences in strength, degree, and participation coefficient comparing anxiety subgroup 0, subgroup 1, and healthy controls to each other.Fig.4D-F show regional differences in strength, degree, and participation coefficient comparing mania subgroup 0, subgroup 1, subgroup 2, subgroup 3, and healthy controls to each other.The blue colour represents the F statistic, orange indicates whether a region showed a significant difference in a network parameter across subgroups after multiple testing correction (FDR) with significance level P FDR < 0.05.

Table 1
Demographic characteristics of the study population.

Table 2
Multiple linear regression results for anxiety and mania subgroups for model 3 (adjusted for clinical and technical covariates).
Adriana P. C. Hermans et al.

Table 3
Multiple linear regression results for assortativity in anxiety subgroups for model 3 (adjusted for clinical and technical covariates).