Brain network reconstruction
Structural brain networks were reconstructed from the structural and diffusion magnetic resonance imaging (MRI) data provided by the Human Connectome Project (HCP) database26† Individual nodes and links in the networks were chosen to represent brain regions and anatomical connections between them. The HCP scanning protocol was performed in compliance with the approval of the local institutional review board at Washington University in St. Louis. From the database, the preprocessed structural and diffusion MRI data of 100 unrelated subjects from HCP 1200 Subject Release was used to reconstruct the structural brain networks of 100 healthy young adults. The full details of the preprocessing pipeline have been described in the previous study27† Reconstruction of the structural brain networks was carried out through the following procedure by using the MRtrix3 package28: (1) Generate a tissue-segmented image from structural MRI data for Anatomically-Constrained Tractography (ACT)29† (2) match the parcellated brain image to the Destrieux atlas30, which consists of 164 cortical and subcortical regions; (3) calculate the fiber orientation density (FOD) with diffusion MRI data by multi-shell, multi-tissue constrained spherical deconvolution (MSMT-CSD)31† (4) trace white-matter neural fibers connecting each pair of brain regions by using probabilistic tractography algorithm32† (5) construct the adjacency matrices for each structural brain network by enumerating the number of neural fibers between brain regions. The overall streamline count number was set to 10 million, which is the default value of the MRtrix3.
We binarized each brain network to the presence or absence of inter-regional connectivity, to compare brain networks with various types of binarized real-world complex networks and to solely focus on investigating the contribution of connectivity pattern in determining community architecture. If connectivity exists between region i and j, the value of the corresponding network element is set to 1; otherwise, the value is set to 0. We set each brain network to have a link density of 10%, and further reconstructed a group-averaged brain network (the consensus brain network) by selecting all connections that are present in at least 75% of the 100 subjects. Through this process, each of the reconstructed brain network and the group-averaged brain network are realized in the form of an undirected and unweighted network. The community architecture of brain networks was also determined at different link densities (Fig. S2).
Network collection
For comparison with the structural brain networks, we collected 157 real-world complex networks from 4 disparate fields, ranging from transportation networks (31 networks), technological networks (36 networks), biological networks (30 networks), and social networks (60 networks ). These networks were collected from the Stuart et al†33, and were also used in the form of an undirected and unweighted network for further analysis. Details of each network are summarized in Table S1.
We generated 100 randomized and 100 latticized (ie grid-structured) null networks derived from the group-averaged brain network of 100 subjects, which we described above. Each type of null network was generated using the functions randmio_und_connected or latmio_und_connected in the Brain Connectivity Toolbox34† These functions randomly permute the edges of a network or latticizes a network while preserving the number of nodes, links, and its degree distribution.
We also prepared synthetic modular networks with different community segregation scores. Model networks with equally-sized modules were created by manipulating the probability of within- and between-community edge placement (Code is available at http://github.com/macshine/integration/guimera_model.m)35† The probability of within-community edge placement was varied between 0.7 and 1 with 0.02 intervals, and the probability of between-community edge placement was varied between 0.01 and 0.35 with 0.02 intervals. To investigate a general relationship between the degree of community segregation and overlap, the values of the two probabilities were obtained across the parameter space for generating each synthetic modular network. We used 7 different modular configurations in this study; (the number of nodes N, the number of modules M) = (100, 5), (200, 5), (300, 5), (120, 4), (120, 6), (120, 8), and (120, 10).
Community architecture in complex networks
To characterize the community architecture of human brain networks and compare it with that of other real-world complex networks, we quantified two distinct properties of community architecture with the following metrics: (1) the segregation score, which quantifies the degree of segregation between communities , and (2) the overlap score, which quantifies the degree of overlap between communities. To measure the segregation score, we partitioned a complex network into non-overlapping communities and computed the modularity index using the Louvain community detection algorithm (gamma = 1)5as implemented in the Brain Connectivity Toolbox34† The resultant modularity index was utilized as a segregation score of the network. To measure the overlap score, we partitioned a complex network into overlapping communities using the link community algorithm6as implemented in the Brain Connectivity Toolbox34† Then, we quantified the overlap score by calculating the average number of different communities to which each node belongs. In the process of calculating the overlap score, we ensured that each node in a network belongs to at least one community. We also examined the number of overlapping communities of each network to check whether it maintains its segregated community structure. If the number of overlapping communities of a network is equal to or higher than the number of predefined non-overlapping modules, we determined that the network maintains its segregated community structure. In addition, when comparing network groups with different attributes, such as brain networks and real-world complex networks, the overlap score of each network has been normalized by the average overlap score of its 10 random null networks for group-wise comparisons.
Topological reinforcement model
To analyze the functional benefits of human brain network community architecture, we employed the topological reinforcement (TR) model, which mimics activity-dependent adaptation and synaptic plasticity during a neurodevelopmental process. By promoting connections between nodes that share more common neighbors than others, this model reliably evolves initial networks with a random structure into more modularized networks and supports the role of TR as a contributor to the emergence of modular brain networks. A more detailed description of the TR model is given in the previous study7†
We prepared 50 Erdős–Rényi (ER) random networks without self-connections of size N = 164 nodes and average node degree (lambda) = 16 (equivalent to a density of 10%). At each rewiring step, we randomly selected N/2 nodes that are neither disconnected nor fully connected to all other nodes. Then we added simultaneously one link on each selected node with a non-neighbor that has the highest topological overlap score. The topological overlap score represents the neighborhood similarity of a pair of nodes by counting the number of their common neighbors7† After inserting links, the same number of links were randomly eliminated to maintain the original link density. We set the number of evolution steps to 60 so that sufficient link rewiring can take place for modularity emergence as long as the network is connected. In each rewiring step, we collected resultant networks to analyze the functional benefits of human brain network community architecture.
Graph theoretical analysis
All graph-theoretical analyzes were performed using the Brain Connectivity Toolbox34 (http://sites.google.com/site/bctnet/).
linear threshold model
We employed the linear threshold model to compare communication efficiency in view of network diffusion. In this model, every node of a network has a common threshold parameter within [0, 1] and a certain percentage of nodes are selected as ‘seed nodes’ for initial stimuli. An activated node can have 1/node_degree of influence on each of its neighbors, and the influence a given node receives is calculated by summing the influences from its activated neighbors. If the influence received by a given node exceeds a predefined threshold parameter, it becomes an activated node that can influence its neighbors. This chain reaction is finished until there is no more changes in the number of activated nodes in a network, which is referred to as a steady state (code is available at https://github.com/hhchen1105/networkx_addon/information_propagation/linear_threshold. py)36†
We varied the proportion of seed nodes from 5 to 50% and set a threshold parameter to 0.5. We randomly selected seed nodes and then simulated a diffusion process until the number of activated nodes in a network no longer changes. Communication efficiency in view of network diffusion was quantified by the ratio of activated nodes at a steady state.
Dynamic flexibility of functional brain networks
To examine the dynamic flexibility of functional brain networks of the 100 subjects, we used each subject’s resting-state functional MRI (rfMRI) time-series data (sampling rate = 0.72 s.) that are publicly released through the central HCP Connectome DB database. First, the time-varying inter-regional functional connectivity was estimated by correlation-based sliding-window analysis using the time-series data with a different number of independent component analysis (ICA) dimensionalities (N = 50 or 100). This results in 50- or 100-node time-series of 4,800 time points for each subject. We used a window length of 30 s. and each window was shifted 0.72 s. across the whole scan, resulting in a total of 4759-time windows (W). We generated a 3D matrix of correlation coefficients of size N × N × W, composed of 2D N × N matrices for every time window W. Using the dynamic community detection algorithm5, we identified network modules in each time window and tracked their evolution over time. By estimating the probability that a brain region changes its allegiance to modules between any two consecutive time windows, we computed the dynamic flexibility of each brain region. To investigate a structure–function relationship, we correlated the mean dynamic flexibility of the functional brain networks and the overlap score of the 100 structural brain networks (Fig. 3B and S3).
Behavioral measures of cognitive flexibility
We used three behavioral measures of cognitive tasks that enabled the examination of the three aspects of cognitive flexibility: cognitive inhibition, executive function, and processing speed. These cognitive tasks include (1) the pattern completion task (processing speed), (2) the card sorting task (executive function), and (3) the flanker task (cognitive inhibition)37† The behavioral measures of each task are recorded for each HCP participant. These measures are part of the HCP NIH toolbox, which is a multidimensional set of behavioral measures evaluating diverse aspects of cognitive functions26.37† To investigate a relationship between the community architecture and cognitive flexibility, we correlated the performance score and the overlap score of the 100 structural brain networks in each cognitive task.
Statistical analysis
We performed two tailed t-tests to determine whether the characteristics of the given network groups shows a significant difference in the perspective of community architecture.
