# Arbor-TVB: a novel multi-scale co-simulation framework with a case study on neural-level seizure generation and whole-brain propagation

Thorsten Hater, Juliette Courson, Han Lu, Sandra Diaz-Pier, Thanos Manos

► To cite this version:

Thorsten Hater, Juliette Courson, Han Lu, Sandra Diaz-Pier, Thanos Manos. Arbor-TVB: a novel multi-scale co-simulation framework with a case study on neural-level seizure generation and whole-brain propagation. Frontiers in Computational Neuroscience, 2026, 19, pp.1731161. 10.3389/fn-com.2025.1731161. hal-05489179

HAL Id: hal-05489179

https://hal.science/hal-05489179v1

Submitted on 2 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_2.png" /></figure>

Check for updates

OPEN ACCESS

EDITED BY
Tiago Ribeiro,
National Institutes of Health (NIH), United
States

REVIEWED BY
Gerald Kaushallye Cooray,
Karolinska Institutet (KI), Sweden
Gian Marco Duma,
Eugenio Medea (IRCCS), Italy

*CORRESPONDENCE

Thorsten Hater

<figure><img src="image_3.png" /></figure>

t.hater@fz-juelich.de

Thanos Manos

<figure><img src="image_4.png" /></figure>

thanos.manos@cyu.fr

RECEIVED 23 October 2025

REVISED 15 December 2025

ACCEPTED 17 December 2025

PUBLISHED 02 February 2026

CITATION

Hater T, Courson J, Lu H, Diaz-Pier S and
Manos T (2026) Arbor-TVB: a novel
multi-scale co-simulation framework with a
case study on neural-level seizure generation
and whole-brain propagation.

Front. Comput. Neurosci. 19:1731161.

doi: 10.3389/fncom.2025.1731161

COPYRIGHT

© 2026 Hater, Courson, Lu, Diaz-Pier and
Manos. This is an open-access article
distributed under the terms of the Creative
Commons Attribution License (CC BY). The
use, distribution or reproduction in other
forums is permitted, provided the original
author(s) and the copyright owner(s) are
credited and that the original publication in
this journal is cited, in accordance with
accepted academic practice. No use,
distribution or reproduction is permitted
which does not comply with these terms.

Arbor-TVB: a novel multi-scale
co-simulation framework with a
case study on neural-level
seizure generation and
whole-brain propagation

Thorsten Hater¹*, Juliette Courson²,³, Han Lu¹, Sandra Diaz-Pier¹ and Thanos Manos²*

¹Simulation and Data Lab Neuroscience, Jülich Supercomputing Centre (JSC), Forschungszentrum Jülich GmbH, Jülich, Germany, ²ETIS Lab, ENSEA, CNRS, UMR8051, CY Cergy-Paris University, Cergy, France, ³Department of Computer Science, University of Warwick, Coventry, United Kingdom

Computational neuroscience has traditionally focused on isolated scales, limiting understanding of brain function across multiple levels. While microscopic models capture biophysical details of neurons, macroscopic models describe large-scale network dynamics. Integrating these scales, however, remains a significant challenge. In this study, we present a novel co-simulation framework that bridges these levels by integrating the neural simulator Arbor with The Virtual Brain (TVB) platform. Arbor enables detailed simulations from single-compartment neurons to populations of such cells, while TVB models whole-brain dynamics based on anatomical features and the mean neural activity of a brain region. By linking these simulators for the first time, we provide an example of how to model and investigate the onset of seizures in specific areas and their propagation to the whole brain. This framework employs an MPI intercommunicator for real-time bidirectional interaction, translating between discrete spikes from Arbor and continuous TVB activity. Its fully modular design enables independent model selection for each scale, requiring minimal effort to translate activity across simulators. The novel Arbor-TVB co-simulator allows replacement of TVB nodes with biologically realistic neuron populations, offering insights into seizure propagation and potential intervention strategies. The integration of Arbor and TVB marks a significant advancement in multi-scale modeling, providing a comprehensive computational framework for studying neural disorders and optimizing treatments.

KEYWORDS

Arbor, mouse brain connectome, multi-scale neural models, seizures, The Virtual Brain

# 1 Introduction

The human brain consists of billions of neurons and an equally vast population of non-neuronal cells, intricately organized into layers and regions (Herculano-Houzel, 2009, 2012). Each neuron operates as a highly sophisticated biochemical machinery (West et al., 2002; Augustine et al., 2003; Darnell, 2013; Lu et al., 2025), coordinating signal transmission within an extensive network in health (Reyes, 2003; Barral et al., 2019; Dicks, 2022) and disease (see e.g., Tetzlaff et al., 2025). Ever since the Hodgkin-Huxley model was introduced to describe membrane potential dynamics (Hodgkin and Huxley, 1952c),
computational neuroscience has played a pivotal role in enhancing our understanding of brain function. Yet, due to the immense complexity of the brain and computational constraints, most modeling studies focus on a single scale simulator or rely on standalone simulation codes.

Modeling the large-scale electrical activity of the brain is a complex task. It not only demands familiarity with advanced mathematical methods, but also a solid grasp of the brain's physiology and anatomy. Neural field theory offers a way to study the nonlinear behavior of large groups of neurons at a population level, while still keeping the mathematics manageable. These models give us a strong theoretical framework for understanding key processes in neural tissue, including how the brain transitions between different activity states, such as those seen in sleep or during epileptic events, see e.g., Cook et al. (2022) for a recent review. Moreover, multi-scale computational modeling provides a framework for connecting neural mechanisms with measurements ranging from unit recordings to electroencephalogram (EEG), magnetoencephalography (MEG), and functional magnetic resonance imaging (fMRI). Such models clarify how neural systems compute and interact, and they are essential for integrating empirical findings into a robust theoretical understanding of brain function, see e.g., Deco et al. (2008) and Cooray et al. (2023). Along this direction, recently in Cooray et al. (2025), the authors also studied oscillatory activity in cortical tissue arising from not uniform neural connections. and showed that oscillations can be maintained under a wide range of anisotropic and time-varying connectivity patterns.

Simulations incorporating biophysical properties and neural morphology typically concentrate on individual neurons using simulators such as the NEURON simulator (Carnevale and Hines, 2006). At the microscopic level, studies have explored questions such as how the ion channel kinetics influence neural excitability (see e.g., Gurkiewicz et al., 2011; Suma et al., 2024), how proteins, enzymes, and calcium concentration are distributed among neighboring spines to impact plasticity (see e.g., Luboeinski and Tetzlaff, 2021; Chater et al., 2024), and how signal propagation along axonal fibers relates to neuropathic pain (see e.g., Tigerholm et al., 2014, 2015). Some studies examine how neural morphology—such as dendritic tree growth (see e.g., Yasumatsu et al., 2008) and morphology-dependent plastic interactions (see e.g., Hananeia et al., 2024)—affects function. These studies, while often limited to small patches of the neural membrane, a few dendritic segments, or a small local network, provide valuable approximations of broader neural phenomena.

At the mesoscopic level, researchers simplify neuronal representations using point leaky-integrate-and-fire neurons [based on simulators such as NEST (Gewaltig and Diesmann, 2007) or Brian/Brian2 (Stimberg et al., 2019)], allowing studies on larger networks without explicit neuronal morphology or with some degree of self-customized morphology, using, e.g., NESTML (Linssen et al., 2024). This approach has advanced our understanding of neural heterogeneity (Demirtaş et al., 2019; Nayebi et al., 2021; Gast et al., 2024), self-organization (Zheng et al., 2013; Diaz-Pier et al., 2016; Miner and Triesch, 2016), neural capacity (Emina and Kropff, 2022), energy efficiency (Sacramento et al., 2015), and neural plasticity in disease and health (Manos

et al., 2021; Lu et al., 2024). Most microscopic and mesoscopic models remain theory-driven, using mathematical approximations to infer neural behavior rather than directly establishing model based on large datasets (see e.g., Popovych et al., 2019 for a recent review).

Data-driven modeling has gained traction at the macroscopic level with the rise of open-source brain imaging databases [such as OpenfMRI (Poldrack and Gorgolewski, 2017)]. High-resolution structural and functional data from magnetic resonance imaging (MRI) and diffusion tensor imaging (DTI) enable whole-brain modeling based on real anatomical features. The Virtual Brain (TVB) (Sanz Leon et al., 2013; Sanz-Leon et al., 2015; Ritter et al., 2013) and Virtual Brain Twin (VBT) (Hashemi et al., 2025) platforms, for instance, integrate functional MRI and DTI datasets to build individualized models, using coupled oscillators to represent regional activity. TVB has contributed significantly to the understanding of neurological disorders and serves as a testing ground for therapeutic interventions (see e.g., Stefanovski et al., 2021; Monteverdi et al., 2023; Courson et al., 2024 and references therein), for studying self-organization on macroscale (see e.g., Fousek et al., 2024), consciousness (see e.g., Breyton et al., 2024) and healthy aging (see e.g., Lavanga et al., 2023).

With advances in computing resources and simulation technologies, the integration of models across different scales has become both feasible and essential to strike a balance between retaining detailed information and achieving a broad-scale understanding. Recently, a co-simulation framework was introduced that employs NEST and TVB to bridge mesoscopic and macroscopic modeling. This work has pioneered cross-scale modeling (Kusch et al., 2024) and has demonstrated the benefits of integrating models across spatial levels. A notable application of the virtual deep brain stimulation model (Meier et al., 2022; Shaheen et al., 2022; Wang et al., 2025) demonstrated its utility in multiscale simulations. Similar tools have been made available within the European digital neuroscience platform, EBRAINS (Schirner et al., 2022).

However, integrating microscopic and macroscopic models remains technically challenging. At the core lies the vast amount of information being processed, billions of cells with thousands of connections, and the immense gap in timescales, from microseconds in ion channel dynamics to minutes or hours for plastic changes of the connectome. To solve this challenge, we used the Arbor simulator (Abi Akar et al., 2019) and, more specifically, its most recent next-generation version (Cumming et al., 2024), at the microscopic end. Designed for single-neuron and large-scale network simulations, Arbor leverages GPU resources to enhance computational speed and energy efficiency.

In this Methods paper, we successfully established efficient communication between Arbor and the TVB that respects their different operational time scales and provided a use case example of the cross-scale interaction. To demonstrate a first showcase, we used a mouse brain connectome provided by TVB, where each region represents the mean mass neural activity of a brain area modeled by a macroscopic model. Using the co-simulation interface, we replaced one TVB node with a network of detailed neurons modeled in Arbor. An extended Hodgkin-Huxley-based neuron model (Depannemaecker et al., 2022) was utilized in Arbor
to simulate different neural activity patterns (e.g., spiking, bursting, seizure-like, etc.). By tweaking a single parameter, we showcased that the seizure-like events generated in Arbor propagated to other nodes modeled in TVB. This platform provides users with the freedom to use existing models across scales with minimal additions and enabled future development in building brain digital twins that contain both micro- and macroscopic information for therapeutic applications.

# 2 Materials and methods

## 2.1 The Arbor simulator

Arbor is an open-source library for building simulations of biophysically detailed neuron models (Abi Akar et al., 2019). It provides an alternative to software like NEURON (Carnevale and Hines, 2006), but with a strong emphasis on modern hardware and scalability to large-scale systems (Hines, 1984). Its overall set of capabilities allows Arbor to model neural networks at a level of resolution beyond point models to explore phenomena like dendritic computation. Thus, support for bulk-synchronous parallelism via MPI, shared memory parallelism by utilizing a thread-pool and job system is central to Arbor, and certain cell types—primarily cable cells—can further leverage SIMD and GPU hardware. Arbor is written in C++, though most users interface with it through an intuitive, high-level Python interface built on top of the lower level implementation. The underlying numerical model of Arbor is the cable equation:

$$c \frac{\partial U}{\partial t} = \frac{\partial}{\partial x} \left( \sigma \frac{\partial U}{\partial x} \right) + i \quad (1)$$

where the membrane potential U is computed over the morphological structure of the neural tree; the spatial coordinate x and the derivative are to be understood within this structure (Von Helmholtz, 1850; Thompson and Kelvin, 1855; Hodgkin et al., 1952; Hodgkin and Huxley, 1952a,b; Loligo, 1952; Hodgkin and Huxley, 1952c; Scott, 1975). The parameters c and σ define the membrane capacitance and longitudinal resistance. The trans-membrane current density i models the entirety of ionic and non-ionic currents. In both NEURON and Arbor, these are calculated from user-specified sets of differential equations, potentially varying along the morphology. The equations for i and U are solved in alternation (Lie-Trotter splitting) using a first-order implicit method.

## 2.2 Single neuron and network models in Arbor

We begin by selecting a dynamical model that allows for relatively easy yet realistic simulation of a broad spectrum of neural activity at the single-neuron level, governed by a small set of biophysical parameters. The neural model from Depannemaecker et al. (2022) was formulated for Arbor in the Neuron MODeling Language (NMODL). The following equations form the slow part of the system, describing the evolution of ion concentrations due to voltage-gated channels, active pumps, and buffering by an external

bath, see Figure 1 for a schematic of the dynamics. It describes the ionic exchanges between the intracellular and extracellular spaces (ICS, ECS) of a neuron immersed within an external bath, acting as a potassium buffer of concentration K<sub>bath</sub>. Ions flow between the ICS and ECS through a sodium-potassium pump and the sodium, potassium and chloride voltage-gated channels, driving changes in the internal (K<sub>i</sub>, Na<sub>i</sub>, Cl<sub>i</sub>) and external (K<sub>o</sub>, Na<sub>o</sub>, Cl<sub>o</sub>) ionic concentrations. By gradually increasing the external bath concentration of potassium ions K<sub>bath</sub>, the model sequentially presents these patterns: resting state (RS), spike train (ST), tonic spiking (TS), bursting, seizure-like events (SLE), sustained ictal activity (SIA) and depolarization block (DB), see Figure 2. The fast dynamics of the membrane potential V is modeled in Arbor via the cable equations, see above, which require computing the ion current densities i<sub>X</sub> used in the simulator update as:

$$i_X = g_X(V - E_X) \quad (2)$$

$$E_X = C \cdot \log \left( \frac{X_o}{X_i} \right) \quad (3)$$

with the ion species X = {K, Na, Cl} and a non-ion current density:

$$i_{\text{pump}} = \frac{\rho}{(1 + \exp(10.5 - 0.5Na_i))(1 + \exp(5.5 - K_o))} \quad (4)$$

These currents enter the cable equation Equation 1 as the trans-membrane current i via:

$$i = i_{\text{pump}} + \sum_{X} i_{X}$$

Following the Hodgkin-Huxley model, conductivities g<sub>X</sub> are written as:

$$g_K = g_{0,K}n + g_{l,K} \qquad g_{Na} = g_{0,Na}mh + g_{l,Na} \qquad g_{Cl} = g_{0,Cl} \quad (5)$$

The—internal i and external o—ion concentrations are modeled as:

$$\begin{align}
K_i &= K_{0,i} + \Delta K_i & Na_i &= Na_{0,i} - \Delta K_i & Cl_i &= Cl_{0,i} \nonumber \\
K_o &= K_{0,o} - \beta \Delta K_i + K_g & Na_o &= Na_{0,o} + \beta \Delta K_i & Cl_o &= Cl_{0,o} \tag{6}
\end{align}$$

The variables {ΔK<sub>i</sub>, K<sub>g</sub>} evolve as:

$$\frac{d\Delta K_i}{dt} = -\gamma (i_K - i_{\text{pump}}) \quad (7)$$

$$\frac{dK_g}{dt} = \epsilon (K_{\text{bath}} - K_o) \quad (8)$$

where γ converts the current density i<sub>X</sub> to molar flux, summarizing the effect of the ion pump in Figure 1A and the external buffer. Finally, fast dynamics were reduced and adjusted to mammalian neurons:

$$\frac{dn}{dt} = \frac{1}{\tau} (n - n_{\infty}(V)) \quad (9)$$

$$n_{\infty}(V) = \frac{1}{1 + \exp(-(19 + V)/18)} \quad (10)$$

$$m = m_{\infty}(V) = \frac{1}{1 + \exp(-(2 + V/12))} \quad (11)$$

$$h = h(n) = 1.1 - \frac{1}{1 + \exp(3.2 - 0.8n)} \quad (12)$$
<figure><img src="image_5.png" /><figcaption>FIGURE 1</figcaption></figure>

Biophysical neuron model and Arbor network. (A) For single cell dynamics, three ion concentrations (K, Na, and Cl) are modeled in the cell's interior and a thin shell of its extracellular medium. The latter is, in turn, surrounded by a bath of a fixed potassium concentration $K_{\text{bath}}$. The model simulates changes to the concentration in addition to the current contributions based on three voltage-gated ion channels, an active pump between potassium and sodium, and the buffering effect of the surrounding potassium bath. (B) We choose typical values of $K_{\text{bath}}$ for the single models to generate the tonic spiking and seizure-like event behaviors. In most cases, a fully connected network using exponential synapses with weight $w = 0.5$ is used. As an example, we show here the network instantiation for a network size $N = 5$ and a ratio of SLE to tonic neurons of $f = 0.2$.

based on the observations that the reaction of the sodium gating variable to changes in V is high instantaneous and $h(t) + n(t) = \text{const}$.

The resulting ion channel was added to a basic, spherical, single-compartment neuron. After implementing this biophysical model, we reproduced the firing patterns using the parameters of the reported model (Figures 2A–F), see also Depannemaecker et al. (2022) for more details and motivation for model parameter choices. Note that despite the values given in the original publication, neither the Arbor nor the published reference implementation produces the depolarization block pattern at $K_{\text{bath}} = 20\text{mM}$ but only at around $K_{\text{bath}} = 22.5\text{mM}$. From here, a simple model network was developed, comprising N total cells, with a mixture of tonic $f \cdot N$ and SLE $(1-f) \cdot N$ cells, where both sub-populations are assigned individual values for $K_{\text{bath}}$, sketched in Figure 1. Cells are connected using exponential synapses with an internal weight of $w = 0.5$ chosen to produce an activity similar to Rabuffo et al. (2025) which uses delta synapses.

In addition to the elementary spherical morphology, we also investigated two multi-compartmental neuron models. The first model included a single dendritic segment of 25µm, subdivided into 5µm compartments. The second model extended the dendrite into a randomly generated tree composed of 5µm compartments. To introduce variability across cells, the random number generator was seeded with each cell's unique identifier (see Figure 3 for examples). In future work, these synthetic morphologies will be replaced with reconstructions derived from neural imaging data available in public databases. Supporting such models in Arbor will require only a simple command to load per-cell morphology data from disk.

## 2.3 TVB network model

Following Deco et al. (2013), we used the reduced Wong-Wang model (Wang, 2002) to simulate resting-state activity and to

investigate the dynamics of local brain regions embedded within a large-scale brain network. The mean firing rate $H(x_I)$ and mean synaptic gating variable $S_I$ of region I are described by:

$$ \frac{dS_I}{dt} = -\frac{S_I}{\tau_s} + (1 - S_I)\gamma H(x_I) \quad (13) $$

$$ H(x_I) = \frac{ax_I - b}{1 - \exp(-d(ax_I - b))} \quad (14) $$

$$ x_I = \omega J_N S_I + G J_N \sum_K c_{IK} S_K + I_0, \quad (15) $$

where $x_I$ is the synaptic input to the I-th region, $\omega = 1$ denotes the local excitatory recurrence, $c_{IK}$ is the strength of the structural connection from the local area I to K, and G is a global coupling strength. The parameters are set to the same values as those used in the TVB implementation. $J_N = 0.2609\text{nA}$ is the synaptic coupling of NMDA receptors and $I_0 = 0.33\text{nA}$ is the baseline external input. The kinetic parameters are $\tau_s = 100\text{ms}$ and $\gamma = 0.641$. The parameters of the input-output function H are $a = 0.27\text{nC}^{-1}$, $b = 0.108\text{kHz}$, and $d = 154\text{ms}$. Depending on the tuning of G, the system exhibits a multi-stable regime, with steady states of high and low spiking activity. Here, we set $G = 0.096$.

## 2.4 Co-simulation framework of Arbor and TVB

Both Arbor and TVB offer support for attaching a second simulator to perform co-simulation, potentially at different scales. Co-simulation from TVB's viewpoint is the simpler technology of the two frameworks, since TVB is designed to execute as a single process. TVB allows for exchange of any variable relevant to the region models and any number of variables. The co-simulation partner is encapsulated in one or more TVB regions, called proxy nodes, see Figure 4. These proxies present a conforming interface to TVB, exchanging the salient variables as a table, one row per time-step, one column per variable. As TVB advances in lockstep
<figure><img src="image_6.png" /><figcaption>FIGURE 2
Different neural spiking patterns. (A) Spike Train, $K_{\text{bath}}$ =7.5mM. (B) Tonic Spikes, $K_{\text{bath}}$ =9.5mM. (C) Bursting, $K_{\text{bath}}$ =12.5mM. (D) Seizure-Like Event (SLE), $K_{\text{bath}}$ =17.0mM. (E) Sustained Ictal Activity (SIA), $K_{\text{bath}}$ =17.5mM. (F) Depolarization Block, $K_{\text{bath}}$ =22.5mM. Note that by setting $K_{\text{bath}}$ =4mM, one can obtain Resting State (RS) activity too (result not shown here).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>FIGURE 3
Compartmental neuronal morphologies used for simulations. Each morphology comprises a soma (blue) and a dendritic section (violet). For numerical simulations, dendritic segments are discretized into 5µm compartments, whereas the soma is modeled as a single compartment. Since the cable model neglects extracellular effects, a 1.5-dimensional morphology is employed internally, rendering spatial relationships irrelevant to the model dynamics. When multiple connections converge onto a single cell, synaptic assignments follow a round-robin scheme. (A) Soma only. Equivalent to a point model, the cylindrical segment has a radius *r* and a length of 2*L*, chosen to yield the same surface area as a sphere of radius *r*. Both synapses and spike detectors are attached at the center of the soma. (B) Ball and stick. A straight dendritic segment is added, featuring passive current flow and a synapse attached at a fixed distance from the soma. (C) Random trees. More complex morphologies are generated as random binary trees of depth five. Synapses are placed at a fixed distance from the soma and may be targeted by connections originating from spike detectors, which are positioned at the soma center. See Supplementary Figure S1 for more examples of randomly generated morphologies.</figcaption></figure>

on a global time-step, this is almost identical to normal operation. However, co-simulation introduces the concept of an 'epoch' to TVB, i.e., the length of time that conforms to the smallest delay $\tau_{\min}$ in the set of inter-region connections delays $\tau_{IJ}$, with *I* and *J* referring to two connected regions. These delays are part of

the connectome data used to construct a TVB simulation. In the case that a connectome contains zero-valued delays these must be replaced with a pre-defined finite value. Further, it is required that the time-step evenly divides $\tau_{\min}$. Co-simulation thus can integrate all nodes' state, including the proxy, for one epoch $\tau_{\min}$
without exchanging data. This is correct as an event emanating from any region *I* at time *t* influences any other region *J* at time *t* + τ<sub>IJ</sub> ≥ *t* + τ<sub>min</sub>. Only after an epoch, data need to be exchanged between the proxy and the rest of the TVB regions. A TVB-NEST demonstration has been published to showcase the interaction between a local network of spiking neurons and the whole-brain network dynamics (Kusch et al., 2024).

Arbor has a different design in terms of connectivity. Interaction between physically separate cells is mediated by action potentials, e.g., when the membrane potential triggered by dedicated sources crosses a configurable threshold. Cells are connected by wiring these sources to corresponding sinks like synapses via an abstract connection object comprising a delay and weight, modeling transmission and attenuation via an axon. In contrast to TVB, Arbor is fundamentally a distributed system and internally employs the same approach to decoupling via the minimum network delay as explained above. To initiate co-simulation, Arbor utilizes an additional interface to manage the spike exchange from external connections that are originated from outside but terminate at cells simulated in Arbor. On the technical side, the latter part leverages MPI_Allgather v through an inter-communicator and effects that the concatenation of all spikes sent from all MPI ranks running TVB arrive on all ranks running Arbor and vice versa (Figure 4B). This allows co-simulation in conjunction with arbitrary numbers of ranks on both sides and even in compounds with more than two simulators.

Finally, bi-directional translation between TVB's variable concept and Arbor's representation of action potentials is required. As the former depends on the region models used, we chose to bundle this with the remaining TVB functionality as part of the Arbor proxy node. For the TVB models used in this study, the main variable is the per-region mean activity rate $v_I$ which is conceptually compatible with the concept of spike generation. For each region *I* connected to the proxy node *P*, i.e., with connectome weight $c_{IP} > 0$, a set of synthetic events must be generated such that the mean activity conforms to $v_I$. This is an ambiguous process, even if we prescribe a population (list of cell identifiers) and a per-cell distribution, e.g., a Poisson point process, from which to draw events, which likely must be resolved by ensembles of simulations. In general, this is both model dependent and mathematically intractable, so we leave the general case as a customization point in the framework. For our running example, however, we make the following choice: Event timings for the current step *k* will be drawn from a uniform distribution and dispatched to all cells in the Arbor network. Note that while these events are created at given time a per-connection delay is applied and thus delivery occurs at a later time.

The inverse direction, converting spike events to mean rates, while being well-defined, is still subject to customization. We explore two options here. First, simple running averages, i.e., all spikes that originate within the Arbor network during the current epoch, are collected and sorted into bins of width $\Delta t$. This list is then normalized to the cell count and time step and sent to TVB as the mean activities as a function of time. Although straightforward, this can lead to unrealistically rough activity traces, especially if cell populations are small. Second, as inspired by high-speed calcium imaging experiments (Grewe et al., 2010), a mechanism to track cell

activity through calcium level is implemented as:

$$ \frac{dC_p}{dt} = -\frac{C_p(t)}{\tau} + \beta \sum_{t_{\text{spike},p}} \delta(t - t_{\text{spike},p}), \quad C_p(0) = 0 \quad (16) $$

with per cell *p* having a decay parameter $\tau$ and a weight $\beta$. Computing the activity becomes the average:

$$ v_P(t) = \langle C_p(t) \rangle_{p \in P}, \quad (17) $$

yielding a smooth trace. This method of converting discrete spiking events into a continuous interval variable is also used in a few plasticity models recruiting a negative feedback control mechanism such as synaptic scaling (Van Rossum et al., 2000) and homeostatic structural plasticity (Butz and Van Ooyen, 2013; Diaz-Pier et al., 2016; Lu et al., 2024) models. The choice of the calcium kernel parameters is based on trial-and-error to match the output of Arbor and the magnitude of activity generated by TVB. Figure 5 compares the impact of this choice on the macro-scale network. In small networks and over short timescales defined by the epoch length as shown in the example, spiking activity occurs in noncontinuous bursts, which is dubious in conjunction with the smooth dynamics of the chosen TVB model. Figure 5A shows the propagation of this noncontinuous activity into the TVB regions, while using the Ca-like model (B) provides smooth dynamics in both the Arbor and TVB models. We thus will use the latter in all simulations from here on out. In general, both methods require scaling by the number of cells in the proxy region to arrive at a scale-free activity measure. A local scaling factor $G_A$ is used to convert between the activity of the detailed network and the activity of the region modeled in TVB. In general, $G_A$ needs to be adjusted to the choice of connectome and TVB model, similar to the choice of the global coupling strength $G$ in the RWW model. In this study, $G_A = 100$ is used as it produces seizure-like propagation patterns comparable to those found in similar studies, see e.g., Melozzi et al. (2017) and Courson et al. (2024) and references therein.

## 3 Results

So far, we have described the components of the co-simulation framework. It consists of the following components: a TVB model based on the connectome and node dynamics, a specified set of TVB nodes where the Arbor models will be placed, one or more internally connected network models in Arbor, a defined mechanism for routing events from TVB to individual cells in Arbor, a method for translating Arbor-generated events into TVB variables, and a translation process for converting TVB variables into events originating from synthetic cells. Each of these components serves as a customization point for the user. While reasonable default configurations can be provided for some, others require user-defined specifications to suit specific modeling needs.

The single cell model has been demonstrated to exhibit the necessary range of behaviors. We have also motivated our choice for converting spikes to rates of using a biologically-inspired exponential smoothing filter via a Ca-like activity over
(A)

<figure><img src="image_8.png" /></figure>

(B)

<figure><img src="image_9.png" /><figcaption>FIGURE 4</figcaption></figure>

Arbor-TVB co-simulation schematic and communication pattern. (A) In a TVB simulation of regions *I*, *K*, and *P*, one region *P* will be replaced by a proxy containing a network of detailed cells simulated in Arbor. Regions are connected via the weights of the connectome and produce an activity values based on the chosen region model. When crossing the boundary between TVB and Arbor models, care needs to be taken to convert between discrete action potentials in Arbor to continuous, region-model-specific variables in TVB. (B) Spikes generated by Arbor and TVB—converted from activity values interpreted as mean spiking rates—are exchanged using an MPI intercommunicator and the All-gather primitive. This is equivalent to concatenating all contributions from all Arbor MPI ranks and sending the result to all TVB ranks and vice-versa.

<figure><img src="image_10.png" /><figcaption>FIGURE 5</figcaption></figure>

Impact of conversion method on activity exchange. For an all-to-all connected network of a mixture of 10 SLE cells ($K_{\text{bath}}$ = 17.5mM) and 90 tonic ($K_{\text{bath}}$ = 9.5mM) neurons in Arbor, we plot the membrane potential traces for four tonic and one SLE cell in (A). This simulation is repeated for two activity exchanging methods, either spikes were binned into buckets of width $\Delta t$ to extract instantaneous rates, or the differential equation Equation 16 emulating the change in Calcium concentration of a biological cell after spiking was used (with $\tau$ = 100ms and $\beta$ = $\frac{0.1}{N}$). The resulting activity traces for the Arbor network (B) and selected TVB nodes (out of 98 regions) are displayed in (C–F).

simple binning and fix parameters to $\tau$ = 100ms and $\beta$ = $\frac{0.1}{N}$. This normalization is important, as it produces results that are invariant under changes in the number of cells *N* in the detailed network.

To induce seizure propagation in mice brain, we use the structural connectivity derived from the Allen mouse brain atlas (Oh et al., 2014), also used in Melozzi et al. (2017), to embed the TVB nodes. Following the work presented in Courson et al. (2024),
the proxy node modeling the Arbor population of tonic-spiking and SLE point neurons is set within left the Hippocampus and in particular in the left-field CA1 (l CA1), a region that is prone to generate widespread seizures.

## 3.1 Seizure induction among morphologically detailed cells

As an initial step, we examined a network of morphologically detailed cells (comprising of a single-compartment soma and a random dendritic tree, see bottom-right inset panel of Figure 6) without embedding them into a co-simulation framework. The objective was to assess how interactions among neurons exhibiting distinct firing patterns influence network dynamics within a small population of SLE neurons at the local (Arbor) level. We constructed a network comprising 80 tonic-spiking neurons ($K_{\text{bath}}$ = 9.5mM) and 20 SLE neurons ($K_{\text{bath}}$ = 17.0mM), initially without internal synaptic connectivity. The system was simulated for $T = 5$s, and the resulting membrane potentials are presented in Figure 6 (this initial phase is indicated by the gray-shaded region). During this initial phase, all neurons independently showed their intrinsic firing behavior, consistent with the activity shown in Figure 6. Following this baseline simulation, integration was paused, the network was reconfigured to full connectivity, and the simulation was resumed. The introduction of connectivity produced an immediate response across the network. The SLE neurons progressively transitioned toward tonic-like spiking, mirroring the dominant (80%) tonic-spiking subpopulation. At the same time, tonic-spiking neurons began to exhibit the bursting activity characteristic of the SLE phenotype (Figure 6, $T \ge 5$s). Note that the use of different $K_{\text{bath}}$ values for neurons within the same group is, at this stage, driven by computational considerations and the need to validate the performance of the Arbor-TV B co-simulator. This choice can be adjusted to study more biologically realistic scenarios which involve the interaction of different neural populations or transition of network dynamics, which is further addressed in the Discussion section.

## 3.2 Seizure induction and propagation employing point neurons in Arbor-TV B co-simulator

Next, we embedded a network of detailed neurons as a proxy node in a simulation of neural mass models in TVB. The general setup is similar as before, however, the proxy node now consists of only SLE-type point neurons in a fully connected network. Figure 7A illustrates the evolution of the membrane potential of individual neurons within the Arbor network, with all neurons exhibiting identical dynamical behavior. The pattern of SLE activity is modulated by neuronal coupling.

We run simulations over 20s, and investigate the effect of SLE emergence in the network once all brain areas have reached their steady state. The production of these patterns in the Arbor node generates changes in the firing rates of local TVB brain areas.

Even though changes in activity occur in most brain areas, these fluctuations occur within different ranges.

In Figure 7B, we show the time-series of TVB nodes' firing rate throughout the simulation, here with $N = 100$ in the Arbor detailed network. A transient period is necessary before all nodes reach their steady-state. The emergence of recurrent SLE patterns in the Arbor node triggers a mixture of periodic and seizure-like patterns in the dynamics of the TVB nodes. We highlight traces corresponding to the l CA1 region modeled using detailed cells (orange) and the four regions with the highest activity modeled by the neural mass model. The inset shows the highlighted traces, excluding l CA1, during a single period.

In Figure 8, we present the propagation of SLE originating in l CA1 (star marker), represented as a fully connected network of $N = 10000$ SLE point neurons. Colors on the brain template show the time distance between seizure emergence in the diseased area and seizure arrival in the different brain regions. The systematic detection of SLE onset in other regions of the whole-brain network (if and when present within the simulation period) is performed as follows. First, the baseline firing activity of each region is defined (first time window of the simulation). The non-SLE nodes initially exhibit a relatively stable firing rate, which gradually becomes influenced by the SLE originating from the Arbor node. When connected to the Arbor node, spiking-tonic brain areas are repeatedly recruited into SLE activity patterns, followed by periods of relaxation. For each region, we define the baseline firing rate as the minimum firing rate observed during the relaxation phase immediately preceding any event in the onset region. When a robust increase in firing rate is detected in a region—indicative of a seizure-like event—the corresponding onset timestamp is recorded. The baseline firing rate serves as a reference, and seizure activity is identified based on a sustained increase in firing rate over a defined duration. We define the onset of a seizure as the beginning of the high-amplitude bursting regime, specifically the point at which the firing rate begins to rise consistently after the first small activity peak. To ensure that only significant deviations are classified as seizures, we require the firing rate to exceed the baseline by at least 10% and remain above this threshold for a minimum of 100 ms. Our seizure onset detection follows a similar (though not identical) approach to that used in related studies—for example, Melozzi et al. (2017) and Courson et al. (2024)—where seizure-like or bursting activity in mouse brain models is identified by monitoring a model variable and applying a threshold to detect the onset. In Figure 8, we also depict the activity time-series of four initially non-SLE nodes of the brain network being recruited in the seizure, namely the ventral part of left Lateral Septal Nucleus (l LSV), l CA3, l ENTCl and r ENTCl. In the main panel (mouse brain template), non-colored regions correspond to areas where seizures either did not occur or had relatively weak effects.

The Arbor network was initialized with fixed connection weights ($w = 0.5$). To assess robustness, we also tested weights drawn from a normal distribution ($\mu = 0.5, \sigma = 0.5$), truncated to positive values. Neural activity was averaged over 20 independent realizations (see Supplementary Figures S2, S3). Experiments with inhibitory/excitatory network variants are also currently in progress. Both baseline and seizure-evoked activity differ across brain regions. Due to the symmetrical
<figure><img src="image_11.png" /></figure>

FIGURE 6
Network level effects induced by SLE activity. Simulation of a 100-cell network with 80 tonic-spiking ($K_{bath}$ = 9.5mM) and 20 SLE ($K_{bath}$ = 17.0mM) neurons. Membrane potentials are shown for one SLE (top) and four tonic-spiking cells. Each neuron includes a single-compartment soma and a random dendritic tree. The model was integrated for 5s (shaded region) without internal connections, then switched to a fully connected network, which settled into a new equilibrium dominated by SLE activity. Inset: Example morphology of the first cell.

<figure><img src="image_12.png" /></figure>

FIGURE 7
Multiscale seizure propagation. (A) Membrane potential and raster plot for four neurons in the detailed fully-connected Arbor neural network. Note that the membrane potential time-series are similar across all neurons. (B) Mean firing rate of various TVB local brain areas vs. the Arbor activity. We track a seizure after a transient period, so that all brain areas have reached their baseline activity. The zoomed-in sections show the propagation of the seizure. We highlight and show in detail the four time series with the largest deviation in activity.

<figure><img src="image_13.png" /></figure>
inter-hemispheric connections in the Allen Mouse Brain SC, 1 ENTCI and r ENTCI share the same baseline firing rate. As the SLE pattern emerges in the left hippocampus, 1 ENTCI exhibits higher spiking rates.

## 3.3 Computational performance of the Arbor-TVB co-simulation framework

We next evaluated the performance of the running example on a single Apple M1 (2021) laptop. Arbor was built with MPI and SIMD (Arm Neon/SVE) support, with cells organized into groups of ten to fully exploit SIMD capabilities. The overall runtime consists of four primary components: (1) Arbor model update, (2) Conversion from spikes to rates, (3) TVB model update, and (4) Conversion from rates to spikes.

The Arbor update runs in parallel with the conversions between rates and spikes, as well as the TVB update. During spike exchange, both simulations synchronize, meaning the slower part must wait in the MPI collective, which accounts for the primary time spent in the collective call. Figure 9 illustrates the total runtime of a 10s simulation for the entire model described above, along with the relative contributions, for system sizes ranging from one to 10,000 cells. Notably, at 10,000 cells, nearly all computational time is spent within the Arbor network model. In future experiments, we plan to leverage additional hardware, including GPUs, to accelerate the Arbor side of the simulation. At this scale, TVB and the spike/rate conversions are potential bottlenecks that will require optimization, potentially through TVB's JIT compilation and GPU acceleration. Additionally, further parallelization and porting of the conversion steps to a more performant programming environment remain promising avenues for improvement.

## 4 Discussion

In this work, we presented a co-simulation framework that offers a novel approach to bridging the gap between microscopic (spiking neuron) and macroscopic (mean-field) models. This framework integrates simulators Arbor and TVB within a parallelized MPI environment, enabling a detailed yet computationally feasible representation of neural dynamics across scales. To model large-scale brain network dynamics, we used the mean-field reduced Wong-Wang model to reproduce resting-state dynamics. Simultaneously, a detailed spiking neural network was simulated with Arbor, employing a physiological model of seizures at the neuron level (Figure 2). The spiking activity of the population was then converted into a smooth trace for the proxy node, which was communicated with TVB (Figure 4B). This co-simulation approach successfully captured the interplay between spiking activity and large-scale brain dynamics, where local neuron dynamics generate global activity wave-fronts. At the microscopic scale, we demonstrated that the structure of the detailed neural network influences its activity patterns, thereby affecting the shape of the activity wavefront (Figure 5).

As a proof of concept, and as a technical showcase, we simulated the emergence of seizure-like activity patterns (SLE) in the mouse

hippocampus, using the Allen Mouse Brain Structural Connectivity data. By tuning the single neural parameter, we modeled the target brain area with a small, fully-connected network of (point and detailed-compartmental) neurons SLE, which produces scale-free activity patterns. This network design can be further adapted to investigate more bio-inspired scenarios, such as the study of: (i) ensembles of neurons within a given brain area that belong to different sub-regions (i.e., are relatively distant from one another), whose dynamical activity may differ significantly from that of other sub-groups and (ii) dynamical transition phenomena where the $K_{\text{bath}}$ value is set near a critical threshold—e.g., inducing a transition from regular spiking to bursting activity, see Figures 2A, B. The latter consideration is more general, in the sense that for a different dynamical model of neural activity, one can choose a relevant control parameter that drives transitions in neural activity patterns—for example, from slow spiking (representing healthy activity) to fast bursting (associated with pathological activity) in a given region, when the system operates near a critical transition point (e.g., from spike trains to bursting). Our approach offers a well-understood and easily controlled platform by showcasing its technical usability. By adapting the Arbor nodes with different cell composition and connectivity, we demonstrated its significant potential for more complex cell models untapped.

Although using seizure propagation as a use case, the seizure activity patterns established in the present study may not capture precisely the complex nature of all epileptic seizures. For example, there are other types of seizures that occur in different brain disorders, i.e., acute symptomatic seizures, which are not like those caused in epilepsy. For example, patients with Alzheimer's disease may also experience seizures, which are classified as progressive symptomatic seizures and typically arise from underlying neurodegenerative processes, see e.g., (Mauritz et al., 2022) for a recent relevant review.

The implementation in the present study extended the functionality and application scenarios of Arbor. Arbor has embarked on many types of computational studies as a new-generation simulator. It enables seamless conversion and simulation of single-neuron models from the NEURON simulator and supports simulations of both individual neurons and large-scale networks. Arbor accommodates various plasticity models, including spike-timing-dependent plasticity, calcium-based synaptic tagging and capture, and structural plasticity. It has been used to study synaptic tagging and capture via the built-in diffusion functionality (Luboeinski et al., 2024, under review). Recent developments focus on co-simulation with membrane dynamics and external kernels, enabling dynamic connectivity modifications in a distance-dependent manner. With its high flexibility and scalability, Arbor stands out as a promising platform developed within the EBRAINS initiative to advance cross-scale simulations in computational neuroscience. Arbor is available as part of the EBRAINS software distribution (ESD) on connected HPC centers and the EBRAINS collab via Jupyter Lab. Providing a bridge between morphologically detailed neurons and neural mass models encompassing the full brain spans a gap of scale from sub-micrometer to decimeter. It allows for placing the resolution—using Arbor and detailed models—where needed and using realistic, data driven environments everywhere else via TVB.
<figure><img src="image_14.png" /><figcaption>FIGURE 8</figcaption></figure>

Propagation of a seizure originating in left-field I CA1 area of the mouse brain model. Time distance between seizure emergence in I CA1 (star marker) and spiking rate increase in each brain area. We also depict the firing activity time-series of four initially non-SLE nodes of the brain network being recruited in the seizure, namely the ventral part of left Lateral Septal Nucleus (I Lsv), I CA3, I ENTCI and r ENTCI. The red arrows indicate the beginning of the bursting activity. Non-colored regions (panel with the mouse brain template) correspond to areas where seizures either did not occur or had relatively weak effects. See text for more information regarding the onset detection of seizures in other regions. See text for more details.

Despite its successes, the Arbor-TV B framework has some
limitations. The co-simulation requires careful exploration
and calibration of coupling parameters to ensure meaningful
interactions between Arbor and TVB, which remains a challenge
when generalizing to diverse neural models. The TVB network
we used here is homogeneous, in the sense that all model
parameters are set to be identical. While this is not a highly
realistic assumption—particularly when assigning a dynamical
mean field model to simulate the activity of a brain region—it
is a common choice among researchers when modeling whole-
brain resting-state dynamics, see e.g., Popovych et al. (2021) and
Manos et al. (2023) and references therein. Note that even with
identical initial settings, firing rates vary across regions due to
the influence of long-distance connectivity weight values and the
resulting complex dynamics. Our use-case simulation of seizure
propagation is not yet directly compared to experimental data
from mice or humans. Nevertheless, the current Arbor-TV B
implementation is capable of indirectly capturing biologically
realistic brain dynamics and activity propagation, as reported
for example in Melozzi et al. (2017) and Courson et al. (2024).
The same applies to the choice of calcium kernel parameters.
This kernel is inspired by experimental observations via calcium

imaging experiments (Grewe et al., 2010) but no exact values are
available. Due care should be taken when matching the magnitudes
of both Arbor output and TVB output. Moreover, within the
Arbor-TV B framework, the chosen dynamical model can be
further tuned to generate neural activity—such as BOLD signals
or firing rates—that more closely aligns with neuroimaging data,
thereby enabling the simulation of more realistic brain dynamics.
Specifically, computational costs for large neural networks may
necessitate further optimization in model parallelization and data
handling. A near-term goal would be to incorporate new features,
such as synaptic plasticity, which could offer valuable insights into
how brain networks adapt and reorganize in response to disrupted
activity.

As stated in the Introduction section, it is natural to compare
the Arbor-TV B co-simulation framework presented in the current
Method paper with the established NEST-TV B co-simulation
(Kusch et al., 2024). Both simulators have their own merits
and application scope. NEST has a longer history (Gewaltig
and Diesmann, 2007) and a rich profile of neural models and
plasticity models. Its performance is excellent in mesoscopic
modeling of spiking neural networks. Arbor is younger but has
made it possible to efficiently simulate networks of the highly
<figure><img src="image_15.png" /></figure>

expensive biophysical HH model with multi-compartments (Abi Akar et al., 2019). It also accommodates a variety of plasticity rules, including heterosynaptic dendritic plasticity rules inside the dendritic shaft (Luboeinski et al., 2024). The Arbor-TV B co-simulation framework has thereby inherited those merits of Arbor, while the NEST-TV B co-simulation framework opens the venue for users to freely use NEST functions for co-simulation. Those two tools are complementary for distinct research purposes. Moreover, since Arbor is designed to make the best use of both GPU and CPU, Arbor-TV B can also be easily adapted to make use of the cutting-edge exascale GPU computing resources.

From an epilepsy-seizure perspective, while the framework provides insights into seizure propagation, additional validation against empirical data would enhance its applicability to clinical settings. This co-simulation framework could enable a detailed investigation of the physiological sources of seizures. Understanding the impact of the structure of the diseased area on seizure patterns and propagation would be of great interest (see e.g., Netoff et al., 2004; Garcia-Ramos et al., 2016). Specifically, we expect the inhibition and excitation ratios in the detailed neural network to play a critical role in seizure patterns (see e.g., Engel, 1996; Liu et al., 2020). Moreover, the Arbor-TV B user can implement various neural models and configuration topologies to simulate different brain regions and to computationally investigate diverse dynamic activities or the effects of medical interventions– for example, modeling subthalamic neurons along with synaptic and structural plasticity under stimulation in Parkinson’s disease (Manos et al., 2021; Meier et al., 2022; Shaheen et al., 2022), Alzheimer’s disease (Stefanovski et al., 2019; Manos et al., 2023) or tinnitus (Manos et al., 2018a,b) etc. Evidently, comparison,

parameter tuning, and validation are also feasible using empirical neuroimaging data, however this was not the primary goal of this work. The simulated neural activity generated in each brain region of the connectome can be transformed into a BOLD signal-similar to the built-in functionality of TVB (Sanz Leon et al., 2013; Melozzi et al., 2017), which computes the hemodynamic response function (HRF) kernel (i.e., fMRI activity) using the Balloon-Windkessel model (Friston et al., 2000)—and can ultimately be aligned with neuroimaging time series data. Hence, a framework like Arbor-TV B can be extended to investigate various brain conditions.

## Data availability statement

The codes necessary to reproduce selected results of this work are available on GitHub: https://github.com/arbor-contrib/arbor-tvb-cosim.

## Author contributions

TH: Writing - review & editing, Software, Validation, Resources, Writing - original draft, Formal analysis, Data curation, Methodology, Investigation, Visualization. JC: Investigation, Writing - review & editing, Software, Writing - original draft, Validation, Visualization, Data curation, Methodology, Formal analysis. HL: Writing - review & editing, Investigation, Validation, Writing - original draft, Formal analysis, Methodology, Data curation, Software, Visualization. SD-P: Project administration, Funding acquisition, Validation, Conceptualization, Supervision, Writing - review & editing, Methodology. TM: Validation, Formal
analysis, Methodology, Conceptualization, Supervision, Project administration, Resources, Writing – original draft, Writing – review & editing, Funding acquisition.

## Funding

The author(s) declared that financial support was received for this work and/or its publication. JC was supported by the LABEX MME-DII (ANR-16-IDEX-0008) PhD grant. This work was supported by the NIH 1R21AG087888-01 grant. HL was supported by EBRAINS2.0. EBRAINS 2.0 has received funding from the European Union's Research and Innovation Program Horizon Europe under Grant Agreement No. 101147319. Open access publication funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)—491111487.

## Conflict of interest

The author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The authors TM, HL declared that they were an editorial board member of Frontiers, at the time of submission. This had no impact on the peer review process and the final decision.

## Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

## Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

## Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fncom.2025.1731161/full#supplementary-material

## References

Abi Akar, N., Cumming, B., Karakasis, V., Küsters, A., Klijn, W., Peyser, A., et al. (2019). "Arbor – a morphologically-detailed neural network simulation library for contemporary high-performance computing architectures," in *2019 27th Euromicro International Conference on Parallel, Distributed and Network-Based Processing (PDP)* (New York City, NY: IEEE), 274–282.

Augustine, G. J., Santamaria, F., and Tanaka, K. (2003). Local calcium signaling in neurons. *Neuron* 40, 331–346. doi: 10.1016/S0896-6273(03)00639-1

Barral, J., Wang, X.-J., and Reyes, A. D. (2019). Propagation of temporal and rate signals in cultured multilayer networks. *Nat. Commun.* 10:3969. doi: 10.1038/s41467-019-11851-0

Breyton, M., Fousek, J., Rabuffo, G., Sorrentino, P., Kusch, L., Massimini, M., et al. (2024). Spatiotemporal brain complexity quantifies consciousness outside of perturbation paradigms. *eLife*, 13:RP98920. doi: 10.7554/eLife.98920

Butz, M., and Van Ooyen, A. (2013). A simple rule for dendritic spine and axonal bouton formation can account for cortical reorganization after focal retinal lesions. *PLoS Comput. Biol.* 9:e1003259. doi: 10.1371/journal.pcbi.1003259

Carnevale, N. T., and Hines, M. L. (2006). *The NEURON Book*. Cambridge: Cambridge University Press.

Chater, T. E., Eggl, M. F., Goda, Y., and Tchumatchenko, T. (2024). Competitive processes shape multi-synapse plasticity along dendritic segments. *Nat. Commun.* 15:7572. doi: 10.1038/s41467-024-51919-0

Cook, B. J., Peterson, A. D. H., Woldman, W., and Terry, J. R. (2022). Neural field models: A mathematical overview and unifying framework. *Mathem. Neurosci. Appl.* 2, 1–67. doi: 10.46298/mna.7284

Cooray, G. K., Cooray, V., and Friston, K. J. (2025). Cortical dynamics of neural- connectivity fields. *J. Comput. Neurosci.* 53, 373–391. doi: 10.1007/s10827-025-00903-8

Cooray, G. K., Rosch, R. E., and Friston, K. J. (2023). Global dynamics of neural mass models. *PLoS Comput. Biol.* 19:e1010915. doi: 10.1371/journal.pcbi.1010915

Courson, J., Quoy, M., Timofeeva, Y., and Manos, T. (2024). An exploratory computational analysis in mice brain networks of widespread epileptic seizure onset locations along with potential strategies for effective intervention and propagation control. *Front. Comput. Neurosci.* 18:1360009. doi: 10.3389/fncom.2024.13 60009

Cumming, B., Yates, S., Hater, T., Lu, H., Huisman, B., Wouter, K., et al. (2024). *Arbor v0.10.0*. Zenodo.

Darnell, R. B. (2013). Rna protein interaction in neurons. *Annu. Rev. Neurosci.* 36, 243–270. doi: 10.1146/annurev-neuro-062912-114322

Deco, G., Jirsa, V. K., Robinson, P. A., Breakspear, M., and Friston, K. (2008). The dynamic brain: from spiking neurons to neural masses and cortical fields. *PLoS Comput. Biol.* 4:e1000092. doi: 10.1371/journal.pcbi.100 0092

Deco, G., Ponce-Alvarez, A., Mantini, D., Romani, G. L., Hagmann, P., and Corbetta, M. (2013). Resting-state functional connectivity emerges from structurally and dynamically shaped slow linear fluctuations. *J. Neurosci.* 33, 11239–11252. doi: 10.1523/JNEUROSCI.1091-13.2013

Demirtaş, M., Burt, J. B., Helmer, M., Ji, J. L., Adkinson, B. D., Glasser, M. F., et al. (2019). Hierarchical heterogeneity across human cortex shapes large-scale neural dynamics. *Neuron* 101:1181–1194. doi: 10.1016/j.neuron.2019. 01.017

Depannemaecker, D., Ivanov, A., Lillo, D., Spek, L., Bernard, C., and Jirsa, V. (2022). A unified physiological framework of transitions between seizures, sustained ictal activity and depolarization block at the single neuron level. *J. Comput. Neurosci.* 2022, 1–17. doi: 10.1007/s10827-022-00811-1

Diaz-Pier, S., Naveau, M., Butz-Ostendorf, M., and Morrison, A. (2016). Automatic generation of connectivity for large-scale neuronal network models through structural plasticity. *Front. Neuroanat.* 10:57. doi: 10.3389/fnana.2016.00057

Dicks, L. M. (2022). Gut bacteria and neurotransmitters. *Microorganisms* 10:1838. doi: 10.3390/microorganisms10091838

Emina, F., and Kropff, E. (2022). Selective connectivity enhances storage capacity in attractor models of memory function. *Front. Systems Neurosci.* 104:983147. doi: 10.3389/fnsys.2022.983147

Engel, J. (1996). Excitation and inhibition in epilepsy. *Canadian J. Neurol. Sci.* 23, 167–174. doi: 10.1017/S0317167100038464

Fousek, J., Rabuffo, G., Gudibanda, K., Sheheitli, H., Petkoski, S., and Jirsa, V. (2024). Symmetry breaking organizes the brain's resting state manifold. *Sci. Rep.* 14:31970. doi: 10.1038/s41598-024-83542-w
Friston, K., Mechelli, A., Turner, R., and Price, C. (2000). Nonlinear responses in fmri: The balloon model, volterra kernels, and other hemodynamics. *Neuroimage* 12:466-477. doi: 10.1006/nimg.2000.0630

Garcia-Ramos, C., Lin, J. J., Kellermann, T. S., Bonilha, L., Prabhakaran, V., and Hermann, B. P. (2016). Graph theory and cognition: a complementary avenue for examining neuropsychological status in epilepsy. *Epilepsy Behav.* 64, 329–335. doi: 10.1016/j.yebeh.2016.02.032

Gast, R., Solla, S. A., and Kennedy, A. (2024). Neural heterogeneity controls computations in spiking neural networks. *Proc. Nat. Acad. Sci.* 121:e2311885121. doi: 10.1073/pnas.2311885121

Gewaltig, M.-O., and Diesmann, M. (2007). Nest (neural simulation tool). *Scholarpedia* 2:1430. doi: 10.4249/scholarpedia.1430

Grewe, B. F., Langer, D., Kasper, H., Kampa, B. M., and Helmchen, F. (2010). High-speed *in vivo* calcium imaging reveals neuronal network activity with near-millisecond precision. *Nat. Methods* 7, 399–405. doi: 10.1038/nmeth.1453

Gurkiewicz, M., Korngreen, A., Waxman, S. G., and Lampert, A. (2011). Kinetic modeling of nav1.7 provides insight into erythromelalgia-associated fl449v mutation. *J Neurophysiol.* 105, 1546–1557. doi: 10.1152/jn.00703.2010

Hananeia, N., Ebner, C., Galanis, C., Cuntz, H., Opitz, A., Vlachos, A., et al. (2024). Multi-scale modelling of location-and frequency-dependent synaptic plasticity induced by transcranial magnetic stimulation in the dendrites of pyramidal neurons. *bioRxiv*. doi: 10.1101/2024.07.03.601851

Hashemi, M., Depannemaeker, D., Saggio, M., Triebkorn, P., Rabuffo, G., Fousek, J., et al. (2025). Principles and operation of virtual brain twins. *IEEE Rev Biomed Eng.* 1–29. doi: 10.1109/RBME.2025.3562951

Herculano-Houzel, S. (2009). The human brain in numbers: a linearly scaled-up primate brain. *Front. Hum. Neurosci.* 3:857. doi: 10.3389/neuro.09.031.2009

Herculano-Houzel, S. (2012). The remarkable, yet not extraordinary, human brain as a scaled-up primate brain and its associated cost. *Proc. Nat. Acad. Sci.* 109, 10661–10668. doi: 10.1073/pnas.1201895109

Hines, G. S. (1984). Advances in simulation software technology. *Simulation* 42, 36–36. doi: 10.1177/003754978404200105

Hodgkin, A. L., and Huxley, A. F. (1952a). The components of membrane conductance in the giant axon of ligo!o. *J. Physiol.* 116:473. doi: 10.1113/jphysiol.1952.sp004718

Hodgkin, A. L., and Huxley, A. F. (1952b). Currents carried by sodium and potassium ions through the membrane of the giant axon of ligo!o. *J. Physiol.* 116:449. doi: 10.1113/jphysiol.1952.sp004717

Hodgkin, A. L., and Huxley, A. F. (1952c). A quantitative description of membrane current and its application to conduction and excitation in nerve. *J. Physiol.* 117:500. doi: 10.1113/jphysiol.1952.sp004764

Hodgkin, A. L., Huxley, A. F., and Katz, B. (1952). Measurement of current-voltage relations in the membrane of the giant axon of ligo!o. *J. Physiol.* 116:424. doi: 10.1113/jphysiol.1952.sp004716

Kusch, L., Diaz-Pier, S., Klijn, W., Sontheimer, K., Bernard, C., Morrison, A., et al. (2024). Multiscale co-simulation design pattern for neuroscience applications. *Front. Neuroinform.* 18:1156683. doi: 10.3389/fninf.2024.1156683

Lavanga, M., Stumme, J., Yalcinkaya, B. H., Fousek, J., Jockwitz, C., Sheheitli, H., et al. (2023). The virtual aging brain: Causal inference supports interhemispheric dedifferentiation in healthy aging. *Neuroimage* 283:120403. doi: 10.1016/j.neuroimage.2023.120403

Linssen, C., Babu, P. N., Bouhadjar, Y., Ewert, L., Wybo, W., Lober, M., et al. (2024). *Nestml* 8.0.0. Zenodo.

Liu, Y., Grigorovsky, V., and Bardakjian, B. (2020). Excitation and inhibition balance underlying epileptiform activity. *IEEE Trans. Biomed. Eng.* 67, 2473–2481. doi: 10.1109/TBME.2019.2963430

Loligo, A. (1952). The dual effect of membrane potential on sodium conductance in the giant. *J. Physiol.* 6, 497–506. doi: 10.1113/jphysiol.1952.sp004719

Lu, H., Diaz, S., Lenz, M., and Vlachos, A. (2024). The interplay between homeostatic synaptic scaling and homeostatic structural plasticity maintains the robust firing rate of neural networks. *eLife*. doi: 10.7554/eLife.88376.2

Lu, H., Garg, S., Lenz, M., and Vlachos, A. (2025). Repetitive magnetic stimulation with itbs600 induces persistent structural and functional plasticity in mouse organotypic slice cultures. *bioRxiv*. doi: 10.2139/ssrn.5203072

Luboeinski, J., Schmitt, S., Shafiee, S., Hater, T., Bösch, F., and Tetzlaff, C. (2024). Plastic arbor: a modern simulation framework for synaptic plasticity-from single synapses to networks of morphological neurons. *arXiv* [preprint] *arXiv:2411.16445*. doi: 10.48550/arXiv.2411.16445

Luboeinski, J., and Tetzlaff, C. (2021). Memory consolidation and improvement by synaptic tagging and capture in recurrent neural networks. *Commun. Biol.* 4:275. doi: 10.1038/s42003-021-01778-y

Manos, T., Diaz-Pier, S., Fortel, I., Driscoll, I., Zhan, L., and Leow, A. (2023). Enhanced simulations of whole-brain dynamics using hybrid resting-state structural

connectomes. *Front. Comput. Neurosci.* 17:1295395. doi: 10.3389/fncom.2023.12.9359

Manos, T., Diaz-Pier, S., and Tass, P. A. (2021). Long-term desynchronization by coordinated reset stimulation in a neural network model with synaptic and structural plasticity. *Front. Physiol.* 12:716556. doi: 10.3389/fphys.2021.716556

Manos, T., Zeitler, M., and Tass, P. A. (2018a). How stimulation frequency and intensity impact on the long-lasting effects of coordinated reset stimulation. *PLoS Comput. Biol.* 14(5). doi: 10.1371/journal.pcbi.1006113

Manos, T., Zeitler, M., and Tass, P. A. (2018b). Short-term dosage regimen for stimulation-induced long-lasting desynchronization. *Front. Physiol.* 9:376. doi: 10.3389/fphys.2018.00376

Mauritz, M., Hirsch, L. J., Camfield, P., Chin, R., Nardone, R., Lattanzi, S., et al. (2022). Acute symptomatic seizures: an educational, evidence-based review. *Epileptic Disord.* 24, 26–49. doi: 10.1684/epd.2021.1376

Meier, J. M., Perdikis, D., Blickensdörfer, A., Stefanovski, L., Liu, Q., Maith, O., et al. (2022). Virtual deep brain stimulation: Multiscale co-simulation of a spiking basal ganglia model and a whole-brain mean-field model with the virtual brain. *Exp. Neurol.* 354:114111. doi: 10.1016/j.expneurol.2022.114111

Melozzi, F., Woodman, M., Jirsa, V., and Bernard, C. (2017). The virtual mouse brain: A computational neuroinformatics platform to study whole mouse brain dynamics. *eNeuro*. 4. doi: 10.1523/ENEURO.0111-17.2017

Miner, D., and Triesch, J. (2016). Plasticity-driven self-organization under topological constraints accounts for non-random features of cortical synaptic wiring. *PLoS Comput. Biol.* 12:e1004759. doi: 10.1371/journal.pcbi.1004759

Monteverdi, A., Palesi, F., Schirner, M., Argentino, F., Merante, M., Redolfi, A., et al. (2023). Virtual brain simulations reveal network-specific parameters in neurodegenerative dementias. *Front. Aging Neurosci.* 15:1204134. doi: 10.3389/fnagi.2023.1204134

Nayebi, A., Attinger, A., Campbell, M., Hardcastle, K., Low, I., Mallory, C. S., et al. (2021). "Explaining heterogeneity in medial entorhinal cortex with task-driven neural networks," in *Advances in Neural Information Processing Systems*, eds. M. Ranzato, A. Beygelzimer, Y. Dauphin, P. Liang, and J. W. Vaughan (Red Hook, NY: Curran Associates Inc), 12167–12179.

Netoff, T. I., Clewley, R., Arno, S., Keck, T., and White, J. A. (2004). Epilepsy in small-world networks. *J. Neurosci.* 24, 8075–8083. doi: 10.1523/JNEUROSCI.1509-04.2004

Oh, S. W., Harris, J. A., Ng, L., Winslow, B., Cain, N., Mihalas, S., et al. (2014). A mesoscale connectome of the mouse brain. *Nature* 508, 207–214. doi: 10.1038/nature13186

Poldrack, R. A., and Gorgolewski, K. J. (2017). Openfmri: Open sharing of task fMRI data. *Neuroimage* 144, 259–261. doi: 10.1016/j.neuroimage.2015.05.073

Popovych, O. V., Jung, K., Manos, T., Diaz-Pier, S., Hoffstaedter, F., Schreiber, J., et al. (2021). Inter-subject and inter-parcellation variability of resting-state whole-brain dynamical modeling. *Neuroimage* 236:118201. doi: 10.1016/j.neuroimage.2021.118201

Popovych, O. V., Manos, T., Hoffstaedter, F., and Eickhoff, S. B. (2019). What can computational models contribute to neuroimaging data analytics? *Front. Syst. Neurosci.* 12:68. doi: 10.3389/fnsys.2018.00068

Rabuffo, G., Bandyopadhyay, A., Calabrese, C., Gudibanda, K., Depannemaeker, D., Takarabe, L. M., et al. (2025). Biophysically inspired mean-field model of neuronal populations driven by ion exchange mechanisms. *eLife*. 14:RP104249. doi: 10.7554/eLife.104249.1

Reyes, A. D. (2003). Synchrony-dependent propagation of firing rate in iteratively constructed networks in vitro. *Nat. Neurosci.* 6, 593–599. doi: 10.1038/nn1056

Ritter, P., Schirner, M., McIntosh, A. R., and Jirsa, V. K. (2013). The virtual brain integrates computational modeling and multimodal neuroimaging. *Brain Connect.* 3, 121–145. doi: 10.1089/brain.2012.0120

Sacramento, J., Wichert, A., and van Rossum, M. C. (2015). Energy efficient sparse connectivity from imbalanced synaptic plasticity rules. *PLoS Comput. Biol.* 11:e1004265. doi: 10.1371/journal.pcbi.1004265

Sanz Leon, P., Knock, S. A., Woodman, M. M., Domide, L., Mersmann, J., McIntosh, A. R., et al. (2013). The virtual brain: a simulator of primate brain network dynamics. *Front. Neuroinform.* 7:10. doi: 10.3389/fninf.2013.00010

Sanz-Leon, P., Knock, S. A., Spiegler, A., and Jirsa, V. K. (2015). Mathematical framework for large-scale brain network modeling in the virtual brain. *Neuroimage* 111:385–430. doi: 10.1016/j.neuroimage.2015.01.002

Schirner, M., Domide, L., Perdikis, D., Triebkorn, P., Stefanovski, L., Pai, R., et al. (2022). Brain simulation as a cloud service: the virtual brain on ebrains. *Neuroimage* 251:118973. doi: 10.1016/j.neuroimage.2022.118973

Scott, A. C. (1975). The electrophysics of a nerve fiber. *Rev. Mod. Phys.* 47:487. doi: 10.1103/RevModPhys.47.487

Shaheen, H., Pal, S., and Melnik, R. (2022). Multiscale co-simulation of deep brain stimulation with brain networks in neurodegenerative disorders. *Brain Multiphysics* 3:100058. doi: 10.1016/j.brain.2022.100058
Stefanovski, L., Meier, J. M., Pai, R. K., Triebkorn, P., Lett, T., Martin, L., et al. (2021). Bridging scales in alzheimer's disease: biological framework for brain simulation with the virtual brain. *Front. Neuroinform.* 15:630172. doi: 10.3389/fninf.2021.630172

Stefanovski, L., Triebkorn, P., Spiegler, A., Diaz-Cortes, M.-A., Solodkin, A., Jirsa, V., et al. (2019). Linking molecular pathways and large-scale computational modeling to assess candidate disease mechanisms and pharmacodynamics in alzheimer's disease. *Front. Comp. Neurosci.* 13:2019. doi: 10.3389/fncom.2019.00054

Stimberg, M., Brette, R., and Goodman, D. F. (2019). Brian 2, an intuitive and efficient neural simulator. *eLife* 8:e47314. doi: 10.7554/eLife.47314

Suma, A., Sigg, D., Gallagher, S., Gonnella, G., and Carnevale, V. (2024). Ion channels in critical membranes: clustering, cooperativity, and memory effects. *PRX Life* 2:013007. doi: 10.1103/PRXLife.2.013007

Tetzlaff, S. K., Reyhan, E., Layer, N., Bengtson, C. P., Heuer, A., Schroers, J., et al. (2025). Characterizing and targeting glioblastoma neuron-tumor networks with retrograde tracing. *Cell* 188, 390–411. doi: 10.1016/j.cell.2024.11.002

Thompson, W., and Kelvin, L. (1855). On the theory of the electric telegraph. *Proc. Royal Soc.* 7, 382–399. doi: 10.1098/rspl.1854.0093

Tigerholm, J., Petersson, M. E., Obreja, O., Eberhardt, E., Namer, B., Weidner, C., et al. (2015). C-fiber recovery cycle supernormality depends on ion concentration and ion channel permeability. *Biophys. J.* 108, 1057–1071. doi: 10.1016/j.bpj.2014.12.034

Tigerholm, J., Petersson, M. E., Obreja, O., Lampert, A., Carr, R., Schmelz, M., et al. (2014). Modeling activity-dependent changes of axonal spike conduction in

primary afferent c-nociceptors. *J. Neurophysiol.* 111, 1721–1735. doi: 10.1152/jn.0077.2012

Van Rossum, M. C., Bi, G. Q., and Turrigiano, G. G. (2000). Stable hebbian learning from spike timing-dependent plasticity. *J. Neurosci.* 20, 8812–8821. doi: 10.1523/JNEUROSCI.20-23-08812.2000

Von Helmholtz, H. (1850). Messungen über den zeitlichen Verlauf der Zuckung animalischer Muskeln und die Fortpflanzungsgeschwindigkeit der Reizung in den Nerven [Measurements of the time course of contraction in animal muscles and the propagation speed of excitation in the nerves]. *Archiv für Anatomie, Physiologie und wissenschaftliche Medicin* 17, 176–364. German.

Wang, H. E., Dollomaja, B., Triebkorn, P., Duma, G. M., Williamson, A., Makhalova, J., et al. (2025). Virtual brain twins for stimulation in epilepsy. *Nat. Comput. Sci.* 5, 754–768. doi: 10.1038/s43588-025-00841-6

Wang, X.-J. (2002). Probabilistic decision making by slow reverberation in cortical circuits. *Neuron* 36, 955–968. doi: 10.1016/S0896-6273(02)01092-9

West, A. E., Griffith, E. C., and Greenberg, M. E. (2002). Regulation of transcription factors by neuronal activity. *Nat. Rev. Neurosci.* 3, 921–931. doi: 10.1038/nrn987

Yasumatsu, N., Matsuzaki, M., Miyazaki, T., Noguchi, J., and Kasai, H. (2008). Principles of long-term dynamics of dendritic spines. *J. Neurosci.* 28, 13592–13608. doi: 10.1523/JNEUROSCI.0603-08.2008

Zheng, P., Dimitrakakis, C., and Triesch, J. (2013). Network self-organization explains the statistics and dynamics of synaptic connection strengths in cortex. *PLoS Comput. Biol.* 9:e1002848. doi: 10.1371/journal.pcbi.1002848