Figures and data
Schematic illustration of brain network model simulations on central processing units and graphics processing units.
Each cell represents a network node at a given simulation time in a given simulation. Dark-colored cells indicate completed calculations, and arrows show the order of execution. Created with BioRender.com.
Optimization of brain network models.
A brain network model describes the dynamics of a network of nodes connected via a structural connectome. As an example, the reduced Wong-Wang model is shown, with three free parameters: G, wp, and JN. The goal of optimization is to fit model parameters to empirical data, such as blood-oxygen-level-dependent (BOLD) signal. Candidate parameter sets are sampled, and parallel simulations are run for each set. The simulated BOLD signals are compared to empirical data using evaluation metrics such as the correlation between simulated and empirical functional connectivity (FC) matrices (FCcorr) and the Kolmogorov-Smirnov distance between simulated and empirical functional connectivity dynamics (FCD) matrices (FCDKS). Additional penalty terms may also be included (e.g., the feedback inhibition control (FIC) penalty in the reduced Wong-Wang model). The evaluation metrics are combined to define the cost function of each simulation. In evolutionary optimization, new parameter sets are iteratively sampled based on cost functions from previous iterations, whereas in grid search, all parameter combinations are simulated in a single run on a predefined parameter grid. After completion of the optimization procedure, the optimal simulation is defined as the one with the lowest cost function.
Grid search of the homogeneous model.
Three-dimensional grids show the variation of cost function and its components (A) as well as the average state variables (B) as a function of the parameters G, wp, and JN. Each dimension is sampled at 22 evenly spaced points. Red spheres mark the optimum, i.e., the simulation with the lowest cost. In A, the cost function, FICpenalty and FCDKS are negated so that in all grids brighter colors indicate better fits.
Evolutionary optimization of the homogeneous model.
(A) Three-dimensional visualization of negative cost values across all parameter combinations evaluated during CMA-ES optimization. (B) Three-dimensional visualization of parameter combinations sampled in eight selected generations, showing gradual convergence of the optimizer toward the optimal region of the parameter space. (C) Evolution of negative cost across particles (dots) and generations. The red line indicates the median negative cost per generation.
Evolutionary optimization of heterogeneous models.
(A) In the map-based heterogeneous model, regional parameters 
Model fit, reliability, and heritability of individualized BNMs.
(A) Comparison of negative cost and its components across homogeneous, map-based heterogeneous, and node-based heterogeneous models. Dots denote the optimal simulations of each subject and model. Horizontal lines indicate significant pairwise differences between models, with the arrow direction showing the better fit. (B) Left : Intraclass correlation coefficients (ICCs) for empirical (red) and simulated (grey) features across day 1 and day 2 sessions. Right : Kolmogorov-Smirnov (KS) distance of functional connectivity dynamics (FCD) distributions between sessions for empirical (red) and simulated (grey) data. (C) Heritability (h2) of empirical (red) and simulated (grey) features based on the averaged data from both sessions. In B and C, violin plots represent the distribution of values across connectome edges or network nodes, with black lines indicating medians and quartiles. Diamonds indicate point estimates.
Scaling of GPU and CPU compute time with the number of simulations and network size.
(A) Compute time as a function of the number of simulations (log-log scale). Solid lines show observed values and dotted lines indicate linear projections. (B) Speed-up relative to a single CPU thread as a function of the number of simulations (log scale on the x-axis). (C) Compute time as a function of network size (log-log scale). The grey band indicates the range of linear growth. GPU compute times increase linearly with the number of nodes, whereas CPU compute times increase supralinearly. (D) Comparison of compute times for simulations with 100 versus 10,000 nodes (log scale on the y-axis). Arrows indicate the fold-increase in compute time between the two network sizes.
Program structure.
Using the Python interface, the user defines a BNMProblem and assigns it to an Optimizer (of type GridOptimizer or PymooOptimizer). The BNMProblem is linked to a model-specific SimGroup (e.g., rWWSimGroup). SimGroup.run() calls the core run_simulations, which instantiates the corresponding C++ model (e.g., rWWModel), which runs the simulation on GPU/CPU. The core returns simulated BOLD and, when requested, FC and FCD, to the SimGroup instance, which scores them against empirical data via SimGroup.score(). The resulting cost function components are passed back to the BNMProblem and its Optimizer. In evolutionary optimization, these values guide the sampling of subsequent generations, with iterations continuing until the maximum generation count or an early termination criterion is reached. For grid search, a single exhaustive iteration is performed. In both cases, once the optimization completes, the optimal parameters are identified and returned to the user.
