The effect of using group-averaged or individualized brain parcellations when investigating connectome dysfunction in psychosis

Psychosis is a neuropsychiatric condition that has long been thought to arise from aberrant neural connectivity, or dysconnectivity, between neuronal populations (Andreasen et al., 1998; Baker et al., 2019; Fornito et al., 2012; Nogovitsyn et al., 2022). Such dysconnectivity is often studied using a network-based approach (Fornito et al., 2016), with the brains of individuals being modeled as a collection of nodes, representing discrete brain regions, connected by edges, representing interregional structural connectivity or functional coupling (FC). This approach has revealed extensive FC disruptions in psychosis patients; these disruptions are often characterized by a global decrease in FC upon which is superimposed more network-specific increases and decreases (Argyelan et al., 2014; Baker et al., 2019; Chopra et al., 2021; Fornito et al., 2012; Hummer et al., 2020; T. Li et al., 2017; Narr & Leaver, 2015; Nogovitsyn et al., 2022; Tu et al., 2013). However, the reported findings have been inconsistent, with reports of increased and decreased FC sometimes found within the same network in different samples (Lynall et al., 2010; Moran et al., 2013; Whitfield-Gabrieli et al., 2009; Woodward et al., 2011).

Some of these inconsistencies may be explained by methodological differences in defining the nodes (brain regions of interest, ROIs) of the constructed brain networks, which is a fundamental step in network analysis that could affect the validity and interpretation of subsequent results (Fornito et al., 2010, 2016; Zalesky, Fornito, Harding, et al., 2010). Each node should ideally represent a functionally specialized area with homogenous activity (Eickhoff, Constable, & Yeo, 2018; Eickhoff, Yeo, & Genon, 2018), but there is no consensus on the optimal way of parcellating the brain, meaning that investigators must rely on various heuristic methods (Eickhoff, Constable, & Yeo, 2018; Eickhoff, Yeo, & Genon, 2018).

The vast majority of studies in patients with psychosis have used a one-size-fits-all, group-based approach in defining distinct ROIs. A parcellation using this approach is often defined in a standardized coordinate space based on a sample average and then mapped to individual participants via a spatial normalization procedure (Eickhoff, Yeo, & Genon, 2018). This approach fails to consider interindividual variability in functional and anatomical brain organization (Amunts et al., 2005; Mueller et al., 2013). Investigation of such variability with resting-state fMRI (rsfMRI) has shown that, although most cortical areas can indeed be robustly identified in every individual, their sizes and shapes vary across the population, especially when using more fine-grained parcellation methods (Gordon et al., 2017). Furthermore, the topographical locations of specific areas tend to shift between individuals, sometimes across anatomical landmarks such as sulci and gyri (Gordon et al., 2017), which are often used as reference points in many standard parcellations (Fornito et al., 2016).

To better accommodate this individual variability, approaches have been developed to derive individualized parcellations at the level of either canonical functional networks (S. Li et al., 2016; Yeo et al., 2011) or cortical regions (Gordon et al., 2017; Kong et al., 2021). These approaches have revealed that individual variability can considerably impact network analyses. For instance, regions assigned to one network in individual parcellations are often assigned to a different network in the group average (Bijsterbosch et al., 2018), which could impact FC analysis. The use of individually tailored parcellations yields more functionally homogeneous regions (Chong et al., 2017; Kong et al., 2021) and can improve predictions of behavior from FC (Kong et al., 2019). Indeed, in healthy samples, individual differences in the locations of functional regions, as represented by individualized parcellation, affect predictions of fluid intelligence (Kong et al., 2019), life satisfaction (Bijsterbosch et al., 2018), participant sex (Salehi et al., 2018), and performance in reading and working memory tasks (Kong et al., 2021). Moreover, some estimates indicate that up to 62% of variance in network edge strength (i.e., FC values) can be explained by the spatial variability of defined regions (Bijsterbosch et al., 2018). These findings suggest that clinically important relationships may be masked when using a group-based parcellation. On the other hand, these approaches present several challenges, such as establishing a correspondence between similar regions in different people and accounting for differences in region size.

A particularly salient point in clinical studies, such as those of schizophrenia, is that standard brain atlases have been derived from healthy participants, which may not adequately capture the characteristic properties in the brain organization of patients (Glasser et al., 2016; Schaefer et al., 2018). Patient-specific individual variability in functional organization can influence the results of brain network analyses. Indeed, one study has found that slight displacements of a seed region in the thalamus can lead to significant differences in disorder-related dysconnectivity (Welsh et al., 2010), emphasizing the importance of a valid and consistent node definition.

One strategy to develop individualized parcellations is to adjust the borders of a group-based template for each individual participant according to predefined functional criteria. For instance, Chong et al. (2017) developed a Bayesian algorithm (called group prior individualized parcellation, GPIP) that uses a group-based template as a prior to find an optimal corresponding parcellation on individual brains using individual FC data. The group-based prior ensures that the same regions are mapped in each individual, while updates to the individualized prior account for variability in the shape and size of each parcellated region. Chong et al. (2017) have shown that this method yields parcellated regions with increased intraregional functional homogeneity and reduced variance in connectivity strength between individuals. Here, we used this approach to compare FC disruptions observed in people with early psychosis using analyses that rely on either a group-based or an individualized parcellation. The parcellation algorithm (Chong et al., 2017) allowed us to match all brain regions across participants while accounting for individual variability. Our analyses were conducted using the high-quality, open-access data provided by the Human Connectome Project for Early Psychosis (HCP-EP) resource (Glasser et al., 2013; HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021). We tested two competing hypotheses of how individual variability contributes to apparent FC disruptions in psychosis. Under one hypothesis, a failure to consider individual variability may lead to erroneous regional parcellations, adding noise to the analyses and reducing statistical power for detecting valid group differences. In this case, we expect to see fewer differences between patients and controls when using the group-based parcellation compared with individualized parcellation. Alternatively, FC differences between groups may be largely driven by variations in the underlying organization of each individual’s brain, rather than reflecting specific differences in FC. In this case, we expect to see more differences using the group-based parcellation.

Here, we present results obtained using group-level cortical parcellations provided by Schaefer et al. (2018) as the basis for our analysis, focusing on the 100-region parcellation (s100). To ensure that our results are robust to the number of regions, we repeated our analysis using the 200-region variant (s200) and after applying global signal regression (GSR). Results obtained using the s200 atlas, and results for both atlases after GSR, can be found in the Supporting Information and are largely consistent with the primary results reported in the following sections.

Figure 1 shows examples of individualized parcellations generated for three individuals compared with the original group-based s100 atlas. The individualized parcellation algorithm preserved the same regions for every individual but shifted their borders and changed their shapes and sizes to accommodate for individualized variations in brain organization. Indeed, on average, 42.56% (SD = 2.37) of vertices were reallocated to a different region as a result of the individualized parcellation algorithm, highlighting the considerable variability of cortical functional organization between individuals. Figure 2A shows the proportion of vertices that were relabeled in controls M (SD) = 43.28% (2.34) and in patients M (SD) = 42.20% (2.31). The difference between the two groups was small but statistically significant, following permutation testing (p = 0.004, Hedges’s g = 0.465). However, at a regional level (Figure 2B), no parcel showed significant differences in the number of vertices relabeled between patients and controls (i.e., all p_FDR > 0.05, corrected with the Benjamini and Hochberg method).

We next compared the average functional homogeneity of the group-based and individualized parcellations. Functional homogeneity was measured out of sample, on functional scans from Run 2 with parcellations generated for scans from Run 1. In controls, the mean homogeneity was 0.364 (SD = 0.09) and 0.372 (SD = 0.08) for the group-based and individualized parcellations, respectively. In patients, the mean homogeneity was 0.297 (SD = 0.06) and 0.304 (SD = 0.06) for the group-based and individualized parcellations, respectively (Figure 2C). A two-way mixed ANOVA revealed that mean homogeneity was higher for the individualized parcellation (F(149) = 54.81, p < 0.0001) and higher in controls compared with patients (F(149) = 30.91, p < 0.0001), with no interaction between parcellation type and diagnostic group (F(149) = 0, p = 0.898). Post hoc analysis showed that individualized parcellation resulted in greater homogeneity scores in patients (t(103) = 5.64, p < 0.0001) and controls (t(46) = 2.90, p = 0.006). When comparing homogeneity scores for individual parcels (Figures 2D and 2E), 55 out of 85 regions showed significant differences in homogeneity between parcellation approaches (i.e., p_FDR < 0.05, corrected with the Benjamini and Hochberg method). Moreover, both methods showed high reliability of homogeneity estimates, as measured with the intraclass correlation coefficient (McGraw & Wong, 1996) (r_group-based = 0.842, p < 0.0001 and r_{individualized} = 0.862, p < 0.0001). To quantify functional distinctions between parcels, we computed the mean Pearson’s correlation of activity between each pair of vertices that were not allocated to the same region. We found that the individualized parcellation (M_corr (SD) = 0.100 (0.066)) delineates parcels that are slightly more functionally distinct than those in the group-based parcellation (M_corr (SD) = 0.102 (0.066)). Although small, the difference was statically significant (t(165) = 14.0, p < 0.001).

Homogeneity scores results were similar for the s200 atlas with and without GSR (Supporting Information Figures 2 and 3). For the s100 atlas with GSR, differences in homogeneity between groups and parcellation approach were similar to the main results. However, there was a significant interaction effect between parcellation type and diagnosis (F(148) = 4.68, p = 0.032; see Supporting Information Figure 1), such that homogeneity scores in patients were more impacted by individualized parcellation than in controls. This result suggests that, at this particular resolution, parcellation type may differentially affect FC estimates in patients and controls only following the application of GSR. The reasons for this sensitivity to parcellation scale and GSR are unclear.

Following exclusion of regions with poor signal (see Methods), the final networks examined comprised 85 regions. The FC matrices resulting from both parcellation methods were positively correlated, with correlations ranging between 0.679 and 0.898 (median = 0.794) across participants (Supporting Information Figure 4a), indicating that the results obtained with individualized and group parcellations are generally similar, although far from identical.

Figure 3A shows the distribution of t statistics across edges, comparing FC between patients and controls estimated using either the group-based or the individualized parcellation. Both distributions have predominantly positive values, consistent with evidence of widespread hypoconnectivity in patients compared with controls. The distribution for the group-based approach is shifted further to the right, indicating that larger group differences are detected with this method, on average. The difference in the means of the distributions was statistically significant, as calculated with a Wilcoxon signed-rank test (Z = 24.053, p < 0.0001). Figure 4 of the Supporting Information shows that most FC edges were positively valued; as such, the significant FC reductions observed in patients result from patients having lower positive FC rather than patients having stronger negative FC. Given the higher functional homogeneity of the individualized parcellation, this result suggests that the group-based parcellation overstates FC differences between patients and controls. Similar results were obtained when looking at the effect size of the differences in edge strength between patients and controls (Supporting Information Figure 4), with the group-based parcellation yielding higher effect size estimates than individualized parcellation, on average (p < 0.0001).

The t matrices obtained using the group-based and individualized parcellations were positively correlated (r = 0.76, p < 0.0001), suggesting that the two approaches show largely similar between-group FC differences. The effects of parcellation type were consistent across the full extent of the t distributions, as indicated by the shift function, which compares differences between distributions at each decile. This analysis showed a significantly higher value in every decile of the group-based parcellation, when compared with the individualized parcellation, with the 95% confidence interval never crossing zero (Figure 3B). There was, however, a more pronounced effect of parcellation type on edges associated with larger case-control differences in FC relative to those with smaller case-control differences, as can be seen by the greater shift observed in the right tail of the distribution relative to the left (Figure 3B). This result implies that variations in parcellation type are more likely to influence the edges that are significantly different between patients and controls. Furthermore, results obtained using the s200 parcellations are in agreement with results obtained from the s100 parcellation (see Supporting Information Figure 2). Following GSR, at both parcellation scales, the mean t values were similar, but the t distribution for the individualized parcellation was narrower than for the group-based parcellation. The shift function showed that significant differences between parcellation approaches were mainly for edges with positive t values (see Supporting Information Figures 1 and 3).

We used the network-based statistic (NBS) for inference on the edge-specific t statistics (Zalesky, Fornito, & Bullmore, 2010). The NBS identified a single connected component with significant FC differences between patients and controls using both the group-based (p < 0.0001) and individualized (p < 0.0001) parcellations, for all primary test statistics thresholds tested. Out of 3,570 possible connections, for a primary threshold equivalent to a p value = 0.05, the group-based and individualized parcellations resulted in components comprising 2,877 edges and 2,672 edges, respectively (Figures 4A–B). Thus, the group-based approach implicated approximately 7.7% more dysconnected edges. The binary edge matrices defining these components were moderately correlated (r_phi = 0.548, p < 0.0001) and both components had a total of 571 edges that differed from each other. There was also some variation in the regional affiliation of the edges. For example, Figures 4C–D show that the insula has a high dysconnectivity degree in both group-based and individualized parcellations, but that the former approach implicates more insula subregions. Furthermore, the right medial prefrontal cortex shows a low degree in the individualized parcellation but not in the group-based parcellation. The NBS was repeated with a primary test statistics threshold equivalent to p values = 0.01 and 0.001. For p = 0.01, the component for individualized parcellation comprised 1,786 edges and for group-based parcellation, 2,120. For p = 0.001, the component for individualized parcellation comprised 775 edges and for group-based, 1,257 edges. Note that for all edges in these NBS networks, patients showed reduced FC compared with controls.

A challenge of using individualized parcellations is that the ROIs can vary in size across individuals, which may bias estimates of FC differences between groups. We therefore examined changes in parcel size resulting from the individualization algorithm, as quantified by the number of vertices in each parcel. On average, parcels changed by 50.7 (SD = 45.2) vertices across patients and 52.0 (SD = 45.3) across controls, with no significant difference between the two groups, according to permutation testing (p = 0.104) (Supporting Information Figure 8a). There was also no significant difference in size difference between patients and controls for any of the parcels, when corrected for multiple comparisons following permutation statistics (i.e., all p_FDR > 0.05). Four parcels had different sizes between patients and controls, without correction for multiple comparisons (visual network parcel 9 of the left hemisphere, p = 0.023; somatomotor network parcel 1 of the left hemisphere, p = 0.026; limbic network parcel 1 in the orbital frontal cortex of the left hemisphere, p = 0.039; limbic network parcel 1 in the orbital frontal cortex of the right hemisphere, p = 0.048). We next correlated the differences in parcel size in individualized parcellation between patients and controls with differences in node degree within the NBS network and mean edge dysconnectivity, given by the mean t value of edges attached to each node for the case-control comparison (Supporting Information Figures 8b–c). Neither correlation was significant (r = 0.148, p_spin = 0.104, and r = 0.133, p_spin = 0.127, respectively), suggesting that parcel size did not impact FC differences between patients and controls in the individualized parcellation.

Having demonstrated that the choice of a parcellation strategy can influence both edge- and region-level inferences about FC disruptions in psychosis, we next examined whether parcellation type affects the specific networks that are considered to be dysfunctional. We therefore examined the proportion of edges within the NBS network that fell within and between each of 7 canonical functional networks (Yeo et al., 2011). Considering the raw number of affected edges across both parcellation approaches, the control network was the most impacted in patients with psychosis, with over 1,100 dysconnected edges, particularly those linking the control and somatomotor networks (Figures 4E–F). By comparison, normalized counts, which are adjusted for the total number of possible edges within or between pairs of networks, suggested a more equal and widespread distribution of FC disruptions across networks. Both the raw count (r = 0.983, p < 0.0001) and the normalized matrices (r = 0.802, p < 0.0001) were strongly correlated across the two parcellation methods. These findings indicate that while parcellation method can influence the specific edges that are identified as dysconnected, these edges generally fall within or between the same canonical networks.

Several studies have reported functional brain dysconnectivity in psychosis. A fundamental step in such analyses involves defining a priori ROIs to serve as nodes in the network analysis, which are typically derived from standard parcellation atlases generated from a population or group average template. Here, we asked whether the failure of such an approach to account for individual differences in brain functional organization can bias estimates of case-control differences in FC. Standard methods could result in either an underestimation of the extent of network dysfunction (owing to noisy FC estimation caused by inaccurate ROI delineations) or an inflated estimate of the dysfunction (owing to FC differences being attributable to ROI misalignment), compared with when accounting for individual differences in functional organization of the brain. Our findings indicate that group-based parcellations might inflate estimates of FC differences in psychosis, especially at the edge level. Moreover, the use of individualized parcellations, while yielding a generally consistent pattern of findings, leads to some different conclusions about the specific edges and regions most affected by the disorder, although inferences at the network level were robust to parcellation variations. Together, our findings suggest that the use of individualized parcellations can impact findings of brain dysconnectivity in psychosis and, by extension, other disorders.

The individualized parcellations resulted in nearly half (over 40%) of vertices being assigned to regions that differed from the group-based atlas, as per prior work (Chong et al., 2017). This finding reiterates how group-based parcellations can result in a substantial misspecification of regional borders in individuals and highlights the high degree of variance present in the topographical organization of functional areas. Despite the high percentage of vertices relabeled, FC matrices generated by both parcellations were highly correlated, indicating that the overall FC patterns seen with group-based parcellation are maintained with the individualized parcellation. Note that with GPIP, correspondence between regions is determined based on similarity in FC profiles rather than on spatial location. As such, corresponding regions can shift their spatial location from person to person (see Figure 1).

The higher functional homogeneity of the individualized parcellations supports its improved validity, although the increment was small (2.4%), which is consistent with past reports (Kong et al., 2021; Y. Li et al., 2022), increased homogeneity was seen in the majority of parcels. Regional homogeneity was also marginally (2.3%) higher in controls compared with patients. This differential improvement in homogeneity was expected, as the starting point for the GPIP algorithm was the Schaefer atlas (Schaefer et al., 2018), which was derived from a sample of people with no psychiatric disorders. Defining an initial group atlas in patients would better account for differences in cortical functional organization caused by psychosis. However, it would complicate comparisons between groups because of the requirement to have consistently defined nodes in both patients and controls, which is one of the challenges of using individualized parcellation. The interaction effect between diagnosis and parcellation approach was not significant in most cases (apart from s100 with GSR). This result indicates that individualized parcellations led to a similar improvement in patients and controls. Since most case-control studies use data obtained in healthy individuals to establish a normative benchmark for measures acquired in patients (Chopra et al., 2021; Nabulsi et al., 2020; Nogovitsyn et al., 2022), we relied on the Schaefer parcellation in our analysis. Future work could develop methods to better capture variations in functional organization associated with psychosis.

Widespread decreases in FC in patients with psychosis were identified using both parcellation approaches, highlighting that the dominant effect of both parcellations is generally similar. However, the magnitude of the differences in FC was greater in the group-based parcellation compared with individualized parcellation. Notably, the shift function analysis indicated that differences between the two parcellation approaches were greater for edges associated with large case-control differences. These edges are precisely the ones that are most likely to be declared to be statistically significant following the application of some thresholding procedure. Accordingly, comparison of NBS results revealed a 7.7% reduction in the size of the dysfunctional component identified using the group-based parcellation. Given the higher functional homogeneity, and thus validity, of the individualized parcellation, these results support the hypothesis that at least part of the group differences identified in past studies in psychosis samples do not reflect actual differences in interregional FC but instead result from inaccurate ROI boundaries caused by a failure to account for individual differences in functional organization. These findings imply that individualized parcellations can yield different estimates of FC differences in case-control studies, especially when investigating FC changes at an edge or node level.

While widespread decreases in FC were apparent in patients with psychosis using both parcellation methods, the specific edges affected varied considerably. The NBS components of both group-based and individualized parcellations showed differences in 571 edges (i.e., 19.8% of the total identified with the group-based parcellation). Examining the regions most affected by quantifying the node degrees of the NBS components resulted in broadly similar patterns, but there were some notable differences in location. For example, the right medial frontal region accounts for 1.7% of dysconnectivity in the group-based and 2.3% in the individualized parcellation. The left insula accounts for a slightly smaller percentage (6.5%) of dysconnectivity in the group-based than in the individualized parcellation (6.7%). These findings suggest that conclusions about the specific edges and regions affected by psychosis can vary depending on the parcellation method used. In contrast, inferences at the network level were largely consistent across the two parcellation approaches, indicating that coarse-grained localizations of FC differences are robust to this methodological choice. This could be attributed to network-level inference effectively reducing the dimensionality of the analysis, minimizing the nuances of more fine-grained individual variations. Therefore, studies looking at group differences in FC at a coarse network level might not be impacted by the use of individualized versus group-based parcellations.

To minimize the computational cost, we used fsaverage5, a surface mesh with a relatively low number of vertices. Since GPIP parameters depend on the number of vertices of the mesh, future work could investigate the impact of different surface mesh resolutions and whether the differences observed here apply at different mesh resolution.

To facilitate comparison between subjects, the individualized parcellation algorithm maintains the same number of regions for every subject and uses a parcellation derived in healthy individuals as a starting point. This approach may mask differences in cortical organization in patients, where regions may merge or split, resulting in a different number of ROIs. However, generating separate parcellations in each group complicates comparisons between groups. Resolving this challenge remains an open problem for the field. Moreover, we only looked at cortical regions, owing to the lack of methods available for individual parcellation of subcortical structures.

A proportion of patients in our sample were medicated, and recent evidence has shown that antipsychotic medication can impact FC, even after only 3 months of use (Chopra et al., 2021). However, given that most samples examined in past research are also medicated, our sample is directly comparable to the broader literature. Similarly, the study included more patients than controls and future work could benefit from a balanced sample size. We also emphasize that this study is not focused on identifying the specific nature of FC disturbances associated with psychosis but instead concentrates on how parcellation type affects FC differences in the same patients. In this context, medication exposure was constant across our main contrast of interest (parcellation type), meaning that it cannot explain the differences that we focus on here. The same reasoning applies to the clinical heterogeneity of the patient sample, which composed people diagnosed with both affective and nonaffective psychoses. Future work could use individualized parcellations to delineate FC differences more precisely between distinct patient subgroups.

We have focused here on how the use of individualized versus group-based parcellations affects group differences in FC. A separate question concerns whether parcellation type also affects correlations with behavioral or clinical variables. Several studies have shown that individualized parcellations yield FC estimates that are marginally more correlated with various forms of behavior, including psychopathological ratings (Bijsterbosch et al., 2018; Kong et al., 2019, 2021). A useful direction for future work could involve investigating whether individualized parcellation improves prediction of clinically meaningful outcomes.

Our findings indicate that traditional reliance on group-based parcellations may inflate case-control differences in FC at a fine-grained level. The use of individualized parcellations can yield a more conservative understanding of brain network disruptions in psychotic and possibly other disorders. However, it does not greatly impact case-control differences in network-level analyses.

All data for this study were collected as part of the HCP-EP study, which is an open-access collection aiming to generate high-quality imaging data in early psychosis patients and healthy controls (HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021). This study includes high-resolution structural and functional MRI data from 121 patients with early psychosis (74 males) and 57 healthy individuals (37 males). Demographic information is provided in Table 1. Data collection by HCP-EP has been approved by the Partners Healthcare Human Research Committee/IRB and complies with the regulations set forth by the Declaration of Helsinki (Lewandowski et al., 2020).

The patient group was composed of outpatients with psychosis, meeting criteria for affective or nonaffective psychosis, according to the DSM-5, who were within the first 5 years of onset of symptoms. Patients were recruited by four clinical programs: Beth Israel Deaconess Medical Center (BMH)–Massachusetts Mental Health Center (BIDMC-MMHC), Prevention of and Recovery from Early Psychosis (PREP) Program; Indiana University Psychotic Disorders Program, Prevention and Recovery for Early Psychosis (PARC); McLean Hospital, McLean On Track; and Massachusetts General Hospital, First Episode and Early Psychosis Program (FEPP) (HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021). Imaging took place in three of these sites.

The control group included volunteers who did not present with anxiety disorders and/or psychotic disorders, had no first-degree relative with schizophrenia spectrum disorder, were not taking psychiatric medication at the time of the study, and had never been hospitalized for psychiatric reasons. All participants were aged between 16 and 35 years old (mean = 23, SD = ±3.9) at the time of the study (Table 1). A total of 11 subjects were excluded because of poor data quality, as detailed below, leaving a final sample of 55 (36 male) controls and 112 (67 male) patients.

The participants recruited from four locations were scanned at three sites: BMH; Indiana University; and McLean Hospital, using Siemens MAGNETOM Prisma 3T scanners. The acquisition parameters between the three sites were harmonized and followed the widely used HCP protocol (Demro et al., 2021; HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021). The project collected whole-brain T1-weighted MRI (T1w), T2-weighted MRI (T2w), diffusion MRI, spin echo field maps with anterior to posterior and posterior to anterior phase encoding directions, and consisted of four resting-state functional MRI (rsfMRI) sessions. The current study uses the T1w and T2w images, the spin echo field maps, and the first two runs of the rsfMRI scans. A 32-channel head coil was used at BMH and Indiana University. A 64-channel head and neck coil, with neck channels turned off, was used at McLean Hospital. Real-time image reconstruction and processing was performed for quality control, and scans with detectable problems were repeated (HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021).

Acquisition parameters followed HCP standards. T1w images were obtained using a magnetization-prepared rapid gradient-echo (MPRAGE), with 0.8 mm isotropic spatial resolution echo time (TE) = 2.22 ms, repetition time (TR) = 2,400 ms, and field of view (FoV) = 256 mm. T2w images were acquired following a 3D-SPACE sequence, with 0.8 mm isotropic spatial resolution, TE = 563 ms, TR = 33,200 ms, and FoV = 256 mm (HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021).

The present study mainly utilized the first rsfMRI run (with anterior to posterior phase encoding). The second run (with posterior to anterior phase encoding) was used to validate the parcellation with out-of-sample analysis of within-parcel homogeneity. Scans were acquired for a length of 6.5 min, resulting in a total of 420 volumes; the first 10 volumes were removed prior to the dataset release. Images have an isotropic spatial resolution of 2 mm, TE = 37 ms, TR = 800 ms, and FoV = 208 mm. A multiband acceleration factor of 8 was used to improve spatial and temporal resolution (HCP Early Psychosis 1.1 Data Release: Reference Manual, 2021).

All analyses were done on the MASSIVE high-performance computing facility (Goscinski et al., 2014). Raw structural and functional images were first visually inspected for large artifacts and distortions. Images were then put through an automated quality-control pipeline (MRIQC; Esteban et al., 2017), which computes 15 image quality metrics for each scan with the purposes of identifying outliers warranting closer inspection. At this stage, three subjects were excluded for missing or unusable structural images.

Head motion is a major source of noise in fMRI signals. Its effects remain present even after volume realignment and can introduce systematic bias in case-control studies when not strictly controlled (Parkes et al., 2018; Power et al., 2012). Head motion during the fMRI scan was estimated using framewise displacement (FD), which is a summary measure of the movement of the head from one volume to the next (Parkes et al., 2018). For each scan, FD was calculated according to the method described by Jenkinson et al. (2002) and the resulting trace was band-pass filtered and downsampled to account for the high sampling rate of the multiband fMRI acquisition (Power et al., 2019). Subjects were excluded if they met at least one of the following stringent exclusion criteria: scans had a mean filtered FD greater than 0.25 mm; more than 20% of frames were displaced by more than 0.2 mm; or any frame was displaced by more than 5 mm. These criteria have previously been shown to effectively mitigate motion-related contamination in fMRI connectivity analyses (Parkes et al., 2018). In total, 11 subjects (2 controls) were excluded for excessive head movement in the scanner.

T1w images were processed using FreeSurfer version 6.0.1 (Dale et al., 1999) to generate cortical surface models for each participant. Surfaces were visually examined for inaccuracies and distortions. The fMRI data were processed according to the minimal preprocessing pipeline for HCP data (Glasser et al., 2013). The pipeline adapts steps from FMRIB Software Library (FSL) and FreeSurfer to account for greater spatial and temporal resolution and HCP-like distortions resulting from acquisition choices such as multiband acceleration (Glasser et al., 2013). Briefly, images were skull stripped by the brain extraction tool (Smith, 2002) of FSL, which removes non-brain matter from the image. Skull-stripped T1w, T2w, and fMRI were aligned using FMRIB’s linear image registration tool (FLIRT; Jenkinson et al., 2002; Jenkinson & Smith, 2001). Spin echo EPI field maps with opposite phase-encoding directions were used to estimate spatial distortion caused by magnetic field inhomogeneities, with corrections applied using FSL’s topup tool (Andersson et al., 2003) and FLIRT. This process was fine-tuned and optimized using FreeSurfer’s BBRegister (Greve & Fischl, 2009). Furthermore, bias field correction was performed on structural images to remove gradients of voxel intensity differences, following the HCP pipeline (Glasser et al., 2013). The fMRI volumes were realigned to the first volume for each participant using FLIRT. The fMRI data were then coregistered to their structural image, and the structural image was nonlinearly normalized into standard Montreal Neurological Institute (MNI) ICBM152 space (Grabner et al., 2006) using FLIRT and FMRIB’s nonlinear image registration tool (FNIRT) (Andersson et al., 2010). The resulting transform was then applied to the functional data.

The functional data were denoised using independent component analysis–based X-noiseifier (FIX), which decomposes the data into spatially independent components and uses machine learning to label each resulting component as either signal or noise (Griffanti et al., 2014; Salimi-Khorshidi et al., 2014). The preprocessed fMRI time series were then regressed against the estimated noise component signals and the residuals were retained for further analysis. Component decomposition was performed using multivariate exploratory linear optimized decomposition into independent components (MELODIC) (Griffanti et al., 2014; Salimi-Khorshidi et al., 2014). HCP’s training set—HCP_hp2000, which includes pretrained weights to classify independent components—was used as the training set for the algorithm. A temporal high-pass filter (2000s full width half maximum) was applied to remove low-frequency signal drifts, as recommended by the HCP preprocessing guideline (Glasser et al., 2013). Following HCP’s guidelines (Glasser et al., 2013), a lenient threshold component labelling in FIX was used (th = 10), regressing out the noise components while controlling for the signal components. The accuracy of the labels was manually verified. The analyses were repeated after applying global signal regression (GSR), which removes widespread signal fluctuations associated with respiratory variations (Aquino et al., 2020; Power et al., 2017) (see the Supporting Information).

The processed images in MNI volume space were resampled to each individual’s cortical surface, as generated by FreeSurfer, and then registered to the fsaverage5 template using a surface-based registration algorithm (Dale et al., 1999; Fischl, 2012). Fsaverage5 is a standard template generated by FreeSurfer; the resulting surface mesh comprises a total of 20,484 vertices.

We used group parcellations provided by Schaefer et al. (2018) as the basis for our analysis, as this parcellation is widely used and has shown superior functional homogeneity compared with other leading approaches (Schaefer et al., 2018). Our study focused on the 100-region parcellation, organized into seven networks (s100), but we repeated the analyses using the 200-region variant to check the robustness of the results (see the Supporting Information). Regions were screened for low BOLD signal intensity, with a method adapted from Brown et al. (2019). Specifically, we found the elbow of the BOLD signal distribution, given by the largest decrease in pairwise differences of the mean BOLD signal of each region. This was used as a cutoff for signal dropout, and regions with lower signal than the cutoff were considered to have signal dropout. Regions that were found to have signal dropout in over 5% of subjects were excluded before analysis. For the s100 atlas, 15 regions were excluded; for the s200 atlas, 16 regions were excluded from further analysis.

To derive individually tailored parcellations, we used the group prior individualized parcellation (GPIP) model (Chong et al., 2017), which relies on a Bayesian formulation with two priors: one based on group FC and one that drives individualized parcel boundaries. The former uses a group sparsity constraint to represent FC between parcels, which allows the model to maintain comparability between subjects. The latter uses a Markov random field in the form of a Potts model to label the set of parcels and maximize the FC homogeneity within each parcel based on individual data. This model allows for comparability between subjects, as it maintains the same areas and labels for every individual while capturing the variability in the shape and size of each parcel to best estimate each subject’s functional regions. Individualized parcel borders were optimized across 20 iterations, starting with the group-based Schaefer atlas and iteratively alternating between updating individual borders and the group FC prior. Further details are provided in Chong et al. (2017). The algorithm was applied to patients and controls together.

For both group-based and individualized parcellations, mean time series were extracted for each region in the s100 and s200 atlases using each individual’s spatially normalized and denoised functional data. Product-moment correlations were then estimated for every pair of regional time series to generate FC matrices. We consider only cortical areas here as, to our knowledge, methods for developing individualized parcellations for subcortical and cerebellar regions have not yet been developed.

To investigate the differences in parcels between the two parcellation approaches, we computed how many vertices were reassigned to a different parcel after applying GPIP. We then compared the number of vertices relabeled between patients and controls at ROI and whole-brain levels. All between-group statistical analyses were evaluated using permutation-based inference, with 5,000 permutations, unless otherwise indicated. Statistically significant effects for ROI-level analysis were identified using an FDR-corrected (Benjamini & Hochberg, 1995) threshold of p_FDR < 0.05, two-tailed.

To measure intrasubject reliability, we also computed homogeneity scores in the first run and compared these results between parcellation approaches, using the intraclass correlation coefficient.

We assessed how parcellation type influences FC differences between patients with psychosis and healthy controls in three ways. First, we examined the distribution of unthresholded t statistics obtained at each edge using a general linear model to quantify mean differences between patient and control groups. This and all subsequent analyses are controlled for age, sex, test site, and mean FD. The contrast was specified such that a larger t statistic indicated lower FC in patients, compared with controls. To compare the similarity of the symmetric t matrices, we vectorized their upper triangles and computed their Spearman correlation. The effect of parcellation type was evaluated using a shift function test on these distributions (Rousselet et al., 2017) to evaluate whether differences between parcellations were restricted to specific quantiles of the t statistic distributions (rather than just comparing the means of these distributions). The shift function computes the difference in value of the nine deciles of the distributions. For inference, it computes the 95% confidence interval associated with each decile difference, based on a bootstrap estimation of the standard error of each decile, controlling for multiple comparisons, via the Hochberg’s method. This analysis thus allowed us to determine whether parcellation type preferentially affected results for edges that showed small, moderate, or large case-control differences.

Second, we compared thresholded results obtained with the network-based statistic (Zalesky, Fornito, & Bullmore, 2010). NBS is an adaptation of cluster-based statistics for network data. A primary threshold of p = 0.05, uncorrected, was applied to the matrix of t statistics obtained using the general linear model described above. Results were repeated with a threshold p = 0.01 and 0.001. The sizes of the connected components of the resulting network (in terms of number of edges) were then estimated. In this context, the connected components represent sets of nodes through which a path can be found via supra-threshold edges. The group labels (patients and controls) were permuted 5,000 times and the previous steps were repeated. At each step, the size of the largest connected component was retained, resulting in an empirical distribution of maximal component sizes under the null hypothesis. The fraction of null values that exceeded the observed component sizes corresponds to a family-wise corrected p value for each component. The test was repeated with different family-wise error rate corrected p values = 0.05, 0.01, and 0.001, all resulting in the same connected component. By performing inference at the level of connected components rather than individual edges, the NBS results in greater statistical power than do traditional mass univariate thresholding methods (Zalesky, Fornito, & Bullmore, 2010). This analysis was repeated for each parcellation type (i.e., group-based and individualized) and scale (i.e., s100 and s200). Differences between significant component sizes observed using the two parcellation methods were then estimated and evaluated with respect to the differences between null component sizes estimated for the two approaches.

We calculated changes in parcel size between parcellation approaches for patients and controls. We compared parcel size difference with a two-sample t test between patients and controls. To understand how parcel size impacted FC measures, we calculated the Spearman’s rho correlation between the t values for parcel size and mean dysconnectivity per parcel and degree of dysconnectivity. The p values were calculated with a spin permutation test, with 5,000 permutations (Alexander-Bloch et al., 2018).

Finally, we examined how parcellation type affects case-control differences at the level of seven canonical networks. We considered the control network; the default mode network; the dorsal attention network; the limbic network; the salience/ventral attention network; the somatomotor network; and the visual network using the seven Yeo network assignments associated with the s100 and s200 atlases (Yeo et al., 2011). Specifically, we quantified the number of edges in the significant NBS component that fell within and between these seven networks. We examined both raw edge counts and counts normalized for the size of each network/network pair and quantified the correlation between the resulting network-level matrices obtained for each parcellation type.

Code used for analysis and image generation can be found online at https://github.com/NSBLab/individualised_parc_psychosis, and code for individualized parcellation can be acquired online at https://neuroimageusc.github.io/GPIP.

Supporting information for this article is available at https://doi.org/10.1162/netn_a_00329.

Priscila Thalenberg Levi: Data curation; Formal analysis; Investigation; Methodology; Software; Visualization; Writing – original draft; Writing – review & editing. Sidhant Chopra: Data curation; Methodology; Software; Supervision; Writing – review & editing. James C. Pang: Software; Supervision; Writing – review & editing. Alexander Holmes: Data curation; Writing – review & editing. Mehul Gajwani: Software. Tyler A. Sassenberg: Software; Writing – review & editing. Colin G. DeYoung: Software; Writing – review & editing. Alex Fornito: Conceptualization; Funding acquisition; Project administration; Resources; Supervision; Writing – review & editing.

Alex Fornito, Sylvia and Charles Viertel Charitable Foundation (https://dx.doi.org/10.13039/100008717), Award ID: 1197431. Alex Fornito, National Health and Medical Research Council, Award ID: 1146292.

Functional coupling (FC):

Coactivation of two brain regions, given by the correlation of the BOLD time series of a pair of regions.

Parcellation:

A method that divides the brain into parcels, or nodes.

Individualized parcellation:

A parcellation approach that adapts each individual’s ROI borders to account for individual differences in brain organization.

Group-based parcellation:

A one-size-fits-all approach to parcellation that is usually based on a group average.

Brain atlas:

A map of the brain that guides parcellation and can be fit into individuals’ brains based on structural landmarks, such as sulci and gyri.

Schizophrenia:

A neurodevelopmental disorder characterized by positive symptoms, such as psychosis, and negative symptoms, such as cognitive impairment.

BOLD signal:

Blood oxygen level–dependent signal, used to infer neural activity, measured by fMRI.