Construction and Analysis of a New Resting-State Whole-Brain Network Model

Background: Mathematical modeling and computer simulation are important methods for understanding complex neural systems. The whole-brain network model can help people understand the neurophysiological mechanisms of brain cognition and functional diseases of the brain. Methods: In this study, we constructed a resting-state whole-brain network model (WBNM) by using the Wendling neural mass model as the node and a real structural connectivity matrix as the edge of the network. By analyzing the correlation between the simulated functional connectivity matrix in the resting state and the empirical functional connectivity matrix, an optimal global coupling coefficient was obtained. Then, the waveforms and spectra of simulated EEG signals and four commonly used measures from graph theory and small-world network properties of simulated brain networks under different thresholds were analyzed. Results: The results showed that the correlation coefficient of the functional connectivity matrix of the simulated WBNM and empirical brain networks could reach a maximum value of 0.676 when the global coupling coefficient was set to 20.3. The simulated EEG signals showed rich waveform and frequency-band characteristics. The commonly used graph-theoretical measures and small-world properties of the constructed WBNM were similar to those of empirical brain networks. When the threshold was set to 0.22, the maximum correlation between the simulated WBNM and empirical brain networks was 0.709. Conclusions: The constructed resting-state WBNM is similar to a real brain network to a certain extent and can be used to study the neurophysiological mechanisms of complex brain networks.


Introduction
Mathematical modeling and computer simulation are the main methods for understanding neural mechanisms, especially in the field of computational neural systems.A whole-brain network model (WBNM) is a network model that is composed of coupled brain regions; these are usually represented by a neural mass model (NMM), which represents the average activity of a brain region on a mesoscopic scale.Their interconnection is usually based on biologically informed data on the brain's structural connectivity [1].In this study, our aim was to use mathematical modeling methods to construct a WBNM with neurophysiological significance.It is similar to real human brains in terms of brain rhythm, structural connectivity, functional connectivity, and brain network characteristics, and the complexity of the model is also considered.This enables researchers to study the neurophysiological mechanisms of functional diseases of the brain and brain modulations based on this model.The human brain is a multiscale complex system composed of many dynamic units that interact to generate various behavioral functions [2].From a small-scale perspective, the entire brain is composed of numerous types of neurons that are connected to one another, and electrical pulses propagate through this vast distributed network.From a large-scale perspective, the entire brain can be seen as a system network with complex dynamic characteristics.At a mesoscopic scale, that is, at the level of the neural masses or brain regions, the collective activity often exhibits rhythmic patterns [3].Neural oscillations promote time coordination in distant regions, leading to the phenomenon of functional interactions, which play a crucial role in interregional communication and information transmission between distributed brain regions [4][5][6].Therefore, a WBNM can help researchers study the functional connectivity of EEG signals and the potential mechanisms of the brain from both mesoscopic and macroscopic perspectives.
In order to construct a WBNM, we need to select a model as a node.An NMM consists of a series of differential equations and provides a new way to study the brain's rhythmic activity and the neural mechanisms of functional diseases of the brain [7].A key strength of NMMs over alternative statistical mathematical frameworks lies in their foundation on anatomical and neurophysiological principles with state variables and parameters that possess well-defined physiological interpretations.This feature provides a substantial correspondence between the mathematical model and the actual EEG signal, and it makes the use of NMMs more reliable when analyzing actual EEG signals [8][9][10][11][12].Freeman [13] initiated the "neural mass" concept by aggregating the activities of neuronal groups into a lumped parameter framework, while its mathematical underpinnings were elaborated by Wilson and Cowan [14].Subsequently, Lopes da Silva et al. [15] introduced the alpha rhythm model (ARM), which integrated two thalamic neuronal populations to mimic thalamic alpha rhythms, underscoring the thalamus's critical function in conveying sensory information to the cerebral cortex [16].On this foundation, Jansen and Rit [17] developed the more sophisticated cortical NMM, which was termed the J&R model and was aimed at modeling the visual cortex's responses to visual stimuli.This model is distinguished by its ability to produce alpha waves using three distinct neuronal populations: pyramidal cells and both excitatory and inhibitory interneurons.Further extending this framework, Wendling et al. [18] incorporated an additional fast inhibitory interneuron class into the J&R model, categorizing the inhibitory neurons into slow and fast subgroups for enhanced model fidelity.In this study, they chose the Wendling model for the nodes of their WBNM.The following were the main reasons: Compared to simple models that only include excitation inhibition circuits, Wendling models have rich dynamic characteristics and can generate six common EEG signals, including background noise, typical alpha waves, and epileptic spike waves; the Wendling model has been widely used in the study of the cerebral cortex, especially for seizures and the transmission mechanism of epilepsy [19][20][21]; the structure of the Wendling model is relatively simple, and the computational cost of constructing a WBNM is tolerable.Ursino [22] developed a new model by modifying the structure of the fast inhibition circuits of the Wendling model, and it was able to generate gamma waves but also had a correspondingly increased.Therefore, its application is not as widespread as that of the Wendling model.
Anatomical studies have divided the brain into different regions.To simulate the coupling between brain regions, structural connectivity data are typically used to construct WBNMs.In recent years, researchers have conducted some valuable studies on WBNMs and their modeling methods [1,23] and used WBNMs to understand the working mechanisms of the brain [2,[24][25][26][27], as well as the neural mechanisms and treatments of functional diseases of the brain [28][29][30][31].Neurolib is a comprehensive Python framework designed for whole-brain modeling, and it allows researchers to simulate the mesoscopic activity of brain regions through neural mass models [1].By incorporating structural and functional brain data, including fMRI-derived BOLD signals, it enables precise calibration of models against empirical observations.The framework facilitates the in-depth analysis and optimization of brain models by providing advanced tools for exploring parameter Brain Sci.2024, 14, 240 3 of 21 spaces and fitting models to multimodal empirical data.The concept of a Virtual Brain Platform was introduced in [23].This work established a neuroinformatics platform based on real brain connectivity that can simulate a whole-brain network.Several widely used neurodynamic models are embedded in this platform, and they serve as the basis for inferring neurophysiological mechanisms at different brain scales and ultimately generating macroscopic neuroimaging signals.
Papadopoulos et al. [2] built a large-scale brain model wherein 82 brain areas were represented as Wilson-Cowan neural masses.They investigated how the stimulation of a particular brain area affected the network dynamics.The research results indicated that the network exhibited different responses to regional perturbations based on the baseline state of the system.Through advanced fMRI and diffusion spectrum imaging tractography coupled with computational modeling, Honey et al. [24] uncovered that resting-state functional connectivity in the human cerebral cortex is intricately modulated by its anatomical structure, despite the presence of connections between structurally unlinked regions.These results challenged the feasibility of inferring structural connectivity from functional connectivity, highlighting the role of indirect connections and spatial considerations in shaping functional network dynamics.This study bridges the gap between structural and functional brain networks, offering insights into their complex interplay.
Schirner et al. [25] developed a groundbreaking model that merged EEG-driven neural dynamics with connectome data; they replicated resting-state fMRI signatures and decoded complex neurophysiological relationships.This method bridged multiple observational scales-from oscillatory activity to network topology-advancing our understanding of brain function.Their findings provide a comprehensive framework for interpreting non-invasive brain measurements.Giannakakis et al. [26] leveraged the Wilson-Cowan framework and synaptic plasticity to simulate brain dynamics in real time, enabling studies of long-term changes due to injury or development.This represents a crucial advancement in understanding brain plasticity and recovery.Endo et al. [27] used the Larter-Breakspear model, and they correlated fMRI and EEG/MEG data to uncover the structural influence on resting-state brain dynamics.Their findings pinpoint synchronization patterns that bridge fast EEG microstates with slow fMRI fluctuations, enhancing our grasp of neural connectivity's role in brain function.
Hann et al. [28] constructed a WBNM with interconnected neural masses, where each neural node was based on the Arm model.The authors used this WBNM to serve as a virtual test ground for developing effective interventions for Alzheimer's disease.They explored six intervention strategies that adjusted the excitatory or inhibitory levels of neurons, aiming to identify potential methods for delaying the degradation of brain network connectivity.Their findings underscore that selectively boosting excitatory neuron activity emerges as the most efficacious intervention, markedly enhancing network integrity and, thereby, illuminating a promising pathway for therapeutic development aimed at sustaining or rejuvenating the neuronal network in AD.Taylor et al. [29] constructed an epilepsy WBNM based on Benjamin's epilepsy model.They integrated dynamic computational models with static structural connectivity data from diffusion-weighted MRI, offering a novel approach to mapping epilepsy's complex neurodynamics.Hutchings et al. [30] utilized a whole-brain epilepsy model tailored to individual structural connectivity profiles from DTI scans of 22 left TLE patients and 39 controls to explore seizure dynamics and surgical outcomes.Kunze et al. [31] applied transcranial direct-current stimulation (tDCS) across a network model of 74 cerebral areas, exploring its effects on brain dynamics and structural connectivity while being guided by the human connectome.
In this study, our objective was to develop a WBNM that more accurately reflects the functional connectivity and structural network organization of the human brain.The Wendling model was used as the node for the WBNM because of its rich spectrum and simple structure.In order to construct the edge of the WBNM, that is, the coupling strength matrix, we used actual MRI data to form an empirical structural connectivity matrix and empirical functional connectivity matrix.But only the relative coupling strength matrix could be obtained by calculating the fiber tracts between region pairs (please refer to Section 2 for specific calculation methods), so an overall multiple was needed to adjust the coupling strength, that is, a global coupling coefficient C was used to rescale the structural connectivity strength between regions.We hope that by adjusting the global coupling coefficient C, the WBNM can simulate real brain functional connections to a certain extent.Finally, we analyzed the differences between simulated brain networks and actual brain networks at different threshold levels, along with their small-world network properties, using four graph-theoretical measures: the average degree, characteristic path length, clustering coefficient, and global efficiency.We found that the simulated brain networks and actual brain networks exhibited strong similarities in their network characteristics.

Materials and Methods
Constructing a whole-brain network model (WBNM) based on structural connectivity is a feasible method for simulating whole-brain signals to elucidate the brain's functional activities and network mechanisms.In order to construct a WBNM, each brain region can be represented by an NMM, and interconnections between brain regions can be constructed through excitatory connections of pyramidal neuronal populations.The functional activities of the WBNM depend on the kinetic behavior of each NMM and the interactions between multiple brain regions.

Coupling Strength Matrix of the WBNM
To construct a resting-state WBNM, we needed to obtain a coupling strength matrix (CSM) to simulate the interconnection of the brain regions.Usually, an empirical structural connectivity (SC) matrix can be used as a relative coupling strength matrix.An empirical resting-state functional connectivity (rsFC) matrix was used to fit the simulated rsFC of the WBNM.We selected 15 participants aged between 18 and 31, including 8 females.Participants who self-reported a history of neurological, cognitive, or psychiatric disorders were excluded from the experiment.We collected diffusion-weighted MRI, T1-weighted MRI, and fMRI data.The research was performed in compliance with the Code of Ethics of the World Medical Association (Declaration of Helsinki).Written informed consent was provided by all subjects with an understanding of the study prior to data collection, and the experiment was approved by the local ethics committee in accordance with the institutional guidelines at Charité Hospital Berlin.The human data utilized in this study originated from the preliminary research conducted by Schirner et al. [25], and they were automatically processed using the pipeline for functional MRI (fMRI) and diffusion MRI (dMRI) data.The T1-weighted MRI (TR: 1900 ms, TE: 2.52 ms, TI: 900 ms, FA: 9  Germany).Detailed insights into diffusion-weighted imaging and fMRI preprocessing can be found in the study of Schirner et al. [32].In essence, the study integrated probabilistic white-matter fiber tractography with trajectory clustering within an automated framework and introduced a novel metric for connectivity.This metric weighted the connections between brain regions according to the minimum interface area of gray/white matter, moving away from traditional methods that rely on unweighted tract counts, which are susceptible to anatomical variability.Schirner highlighted that this approach reduces bias and yields a more accurate depiction of brain connectivity.In this study, structural connectivity matrices were extracted to represent the relative strengths of interactions between regions as mediated by white-matter fiber tracts.Empirical rsFC matrices were extracted from the BOLD signals of fMRI data.T1-weighted MRIs were segmented and parcellated into 68 regions according to the Desikan-Killiany brain atlas by using FreeSurfer7.2.0 [33,34].The cortical parcellation is shown in Figure 1A.Table 1 lists the 68 Desikan-Killiany atlas region labels, which are abbreviations of the brain region name, and provides mappings from the individual region to each lobe.The lobes were divided into the frontal, parietal, occipital, and temporal lobes.The Desikan-Killiany atlas comprised 34 cortical regions in each hemisphere.The left hemisphere's brain regions were sequentially numbered from 1 to 34 based on the sequence number of the 34 brain regions, and the prefix LH was added before the labels.The right hemisphere's brain regions were sequentially numbered from 35 to 68 based on the sequence number of the 34 brain regions, with the prefix RH being added before the labels.
Brain Sci.2024, 14, x FOR PEER REVIEW 5 of 23 of interactions between regions as mediated by white-matter fiber tracts.Empirical rsFC matrices were extracted from the BOLD signals of fMRI data.T1-weighted MRIs were segmented and parcellated into 68 regions according to the Desikan-Killiany brain atlas by using FreeSurfer7.2.0 [33,34].The cortical parcellation is shown in Figure 1A.Table 1 lists the 68 Desikan-Killiany atlas region labels, which are abbreviations of the brain region name, and provides mappings from the individual region to each lobe.The lobes were divided into the frontal, parietal, occipital, and temporal lobes.The Desikan-Killiany atlas comprised 34 cortical regions in each hemisphere.The left hemisphere's brain regions were sequentially numbered from 1 to 34 based on the sequence number of the 34 brain regions, and the prefix LH was added before the labels.The right hemisphere's brain regions were sequentially numbered from 35 to 68 based on the sequence number of the 34 brain regions, with the prefix RH being added before the labels.

P(t)
The schematic structure of the Wendling model     The fiber tracts between region pairs are shown in Figure 1B, and they were obtained by applying probabilistic white-matter tractography algorithms to diffusion-weighted MRI.Then, the structural connectivity matrix of each subject could be calculated, and it was defined as the number of fibers within each region pair divided by the maximum number of fibers in all region pairs.In order to eliminate individual differences, we averaged the structural connectivity matrices of all subjects to obtain a group-representative structural connectivity matrix with the elements, i,j = 1,. ..68, as shown in Figure 1C.This empirical structural matrix was used as the relative coupling strength matrix to couple the NMMs to construct the WBNM.As shown in Figure 1C, the strength of the edge connections within the same hemisphere was stronger, while the strength of edge connections between hemispheres was weaker.
Then, we computed the empirical resting-state functional connectivity matrix.From functional MRI, we obtained blood-oxygen-level-dependent (BOLD) signals based on hemodynamics.The BOLD signals had a high spatial resolution with a sampling interval of 1.94 s.The first seven images of each scanning run were discarded so that the BOLD signals could reach a steady state.The representative time series in each region were extracted by averaging the time series of the voxels therein.The parcellation mask was defined according to the Desikan-Killiany atlas.We used the phase-locking value (PLV) to calculate the empirical functional connectivity matrix for each subject.Then, a group-representative empirical rsFC matrix was obtained by averaging the rsFC matrixes of all subjects.
The PLV can quantify phase coherence between two different brain regions.Firstly, the phase θ(t) was extracted with the Hilbert transform of each brain region's time series.Then, we calculated the PLV between brain region i and brain region j with the following formula: where T s is the number of sample time points.The range of ρ ij is between 0 and 1.If the phase difference ∆θ ij (t) = θ i (t) − θ j (t) is constant over a given time window, ρ ij is equal to 1, which indicates strong phase synchronization between the two time series, whereas if the phase differences are uniformly distributed, ρ ij is approximately 0, which indicates that there is almost no phase synchronization.In a whole-brain model, setting the coupling strength matrix between different regions of the model is a crucial and thought-provoking issue.In fact, the empirical structural connectivity matrix only provides the relative coupling strength between regions.In order to determine the value of the coupling strength matrix in the simulation model and make the functional connections of the simulation signal closer to those of an actual brain system, we need a global coupling coefficient C to regulate the overall strength of the long-range coupling between the regions of the WBNM: where CS ij represents the elements of the coupling strength matrix of the WBNM, and C ij represents the elements of the empirical structural connectivity matrix.
In the study, we needed to select an appropriate parameter value for the global coupling coefficient C to achieve the highest correlation between the simulated rsFC matrix and the empirical rsFC matrix.We used the PLV to calculate the simulated rsFC matrix with different global coupling coefficients C and the empirical rsFC matrix for 68 brain regions.Then, by calculating the Pearson correlation coefficient (PCC) between the elements in the lower triangle of the simulated rsFC matrix and the corresponding positions in the empirical rsFC matrix, a similarity measure of functional connections between the simulated brain network and the empirical brain network was obtained in different situations of the global coupling coefficient C. The C value corresponding to the strongest similarity was the most suitable point, and we could obtain the most suitable coupling strength matrix for the WBNM.

The Nodes of the WBNM
The selection of nodes for a WBNM is also important.Figure 1D shows a diagram of the Wendling model's structure.The Wendling model was used to construct each brain region in the WBNM due to its rich dynamic characteristics and relatively simple structure.The Wendling model comprises four main components: excitatory pyramidal neurons, fast inhibitory interneurons, slow inhibitory interneurons, and background activity.The excitatory neurons are responsible for transmitting excitatory signals, while the fast and slow inhibitory neurons provide short-term and long-term network inhibition, respectively [18].These components are interconnected through mathematical equations that simulate the dynamic interactions between neurons and how these interactions collectively influence the brain's electrical activity patterns.Background activity comprises external inputs and internal noise, adding randomness and complexity to the model [18].The dynamic characteristics of each node also affect the entire brain network.Figure 1E illustrates the architectural framework of the WBNM.The left part of Figure 1E uses different-colored spheres to represent the distribution of the 68 brain regions, with red representing frontal, yellow representing parietal, green representing occidental, and blue representing temporary.The size of each sphere represents the strength of the node, indicating its degree of synchronization with other nodes in the network.Further details on this will be discussed later in Section 2.3.
A block diagram of the Wendling model with mutual coupling between different regions is shown in Figure 2. Within the Wendling model, the transformation functions, which are denoted as he, hi, and hg for excitatory, slow inhibitory, and fast inhibitory postsynaptic membrane potentials, respectively, are defined through mathematical formulations.These functions effectively transduce the average pulse density of action potentials into corresponding synaptic potentials, offering a quantitative framework for capturing the nuanced dynamics of glutamatergic excitation, GABA_A-mediated fast inhibition, and GABA_B-mediated slow inhibition within neural networks.S delineates the nonlinear characteristics of neuronal excitability thresholds, mapping input signals to a bounded output range to simulate the activation mechanisms of neurons and their saturation features.C represents a global coupling coefficient that can rescale the structural connectivity strength between regions.C ij represents the structural connectivity strength from region j to region i. y i out represents the output of brain region i and the input to other brain regions.The external inputs from the surrounding environment or external disturbances are encapsulated by the term P(t).These inputs are mathematically modeled as Gaussian white noise, which is characterized by a mean (µ) and variance (σ 2 ).Brain Sci.2024, 14 The linear dynamic transform functions play a crucial role in converting the average pulse density of action potentials into post-synaptic membrane potentials, which cover both excitatory and inhibitory variations: where W ∈ {A, B, G}, w ∈ {a, b, g}, W denotes the average synaptic gains, and w symbolizes the inverse of the average time constants.
The impulse response in Formula (3) corresponds to the following two first-order linear ordinary differential equations: where x(t) represents the input average pulse density of the action potential, and y(t) represents the convolution output's post-synaptic membrane potential.
The S function is a crucial component that describes the response of the average membrane potential to post-synaptic currents.This function is nonlinear, is typically represented as a sigmoid function, and is given by where 2e0 represents the maximum firing rate, i.e., the maximum output frequency of the neuron, v is the average membrane potential, v0 is the half-activation potential, meaning that when v0 = v, the output frequency is e0, and r is a slope parameter that determines the shape of the S function.In summary, node i of the WBNM can be described by the following first-order differential equations: The linear dynamic transform functions play a crucial role in converting the average pulse density of action potentials into post-synaptic membrane potentials, which cover both excitatory and inhibitory variations: where W ∈ {A, B, G}, w ∈ {a, b, g}, W denotes the average synaptic gains, and w symbolizes the inverse of the average time constants.The impulse response in Formula (3) corresponds to the following two first-order linear ordinary differential equations: .
where x(t) represents the input average pulse density of the action potential, and y(t) represents the convolution output's post-synaptic membrane potential.The S function is a crucial component that describes the response of the average membrane potential to post-synaptic currents.This function is nonlinear, is typically represented as a sigmoid function, and is given by where 2e 0 represents the maximum firing rate, i.e., the maximum output frequency of the neuron, v is the average membrane potential, v 0 is the half-activation potential, meaning that when v 0 = v, the output frequency is e 0 , and r is a slope parameter that determines the shape of the S function.In summary, node i of the WBNM can be described by the following first-order differential equations: where y i 0 to y i 9 are the state variables of region i, i = 1,2,. ..68. y j , j = 1,2,. ..68 is the output of region j, which couples with region i.
The output of region i is the post-synaptic potential of the pyramidal cell population, which is defined as The Wendling model was employed for all brain regions with identical parameter values.The standard numerical values and their physiological interpretations are listed in Table 2, as referenced in [18].
The mean and variance of the external input µ = 90, σ 2 = 30 In this study, the WBNM was solved with the fourth order Runge-Kutta algorithm in Matlab 2018b.The sampling rate was set to 1000 Hz, and the simulation length was set to 2 s.We discarded the first 1 s of the simulated signals for each brain region to eliminate the influence of the initial values.Simulations with different parameters were repeated 20 times with different initial values.

Methods of Analyzing the Brain Functional Network Based on Graph Theory
Graph theory analysis plays a pivotal role in understanding the functional connectivity networks of the brain.By conceptualizing the brain as a complex network at a macroscopic level, with nodes representing brain regions and edges representing connection strengths, Brain Sci.2024, 14, 240 10 of 21 the construction of adjacency matrices for graph-theoretical metrics allows us to precisely quantify and characterize the connectivity patterns between different regions of the brain.In this study, graph theory was used to obtain indicators for measuring the network characteristics at the node, edge, local network, and global network levels.
Before calculating the network metrics in graph theory, we first constructed an adjacency matrix using the absolute thresholding method.By setting a threshold, all values in both the simulated rsFC matrix and the empirical rsFC matrix were compared with this threshold.Elements greater than or equal to the threshold were considered as edges in the network (connections exist), while those less than the threshold were regarded as non-connected.The range of threshold values explored was from 0 to 0.5, with a step size of 0.01, resulting in a total of 51 thresholds.The network's graph theory properties were calculated at each threshold value.
The most commonly used measures in graph theory for analyzing complex brain networks are the node strength, node degree, characteristic path length, clustering coefficient, and global efficiency.These metrics comprehensively reflect a network's topological characteristics from various perspectives, including but not limited to those of information transmission efficiency, clustering properties, and key brain regions.They are used to compare the structural and functional similarities and differences between a WBNM and empirical brain networks.These measures were calculated by using the Brain Connectivity Toolbox in this study.In the following, we provide an introduction to these five metrics.
(1) Node strength: Strength represents a fundamental metric for quantifying the connectivity and importance of nodes within weighted networks.The calculation of the strength of node i is defined as where C ′ ij represents the weight of the edge connection between node i and node j, and N denotes the set of nodes that are adjacent to node i. Larger node strengths point to central influential positions with critical roles in maintaining network connectivity and function.In contrast, smaller node strengths signify peripheral positions with less influence on the network's global structure.
(2) Node degree: This degree measures the number of edges that a node shares with other nodes.The calculation of the degree of node i is defined as where a ij ∈ {0, 1} represents whether there is an edge connection between node i and node j.The degree is the most basic network metric, and most other metrics are related to the node degree.The smaller the value of D i , the smaller the role of node i in network information transmission.The larger the value of D i , the greater the role of node i in network information transmission.
The average degree of a network is a key metric in network analysis that quantifies the typical level of connectivity of nodes within a network.It is calculated as the mean value of the degrees of all nodes in the network, providing a simple yet insightful measure of the network's overall connectivity and density.The calculation of the average degree of a network is defined as where N represents the total number of nodes in the network.D = 1 indicates that there are edge connections between any two nodes in the network, while D = 0 indicates that all nodes in the network are isolated nodes without edge connections.
(3) Characteristic path length: The characteristic path length is generally defined as the average of the shortest path lengths between all possible pairs of nodes in a network.The calculation of the characteristic path length of a network is defined as where N represents the total number of nodes in the network, and d ij represents the shortest path length between node i and node j.The maximum number of edges that may exist in an undirected network is N(N − 1)/2.The value of the characteristic path length, L, is inversely related to the speed of information transmission within a network.Specifically, a smaller value of L signifies that the network facilitates quicker information exchange.Conversely, a larger value of L suggests a slower information transmission speed.
(4) Clustering coefficient: The clustering coefficient for a node i describes the ratio of the actual number of edges e i between the neighbors of i to the maximum possible number of edges between them.For node i, the local clustering coefficient CC i is defined as where e i represents the actual number of edges between node i and its neighboring nodes, and k i is the degree of node i, which indicates the number of neighbors.The range of variation in CC i is between 0 and 1.The smaller the value of CC i , the more dispersed the nodes are and the sparser the network is.The larger the value of CC i , the more clustered the nodes are and the closer the network is.
The global clustering coefficient can be understood as the average of the local clustering coefficients of all nodes in the network, and it provides a measure of the network's overall tendency to cluster.Thus, the global clustering coefficient CC can also be calculated through the following formula: where N represents the total number of nodes in the network, and CC i is the local clustering coefficient of node i.When CC = 1, this indicates that any two nodes in the network are adjacent nodes.When CC = 0, this indicates that there are no adjacent nodes between all nodes in the network.
(5) Global efficiency: The global efficiency is a metric used in network analysis to assess the overall efficiency of information transfer within a network.It is based on the shortest path lengths between pairs of nodes in the network, and it reflects the network's ability to process and transmit information on a global scale.The formula for calculating the global efficiency is defined as the average efficiency between all pairs of nodes in the network, and it is mathematically expressed as where N represents the total number of nodes in the network, and d ij represents the shortest path length between node i and node j.The range of variation in E is between 0 and 1.
The smaller the value of E, the lower the efficiency of network information exchange.The higher the value of E, the higher the efficiency of network information exchange.The small-world network index can be used measure whether a network is a smallworld network that lies somewhere between a regular network and a random network.Small-world networks have the characteristics of high clustering coefficients and low characteristic path lengths.The calculation of the small-world network index of a network is defined as follows: where γ represents the clustering coefficient ratio, and λ represents the characteristic path length ratio.CC real and L real represent the clustering coefficients and characteristic path lengths of the network to be analyzed, respectively.CC random and L random represent the clustering coefficients and characteristic path lengths of the random network corresponding to the network to be analyzed, respectively.

The Fitting of the Simulated rsFC Matrix and the Empirical rsFC Matrix
In order to make the functional connections of the WBNM closer to the empirical connections of an actual brain, the Pearson correlation coefficient (PCC) between the simulated rsFC matrix of the WBNM and the empirical rsFC matrix with different global coupling coefficients C and with a step size of 0.1 was calculated and is shown in Figure 3A.The scatter points were marked in blue, and a red solid line was used to fit the scatter plot to represent the trend of correlation changes with different global coupling coefficients, for which the Savitzky-Golay filtering method was adopted.The larger the correlation value, the better the simulated brain functional connection fit the empirical brain functional connection.The results showed that when the global coupling coefficient increased, the correlation showed a trend of first increasing and then decreasing.The correlation increased when C was greater than 70 and then decreased again when C was greater than 75.The correlation reached a maximum value of 0.676 with a significance level of p < 0.001 when the global coupling coefficient was C = 20.3, as marked with a red circle.At this point, the simulated rsFC matrix of the constructed WBNM and the empirical rsFC matrix had the highest similarity, and the WBNM was closest to an actual brain network.Therefore, all subsequent analyses of the WBNM were based on the global coupling coefficient C = 20.3, which is called the optimal global coupling coefficient fitting point.To observe the linear correlation between the simulated and empirical rsFC matrixes at the optimal global coupling coefficient fitting point in detail, a linear regression was calculated, and the results are shown in Figure 3B.The horizontal axis represents the values of the lower triangular elements of the simulated rsFC matrix, while the vertical axis To observe the linear correlation between the simulated and empirical rsFC matrixes at the optimal global coupling coefficient fitting point in detail, a linear regression was calculated, and the results are shown in Figure 3B.The horizontal axis represents the values of the lower triangular elements of the simulated rsFC matrix, while the vertical axis represents the values of the lower triangular elements of the empirical rsFC matrix at the same position.The red line represents the linear regression fitting results of the scatter plot.The PCC between the simulated rsFC matrix and the empirical rsFC matrix at the optimal global coupling coefficient fitting point reached 0.676 and p < 0.001.
The details of the simulated rsFC matrix and the empirical rsFC matrix are shown in Figure 3C.The values of the simulated rsFC matrix and the empirical rsFC matrix were similar in both the intra-hemispheric and inter-hemispheric regions, while in most brain regions, the values of the simulated rsFC matrix were higher than those of the empirical rsFC matrix.

The Characteristics of the Output Signals of the WBNM
In the WBNM, each node represented a brain region.In order to analyze the dynamics of the simulated output signals of each brain region more clearly, we first calculated the node strength of each brain region.Figure 4A shows the node strengths of each brain node, which were basically symmetrical in the left and right hemispheres.The highest node strength corresponded to the 27th brain region (LH.SF), which is marked with red color, and the lowest node strength corresponded to the 31th brain region (LH.FP), which is marked with bright blue.The 51th brain region (RH.POP) showed a moderate node strength, which is marked with orange.The stronger the node strength, the greater the importance of information transmission in the brain network; such nodes are commonly referred to as hub nodes.The smaller the node strength, the less important it is for information transmission in the brain network.
Displaying the simulated time series and corresponding normalized spectra of all 68 brain regions would make the figure very chaotic, so we selected three representative brain regions, as shown in Figure 4B.The simulated output signal of brain region 27 (LH.SF), which had the highest node strength, showed high-amplitude, high-frequency oscillations, with a peak frequency that fell on the alpha band at about 8 Hz; this is marked in red.The output of brain region 31 (LH.FP), which had the lowest node strength, showed characteristics of EEG background noise with a low amplitude and a broadband spectrum; this is marked in bright blue.It can be seen that the coupling strength between the region 31 and other regions was very low in both the simulated rsFC matrix and empirical rsFC matrix in Figure 3C.Therefore, region 31 showed the characteristics of EEG background noise.The output of brain region 51 (RH.POP), which had moderate node strength, is marked in orange.It also showed alpha-wave oscillation, with a lower amplitude and slightly wider spectrum than those of region 27; the dominant frequency was also about 8 Hz.
The baseline potential and dominant frequency of the simulated output signals for all regions are shown in Figure 4C-F.In Figure 4C, the bar plot shows the baseline potential, which was visualized on left, ventral, and right cortical surface MNI templates in Figure 4D by using the BrainNet Viewer toolbox [35].The color of the sphere represents the magnitude of the baseline potential, and the three representative brain regions are marked with red, bright blue, and orange circles.The baseline potentials of the corresponding brain regions in the left and right hemispheres were basically symmetrical.The brain regions with a baseline of larger than 4 mV were almost located in the parietal and occipital regions and are marked with red color, and those with less than 4 mV were located in the frontal and temporal regions and are marked with blue color, as shown in Figure 4D.
In Figure 4E, the bar plot shows the dominant frequency, and a visualization is shown in Figure 4F.The dominant frequencies of the corresponding brain regions in the left and right hemispheres were also basically symmetrical.The brain regions with an alpha dominant frequency (8-12 Hz) were consistent with those with a baseline frequency above Brain Sci.2024, 14, 240 14 of 21 4 mV in the parietal and occipital regions.The values for the other brain regions were less than 8 Hz, and they are marked with blue color in the frontal and temporal regions.
referred to as hub nodes.The smaller the node strength, the less important it is for information transmission in the brain network.Displaying the simulated time series and corresponding normalized spectra of all 68 brain regions would make the figure very chaotic, so we selected three representative brain regions, as shown in Figure 4B.The simulated output signal of brain region 27 (LH.SF), which had the highest node strength, showed high-amplitude, high-frequency oscillations, with a peak frequency that fell on the alpha band at about 8 Hz; this is marked in red.The output of brain region 31 (LH.FP), which had the lowest node strength, showed characteristics of EEG background noise with a low amplitude and a broadband spectrum; this is marked in bright blue.It can be seen that the coupling strength between the region 31 and other regions was very low in both the simulated rsFC matrix and empirical rsFC matrix in Figure 3C.Therefore, region 31 showed the characteristics of EEG background noise.The output of brain region 51 (RH.POP), which had moderate node strength, is marked in orange.It also showed alpha-wave oscillation, with a lower amplitude and slightly wider spectrum than those of region 27; the dominant frequency was also about 8 Hz.
The baseline potential and dominant frequency of the simulated output signals for all regions are shown in Figure 4C-F.In Figure 4C, the bar plot shows the baseline potential, which was visualized on left, ventral, and right cortical surface MNI templates in Figure 4D by using the BrainNet Viewer toolbox [35].The color of the sphere represents the magnitude of the baseline potential, and the three representative brain regions are marked Comparing Figure 4A,C,E, it can be seen that the brain regions with lower baseline potential values and dominant frequencies corresponded to lower node strength values, indicating that these regions had lower participation in network communication, thus representing a background EEG state.

The Characteristics of the Brain Network Based on Graph Theory
A brain network can be represented as a graph of nodes connected to each other through edges, where nodes represent brain regions and edges represent the strength of functional connections between them.Usually, elements with a functional connection value between nodes greater than the threshold are set to 1, and the elements with a value less than the threshold are set to 0. Then, the constructed adjacency matrix can be used to analyze the graph-theoretical measures of brain networks.
Figure 5A shows the variation curves of the numbers of nodes in the simulated WBNM and empirical brain networks under different thresholds (the calculation details of the adjacency matrix are provided in the Section 2).The blue solid line represents the variation curve of the number of nodes in the simulated brain network, while the red solid line represents the empirical brain network.Simulated brain networks and empirical brain networks exhibit isolated nodes without edge connections when the thresholds are greater than 0.23 and 0.24, respectively.Therefore, we chose to observe the four commonly used graph-theoretical measures with thresholds ranging from 0 to 0.24.
(0.24,68) (0.23,68) 0 0.04 0.08 0.12 0.16 0.2 0.24 0 0.04 0.08 0.12 0.16 0.2 0.24 0 0.04 0.08 0.12 0.16 0.2 0.24 To investigate the similarity of the brain network characteristics between the simulated WBNM and empirical brain networks, the four commonly used graph-theoretical measures were calculated under different thresholds, as shown in Figure 5B (the calculation details of the graph theory measures are provided in the Section 2).When the threshold was less than 0.14 for the simulated WBNM and less than 0.1 for the empirical brain network, as marked with the blue and red arrows, respectively, the average degree D was equal to 67, indicating that all nodes in the brain network had edge connections.In the same case, the characteristic path lengths L were equal to 1, which indicated that the shortest path lengths between all nodes in the brain network were one edge, the clustering coefficients CC were equal to 1-indicating that nodes in the network were highly clustered-and the clustering coefficients E were equal to 1-indicating that the information transmission efficiencies between the brain network nodes were the highest.From the analysis of the four graph-theoretical measures, including the average degree, characteristic path length, clustering coefficient, and global efficiency, it was evident that the WBNM and empirical brain networks exhibited strong similarities in their network characteristics.However, at certain thresholds, the simulated network demonstrated a slightly higher efficiency in information transmission compared to that of the empirical network.Beyond a certain threshold, this relationship was reversed.To investigate the similarity of the brain network characteristics between the simulated WBNM and empirical brain networks, the four commonly used graph-theoretical measures were calculated under different thresholds, as shown in Figure 5B (the calculation details of the graph theory measures are provided in the Section 2).When the threshold was less than 0.14 for the simulated WBNM and less than 0.1 for the empirical brain network, as marked with the blue and red arrows, respectively, the average degree D was equal to 67, indicating that all nodes in the brain network had edge connections.In the same case, the characteristic path lengths L were equal to 1, which indicated that the shortest path lengths between all nodes in the brain network were one edge, the clustering coefficients CC were equal to 1-indicating that nodes in the network were highly clustered-and the clustering coefficients E were equal to 1-indicating that the information transmission efficiencies between the brain network nodes were the highest.From the analysis of the four graph-theoretical measures, including the average degree, characteristic path length, clustering coefficient, and global efficiency, it was evident that the WBNM and empirical brain networks exhibited strong similarities in their network characteristics.However, at certain thresholds, the simulated network demonstrated a slightly higher efficiency in information transmission compared to that of the empirical network.Beyond a certain threshold, this relationship was reversed.

Small-World Network Properties
Watts and Strogatz [36] proposed the concept of small-world networks in 1998.Smallworld networks lie somewhere between regular networks and random networks, and they have the characteristics of a high clustering coefficient and a low characteristic path length.Furthermore, Liao et al. [37] demonstrated that complex human brain networks also have small-world network properties.
To measure whether the simulated WBNM and empirical brain network constructed in this experiment had small-world network properties, the clustering coefficients and characteristic path lengths of simulated random brain networks and empirical random brain networks with the same numbers of nodes and edges were calculated under different thresholds, as shown in Figure 6A,B.As the threshold increased, the clustering coefficients CC of all four types of brain networks showed a downward trend, and the clustering coefficients CC of the simulated and empirical brain networks were greater than those of random brain networks.The characteristic path lengths L of the four types of brain networks showed an upward trend, and the characteristic path lengths L of the simulated and empirical brain networks were similar to those of random brain networks.Watts and Strogatz [36] proposed the concept of small-world networks in 1998.Smallworld networks lie somewhere between regular networks and random networks, and they have the characteristics of a high clustering coefficient and a low characteristic path length.Furthermore, Liao et al. [37] demonstrated that complex human brain networks also have small-world network properties.
To measure whether the simulated WBNM and empirical brain network constructed in this experiment had small-world network properties, the clustering coefficients and characteristic path lengths of simulated random brain networks and empirical random brain networks with the same numbers of nodes and edges were calculated under different thresholds, as shown in Figure 6A,B.As the threshold increased, the clustering coefficients CC of all four types of brain networks showed a downward trend, and the clustering coefficients CC of the simulated and empirical brain networks were greater than those of random brain networks.The characteristic path lengths L of the four types of brain networks showed an upward trend, and the characteristic path lengths L of the simulated and empirical brain networks were similar to those of random brain networks.The small-world network index is expressed as the ratio of the clustering coefficient ratio divided by the characteristic path length ratio, and it is used to determine whether a network has small-world properties.A small-world network index σ that is greater than 1 indicates that the network has small-world properties.The small-world network indices of the simulated WBNM and empirical brain networks with different thresholds are shown in Figure 6C.As the threshold increased, both small-world network indices The small-world network index is expressed as the ratio of the clustering coefficient ratio divided by the characteristic path length ratio, and it is used to determine whether a network has small-world properties.A small-world network index σ that is greater than 1 indicates that the network has small-world properties.The small-world network indices of the simulated WBNM and empirical brain networks with different thresholds are shown in Figure 6C.As the threshold increased, both small-world network indices showed an upward trend and were greater than 1, which indicated that both the simulated WBNM and empirical network had small-world network properties.
The Pearson correlation coefficients of the functional connectivity matrix between the simulated WBNM and empirical networks with different thresholds are shown in Figure 6D.As the threshold increased, the correlation first decreased and then increased.At the threshold of 0.22, the correlation achieved its maximum value of 0.709, indicating a high degree of similarity between the functional connectivity matrixes of the simulated WBNM and empirical networks at this threshold.

Discussion
The purpose of this study was to construct a whole-brain network that was close to an actual brain in terms of EEG rhythm, functional connectivity, and brain network characteristics.In this way, researchers can use this whole-brain model to study the neural mechanisms of EEG signal generation, the pathogenic mechanisms of various functional diseases of the brain, and various methods of brain stimulation regulation from a multi-scale perspective.Using computational models of brain dynamics to construct a whole-brain model is not an unprecedented effort.In fact, our work was inspired by previous largescale modeling studies, which elucidated key insights into using various other models to construct an entire brain network.The innovation of this study lies in its introduction of a Wendling model with richer dynamic characteristics to represent the nodes of the brain network, and at the same time, real MRI data and an optimized global coupling coefficient were used to simulate the connections of the brain network, thus constructing a WBNM that was closer to the real brain in terms of brain rhythm, functional connectivity, and network characteristics.
In this study, we successfully constructed a comprehensive brain model to simulate a real whole-brain EEG network, and it was able to balance physiological accuracy and computational processability.Firstly, we partitioned the entire brain based on the Desikan-Killiany brain atlas and then calculated the quantity of white-matter fiber bundles found between different brain regions, which enabled us to derive a relative coupling strength matrix for whole-brain coupling.By coupling 68 Wendling NMMs and optimizing a global coupling coefficient of C = 20.3,we achieved the optimal simulation for the functional connectivity matrix of a real brain network.Then, we analyzed the dynamics of the simulated EEG and employed graph-theoretical metrics to further characterize the attributes of the constructed WBNM, namely, the average degree, characteristic path length, clustering coefficient, global efficiency, and small-world properties.We found that the EEG spectra and network properties of the WBNM constructed in this study were highly similar to those of a real brain network.
In this study, the Wendling model was selected as the node of the whole-brain network, and it made an important contribution to the development of the WBNM.To the best of our knowledge, the utilization of the Wendling model for computational modeling at the wholebrain scale has not been previously explored.When simulating a neural system from a microscopic perspective, it is usually necessary to model a large number of spike discharges of neuronal membrane voltage, which can cause significant computational complexity.It is also difficult to determine how a large number of neurons connect to each other.Therefore, some studies chose to use the "mean field theory" neural mass model as the node of the whole-brain network.An NMM uses a set of differential equations to represent the coupled excitatory and inhibitory neuronal populations, which can represent the true physiological meaning of neuronal populations, such as the population firing rate and postsynaptic current.However, in most previous studies and research involving whole-brain problems, in order to reduce the complexity of WBNMs, the selected models of individual neural units were highly idealized, and the structure of neural group units was overly simplified, as they usually consisted of basic excitatory and inhibitory loops, such as in the Wilson-Cowen model, Later Breakspeed model, Arm model, and Jansen model [1,2,[24][25][26]28].As a result, the whole-brain models that they constructed had certain differences from the actual brain structure.In comparison, the Wendling model can generate six common EEG signals, including background noise, alpha waves, epileptic spike waves, etc., thus enhancing the physiological and physical reality in certain aspects.The signals simulated in the Wendling model have rich dynamic characteristics and wider frequency bands and are closer to real human EEG features [38][39][40][41].Ursino et al. [22] proposed an expanded Wendling model to faithfully replicate actual brain EEG signals, they and obtained different rhythm combinations (β and γ, α and γ, or a wide spectrum).However, the Ursino model is somewhat complex.
Of course, the Wendling model also has some limitations of its own.Although the Wendling model can simulate EEG signals with rich dynamic characteristics, there is still a certain gap with respect to the simulation of real EEG signals because the characteristics of EEG signals are very complex.For example, the Wendling model cannot generate gamma waves and ripple waves.The Wendling model is also a model with very complex bifurcation dynamics.Improper parameter settings and mutual coupling that is too strong may cause the Wendling model to be in an abnormal state.The Wendling model currently has a wealth of research results focused on epilepsy, with relatively few studies on other functional disorders of the brain and cognition.
We selected the Wendling model as the node of the WBNM.Therefore, using the Wendling model as the node of WBNM was beneficial for reproducing a real-world brain neural network while also maintaining a balance between authenticity and computational complexity.In this study, the simulated EEG exhibited rich rhythmic characteristics, including alpha oscillations, a wider spectrum, and the background state (Figure 4).The results indicated that using the Wendling model as a node made the EEG features generated by the WBNM closer to reality.Moreover, by adjusting the parameters of the Wendling model, the dynamic states of different regions could also be adjusted, which can be beneficial when using this WBNM to conduct research on epilepsy transmission and other phenomena.
Another advantage of the WBNM constructed in this study is based on the clear anatomical structure and the use of real MRI signals to form a structural connectivity matrix.The fiber tracts that were calculated from MRI data could represent the relative coupling strength between different regions.However, during the modeling process, the value of the fiber tracts was different from the value of the coupling coefficient between the regions.To keep the relative coupling strength between regions unchanged, we introduced the global coupling coefficient C, which used a multiplicative factor to adjust the specific coupling strength between regions in the model.By calculating the correlation between the functional connectivity matrix of the WBNM and the functional connectivity matrix of the empirical data, the value of C = 20.3 at the strongest correlation point was determined as the optimal global coupling coefficient.By setting the global coupling coefficient to C = 20.3, the rsFC matrices of the simulated WBNM and empirical brain network had strong similarity, and the PCC value reached 0.676 (Figure 3).Therefore, the simulated brain network and empirical brain network had similar graph-theoretical measurements (Figure 5).Both of them have small-world properties, and when the threshold value was 0.22, there was a strong correlation between the simulated rsFC matrix and the empirical rsFC matrix with PCC = 0.709 (Figure 6).All of the above results indicate that the connection characteristics and network characteristics of the WBNM proposed in this study are very similar to those of actual brain networks.In some platforms that built whole-brain models [1,2], different models could be selected as the nodes of the whole-brain network (these models were relatively simple compared to Wendling models), and real structural or functional datasets could be used to form a functional coupling matrix.However, there is no evidence to prove whether the WBNMs that they constructed would generate brain-like EEG waveforms or real-world brain network features.
There are also some limitations in this study.For example, although the flexibility of the code in our constructed WBNM gave it a certain scalability, we still needed to set the global coupling coefficient C to make it fit an actual neural system in terms of EEG waveforms at each node, functional connections, and brain network characteristics.When other subjects' structural connectivity data are used, it is necessary to select a new global coupling coefficient.If we want to change or add brain regions, in the initial step of modeling, which is the calculation of the fiber tracts, we need to recalculate, and the entire model construction process needs to be repeated, which is very time-consuming and laborious.However, we believe that in order to build a more realistic WBNM, this work is necessary.Building a WBNM that simulates the real human brain is not a simple task.Further optimizing the behavior of brain regions requires more complex node models or the optimization of the parameters of each node model.Due to the large number of model parameters, optimizing parameter settings at the whole-brain level is very difficult.However, we can start by reverse-identifying the model parameters of one or several brain regions in order to make the WBNM as close to the actual situation as possible.Our study only constructed a WBNM and did not use this model to study the neurophysiological mechanisms of the brain or functional diseases of the brain.The application aspect is also the direction of our future work.
In future work, our team and other researchers can use this model to perform meaningful work.For example, by adjusting the parameters of the Wendling model, including the excitability, inhibition, and synaptic connectivity coefficient, as well as the structural coupling matrix, we can investigate the mechanisms of reduced brain rhythms, complexity, and functional connectivity in AD and MCI patients.Thus, we can conduct research on the brain functional mechanisms and intervention measures in AD and MCI patients.By adjusting the parameters of the Wendling model, including the excitability, inhibition, and synaptic connectivity coefficient, as well as the structural coupling matrix, we can investigate the brain mechanisms of the reduced brain rhythm, complexity, and functional connectivity in AD and MCI patients.By adding electrical stimulation to appropriate parts of this WBNM, we can investigate the possibility of intervening in cognitive decline.For example, researchers can work on applying direct-current electrical stimulation to this WBNM to investigate the effects of different intensities and positions of electrical stimulation on the waveforms, dynamic behavior, and network characteristics of the whole-brain model's output.Similarly, we can also simulate the spread of epilepsy in the brain and study feasible methods for the closed-loop control of seizures.

Conclusions
In this study, we constructed a resting-state whole-brain network model (WBNM) with 68 brain regions using the Wendling model for the nodes with real MRI data forming the interconnections.By fine-tuning the global coupling coefficient to optimize the correlation between the simulated rsFC matrix and the empirical rsFC matrix, the WBNM had functional connections that were similar to real EEG signals.The EEG signals of network nodes with different node degrees exhibited different rhythms, bandwidths, baseline potentials, and dominant frequencies.Then, we compared the differences between the simulated and empirical brain networks at different thresholds using four widely used graph-theoretical measures.The WBNM that we constructed had network features that were similar to those of empirical neural networks, including the average degree, characteristic path length, clustering coefficient, and global efficiency, as well as small-world network attributes.We suggest that the constructed Wendling WBNM can reproduce real brain networks to a certain extent and can be used to study the neurophysiological mechanisms of cognitive and functional disorders of the brain.

Figure 1 .
Figure 1.Anatomical brain network data and schematic diagram of the whole-brain network model.(A) The cortical parcellation of 68 brain regions.(B) The fiber tracts between region pairs.(C) The empirical structural connectivity matrix.(D) The schematic structure of the Wendling model.(E) The schematic structure of the WBNM.

Figure 1 .
Figure 1.Anatomical brain network data and schematic diagram of the whole-brain network model.(A) The cortical parcellation of 68 brain regions.(B) The fiber tracts between region pairs.(C) The empirical structural connectivity matrix.(D) The schematic structure of the Wendling model.(E) The schematic structure of the WBNM.

Figure 2 .
Figure 2. A block diagram of the whole-brain network model.

Figure 2 .
Figure 2. A block diagram of the whole-brain network model.

Figure 3 .
Figure 3.The WBNM at the optimal global coupling coefficient fitting point.(A) The PCC between the simulated rsFC matrix and empirical rsFC matrix with different global coupling coefficients.(B) Linear regression of the simulated rsFC matrix and empirical rsFC matrix.(C) The simulated rsFC matrix and empirical rsFC matrix.

Figure 3 .
Figure 3.The WBNM at the optimal global coupling coefficient fitting point.(A) The PCC between the simulated rsFC matrix and empirical rsFC matrix with different global coupling coefficients.(B) Linear regression of the simulated rsFC matrix and empirical rsFC matrix.(C) The simulated rsFC matrix and empirical rsFC matrix.

1 BF
27th region,LH.SF 31th region,LH.FP 51th region,RH.POP D The node strength of 68 regions The simulated time series and normalized spectra of three typical brain regions The baseline of 68 regions The visualization of the baseline potential The dominant frequency of 68 regions 10

Figure 4 .
Figure 4.The dynamics of the simulated output signals for all regions.(A) A bar plot of the node strength.(B) The simulated time series and normalized spectra of 27, 31, and 51 brain regions.(C) A bar plot of the baseline potential for all brain regions.(D) A visualization of the baseline potential from the left, ventral, and right angles of view.(E) A bar plot of the dominant frequencies for all brain regions.(F) A visualization of the dominant frequencies from the left, ventral, and right angles of view.

Figure 4 .
Figure 4.The dynamics of the simulated output signals for all regions.(A) A bar plot of the node strength.(B) The simulated time series and normalized spectra of 27, 31, and 51 brain regions.(C) A bar plot of the baseline potential for all brain regions.(D) A visualization of the baseline potential from the left, ventral, and right angles of view.(E) A bar plot of the dominant frequencies for all brain regions.(F) A visualization of the dominant frequencies from the left, ventral, and right angles of view.

Figure 5 .
Figure 5.The characteristics of the simulated WBNM and empirical brain network based on graph theory.(A) The node numbers with different thresholds for the two networks.(B) Four graph-theoretical measures with different thresholds for the two networks.

Figure 5 .
Figure 5.The characteristics of the simulated WBNM and empirical brain network based on graph theory.(A) The node numbers with different thresholds for the two networks.(B) Four graphtheoretical measures with different thresholds for the two networks.

24 Figure 6 .
Figure 6.Analysis of the small-world characteristics of the simulated WBNM and empirical brain network.(A) The clustering coefficients of the simulated WBNM, simulated random WBNM, empirical brain network, and empirical random brain network with different thresholds.(B) The characteristic path lengths of the simulated WBNM, simulated random WBNM, empirical brain network, and empirical random brain network with different thresholds.(C) The small-world network indexes of the simulated WBNM and empirical brain network with different thresholds.(D) The correlation of the functional connectivity matrix between the simulated WBNM and empirical networks with different thresholds.

Figure 6 .
Figure 6.Analysis of the small-world characteristics of the simulated WBNM and empirical brain network.(A) The clustering coefficients of the simulated WBNM, simulated random WBNM, empirical brain network, and empirical random brain network with different thresholds.(B) The characteristic path lengths of the simulated WBNM, simulated random WBNM, empirical brain network, and empirical random brain network with different thresholds.(C) The small-world network indexes of the simulated WBNM and empirical brain network with different thresholds.(D) The correlation of the functional connectivity matrix between the simulated WBNM and empirical networks with different thresholds.

Table 2 .
Standard values and physiological interpretations of the parameters in the Wendling model.