Exploring Neocortical Micro- and Mesocircuitry: Physiology & Experimentation

Neuroscience aims to characterize the dynamics of the brain and understand how they emerge as a product of anatomy and physiology. This is challenging, however, due to the complexity of the brain’s multiscale organization and its study in situ. Large-scale, data-driven, biophysically detailed models (Markram et al., 2015; Billeh et al., 2020) offer a tool that can complement laboratory investigation. This ‘bottom-up’ approach combines detailed constituent models into a 3D model of a given brain volume, thus consolidating disparate, multiscale data sources. Such models aim to provide a general tool for multiscale investigation; access and manipulation of any model component enable predictions about how these components combine to shape emergent dynamics. With increasing biological data, model refinement and validation must be continuous. This necessitates (1) open sourcing of models with high-quality software tools for iterative, community-driven refinement and validation, and (2) rigorous methodology for model parameterization and validation, such that model iterations can be meaningfully compared. Here and in our companion paper (Reimann et al., 2024), we present our methodology for parameterizing the anatomy and physiology of a cortical model, and validating its emergent dynamics. This enabled predictions about how cortical activity is shaped by high-dimensional connectivity motifs in local and mid-range connectivity, and spatial targeting rules by inhibitory subpopulations. Elsewhere, this has allowed the presented model to be used to study the formation of cell assemblies (Ecker et al., 2024b), functional synaptic plasticity (Ecker et al., 2024a), propagation of activity between cortical areas (Bolaños-Puchet and Reimann, 2024), the role of non-random connectivity motifs on network activity (Pokorny et al., 2025) and reliability (Egas Santander et al., 2025), the composition of high-level electrical signals such as the EEG (Tharayil et al., 2025), and how spike sorting biases population codes (Laquitaine et al., 2024). Specifically, we built and validated a model of the entire non-barrel primary somatosensory cortex (nbS1) comprising eight subregions. Whilst our previous data-driven, biophysically detailed model (Markram et al., 2015) provided insights at the scale of a single cortical column (Reimann et al., 2013; Reimann et al., 2017; Reimann et al., 2022; Nolte et al., 2019; Nolte et al., 2020; Newton et al., 2021), the new model is ~140 times larger and to our knowledge offers the first simulations of in vivo-like spontaneous and stimulus-evoked activity in a biophysically detailed cortical model with interregion connectivity. In the companion paper, we introduce the anatomical model (Figure 1, step 1), describing how neuron morphologies were placed within an atlas-based geometry and connected through local and mid-range synapses. Here, we describe our improved techniques to model and validate the electrical properties of neurons and synapses (Figure 1, steps 2 and 3; Reva et al., 2023; Barros-Zulaica et al., 2019), and to compensate for input from missing brain areas (Figure 1, step 4). These improvements enabled enhanced validation of emerging in vivo-like activity (Figure 1, step 5), including the reproduction and extension of five published studies in rodent sensory cortex under a single in vivo-like regime (Figure 1, step 6). Figure 1 Validating that network activity emerges from the same interactions driving in vivo dynamics, requires (1) an approach to tackle the large parameter space without overfitting and (2) comparison of emerging dynamics with laboratory experiments. With respect to (1), we strongly prefer parameterization with directly measured quantities over fitting parameters to yield the correct emerging activity. Additionally, where fitting is applied, we adhere to a principle of compartmentalization of parameters. That is, once a parameter has been parameterized at one biological level, it is no longer a free parameter at a higher level. For example, the maximal conductance of a synapse is fit to biological amplitudes of postsynaptic potentials, but is then never updated in the process of reaching in vivo-like population activity. If a connection is valid on the single-cell level, its contribution at the population level should be equally valid. As a result, the emerging in vivo-like activity is the consequence of only 10 free parameters representing the strength of extrinsic input from other brain regions into 9 layer-specific excitatory and inhibitory populations, and a parameter controlling the noise structure of this extrinsic input. These parameters were fitted using a novel methodology presented here and were in accordance with the mean number of missing synapses for each population (Felleman and Van Essen, 1991; Harris et al., 2019; Gao et al., 2022). With respect to (2), we conducted the following validations: spontaneous activity reproduced layer-wise in vivo firing rates (or a specified proportion of in vivo firing rates to account for in vivo recording bias; Wohrer et al., 2013), varied along a spectrum of asynchronous to synchronous activity, exhibited spatially structured fluctuations, and produced long-tailed firing rate distributions with sub 1 Hz peaks as in vivo (Wohrer et al., 2013; Buzsáki and Mizuseki, 2014). Under the same parameterization, the model reproduced precise millisecond dynamics of layer-wise populations in response to simple stimuli. These initial validations demonstrated that the model was in a more accurate regime compared to Markram et al., 2015—an essential step before testing more complex or larger-scale validations. For example, under the same parameterization we then observed selective propagation of stimulus-evoked activity to downstream areas, and reproduced and extended more complex experiments through accurate modeling of targeted optogenetic stimulation and lesions. Importantly, we highlight where emergent activity shows discrepancies with in vivo activity, to guide future data-driven model refinement. The model generated a number of predictions (Supplementary file 1), including about the role of different layers in driving layer 2/3 stimulus-responses and how inhibitory interneuron types encode contrast, synchronous, and rate-coded information. Additionally, with access to the full structural connectome and tools for precisely editing it (Pokorny et al., 2025), we were able to make predictions about the relationship between structure and function. For example, we predict that an increase in the prevalence of non-random connectivity motifs towards deeper layers leads to a matched increase in spiking correlations, and that subregions more strongly innervated by mid-range connectivity have higher correlated activity locally. Additionally, we generated a new connectome that captured recently characterized spatially-specific targeting rules for different inhibitory neuron types (Schneider-Mizell et al., 2025) in the MICrONS electron microscopy dataset (MICrONS Consortium, 2025), such as increased perisomatic targeting by PV+ neurons, and increased targeting of inhibitory populations by VIP+ neurons. Comparing activity to the original connectome gave predictions about the role of these additional targeting rules. For example, inhibitory populations were more strongly inhibited (increasingly towards central layers) and required more non-local drive to reach the firing rates observed in vivo. Evoked responses increased and decreased in superficial and deeper excitatory populations respectively, suggesting layer-specific roles of the more specific inhibitory targeting. To provide a framework for further studies and integration of experimental data, the full model is made available with simulation tools, as well as a smaller subvolume with the optional new connectome capturing inhibitory targeting rules from electron microscopy (Figure 1, step 7). The detailed modeling approach provides a one-to-one correspondence with most types of experimental data, allowing different datasets to be readily integrated. Due to the incredible speed of discovery in neuroscience, an integrative model will always be lagging behind the latest available data. We believe the solution is to provide a scientifically solid, validated model with clearly characterized strengths and weaknesses, along with the tools to advance or customize it for individual projects. We have therefore also made our tools for building and improving the model openly available: https://www.github.com/BlueBrain. We presented novel methods to build, simulate, and validate models of cortical tissue that correspond directly with 3D digitized brain atlases. As demonstrated, this enables laboratory experiments to have a simulatable in silico counterpart, and vice versa: predictions made by the model automatically have a precise correspondence in biology. By recreating and extending five laboratory experiments under a single model parameterization, we provided strong model validation and demonstrated the model’s natural versatility. To our knowledge the simulations of the full nbS1 represent the first simulations of stable in vivo-like spontaneous and stimulus-evoked activity in a large-scale biophysically detailed model of multiple cortical subregions connected through local and mid-range connectivity. The model generated an initial set of predictions about the relationship between cortical structure and function, and provides a basis for future studies and predictions. The table in Supplementary file 1 describes each quantitative prediction made in this article, and how it might be tested in vivo. Additionally, all assumptions are listed in the table in Supplementary file 8. We made the model and simulation tools available for use, enabling community-driven testing, exploration and refinement. Excitation from non-modeled brain regions is modeled by injecting somatic conductances and is calibrated for specific populations. For the purpose of calibration, we presented a novel technique, which allows rapid and methodical characterization and validation of emergent model activity. The efficiency of the technique allows rapid recalibration after changes to the model such as adding new brain regions or changing anatomical details. Connected divided by unconnected firing rates measured during the process are also a prediction of the degree to which the activity of a population is determined by incoming connections from local vs. extrinsic sources. We found an increase of the importance of local connectivity from superficial to deeper layers, in line with the canonical view of the cortex (Lübke and Feldmeyer, 2007; Feldmeyer, 2012), which places layers 5 and 6 at the end of local information flow. Although abstract relative to the model, the injection amplitudes could be predicted from the number of missing synapses, suggesting it provides a functional representation of their anatomical counterparts. The modular model of extrinsic connectivity can be refined or replaced in the future or when larger fractions of the brain are modeled. The calibration framework could optimize per population parameters for other compensation methods, whilst still offering an interpretable spectrum of firing rate regimes at different levels of PFR. For example, more realistic compensation schemes could be explored which introduce (a) correlations between the inputs received by different neurons and (b) compensation distributed across dendrites, as well as at the soma. We predict that such changes would make spontaneous activity more correlated at the lower spontaneous firing rates which supported in vivo-like responses (PFR:0.1−0.5), which would in turn make stimulus-responses more noise correlated. Our approach to obtaining in vivo-like activity contrasts with a recent approach in a large-scale biophysical model of primary visual cortex of similar scale to our smaller subvolume (Billeh et al., 2020). There, extrinsic inputs were delivered through dendritic synapses instead of somatic injection, which is more anatomically realistic, and intrinsic recurrent weights were adjusted to match extracellularly recorded firing rates. They find a single activity state, instead of the in vivo-compatible spectrum here that allowed us to contrast anesthetized and awake states. The model in general differs from hybrid models, which jointly use biologically detailed models and point neuron models (Billeh et al., 2020; Dura-Bernal et al., 2023), and also from the work of Egger et al., 2020, who modeled activation of L5 pyramidal cells by constraining patterns of synaptic input based on receptive fields and predicted anatomical innervation. Given the complexity of the model, how can we be confident that its activity is in vivo-like? First, we demonstrated correlated spontaneous dynamics in the form of global and local fluctuations and dynamic E-I balance, as found experimentally (Renart et al., 2010). Second, population firing rate distributions were long-tailed with sub 1 Hz peaks, and were similar for spontaneous and stimulus-evoked activity, as reported (Wohrer et al., 2013). Moreover, mean firing rates below 1 Hz are required by metabolic constraints (Attwell and Laughlin, 2001; Lennie, 2003). In response to simple whisker stimuli, response sparsity, spike counts, and the temporal profile and amplitudes of layer-wise populations were similar to in vivo under anesthetized and awake conditions, for different regions of the meta-parameter space. The effect of the spontaneous meta-parameters on evoked responses, also shows that the networks spontaneous state affects stimulus-responses, as observed in vivo (Isbister et al., 2021). Although from the barrel system, these stimuli offer some of the simplest stimulus-response paradigms, enabling principles of neural dynamics and information processing to be studied; in particular, the correspondence between atomic units of sensory information and neural representations. The model also predicts the number of thalamic fibers stimulated for whisker-flick stimuli. Such validations are crucial for complex models, as they provide context for more complex validations. Moreover, in a rich and complex nonlinear system, the source of discrepancies with in vivo activity are simpler to ascertain under simpler protocols. For more complex validations, we reproduced protocols from three laboratory experiments. These experiments used a wide array of techniques, which could be accurately recreated due to the model’s close correspondence with biological tissue. Additionally, we could go beyond the original experiments. For example, we predicted an increased number of neurons would be required to have a decodable effect on postsynaptic activity under in vivo conditions over the in vitro conditions of the original experiment. This can be explained by lower synaptic reliability under in vivo conditions. Also, we predicted a paradoxical effect of optogenetic stimulation on L6 PV+ interneurons, namely a decrease in firing with increased stimulus strength. This is reminiscent of the paradoxical responses found by Mahrach et al., 2020 in the mouse anterior lateral motor cortex (in L5, but not in L2/3) and barrel cortex (no layer distinction), respectively. While Mahrach et al., 2020 conducted their recordings in awake mice not engaged in any behavior, we observed this effect only when drifting grating patterns with high contrast were presented. Nevertheless, consistent with their findings, we found the effect only in deep but not in superficial layers, and only for PV+ interneurons but not for PCs. Our model could therefore be used to improve the understanding of this paradoxical effect in follow-up studies. These examples demonstrate that the approach of modeling entire brain regions can be used to further probe the topics of the original articles and cortical processing. This could be continued in the future using paradigms that take advantage of the multiple regions offered by the model to study processing in cortical hierarchies. Responses to whisker flick stimuli were closer to biology for values of the PFR meta-parameter between 0.1 and 0.5, according with the known inflation of mean firing rates stemming from the bias of extracellular spike sorting techniques towards larger and more active neurons (Olshausen and Field, 2006). Supporting our estimates, patch-clamp experiments show ubiquitous neuronal silence (Crochet and Petersen, 2006), and spontaneous and stimulus-evoked barrel cortex firing rates as much as 10 times lower than for extracellular recordings (0.05–0.15 Hz vs. 0.8–1.5 Hz; Olshausen and Field, 2006). This is also in line with our recent work using the model, which estimated a spike sorting bias corresponding to PFR=0.3 using virtual extracellular electrodes (Laquitaine et al., 2024). Each recording technique comes with challenges and biases (Shoham et al., 2006; Olshausen and Field, 2006; Barth and Poulet, 2012; Wohrer et al., 2013). Calcium imaging and patch-clamp experiments are less biased towards frequently spiking neurons (Wohrer et al., 2013), yet the former can only infer spiking activity and may be biased towards neurons that favor marker expression. Patch-clamp techniques can reduce firing rates, be biased towards larger or more active neurons (Olshausen and Field, 2006) and are limited to recording a few neurons simultaneously, preventing characterizations of single trial population dynamics. An earlier version of the model has been used to simulate cortical extracellular potentials (Reimann et al., 2013). We can now simulate modern extracellular electrodes and compare multi-unit activity and spike sorting results to the in vivo-like ground truth to gain further insights into the potential over-estimation of mean firing rates. During simulations of the entire nbS1, we found emerging spatial inhomogeneities in network activity. Specifically, we observed sharp transitions in firing rates and pairwise correlations at disgranular zones. This highlights that the model supports several, independently acting functional units. Analysis of spatial correlations suggests a radius of approximately 400 μm for a single unit. This property was not present in the previous intra-regional scale model (Markram et al., 2015), but is rather an emergent property of the interregional scale, demonstrating that we achieved our goal of building a model that can be used to study mid-range interactions between regions. Analyzing the source of the inhomogeneities suggested a prominent role of a non-random rich-club structure of mid-range connections at the meso-scale. Within functional units, we found that inhomogeneous increases of spiking correlations were confined to spontaneous packets of L6 activity. In keeping with the philosophy of compartmentalization of parameters and continuous model refinement (see Introduction), it was essential to improve validity at the columnar scale (relative to Markram et al., 2015) as part of demonstrating validity of the full nbS1. Indeed, improved parameterization and validation at smaller scales was essential to parameterizing background input which generated robust nbS1 activity within realistic [Ca2+]o and firing rate ranges. We view this as a major achievement, as it was unknown whether the model would achieve a stable and meaningful regime at the start of our investigation. Whilst we would have liked to go further, our primary goal was to publish a well characterized model as an open resource that others could use to undertake further in-depth studies. In this regard, we are pleased that the parameterization of the nbS1 model has already been used to study EEG signals (Tharayil et al., 2025), as well as propagation of activity between two subregions (Bolaños-Puchet and Reimann, 2024). Investigation, improvement and validation must be continued at all spatial scales in follow-up papers with detailed description, figures and analysis, which cannot be covered in this article. Each new study increases the scope and validity of future investigations. In this way, this model and paper act as a stepping stone towards more complex questions of interest to the community such as perception, task performance, predictive coding and dendritic processing. This was similar for Markram et al., 2015 where the initial paper was followed by more detailed studies. Unlike the Markram et al., 2015 model, the new model can also be exploited by the community and has already been used in a number of follow-up papers studying (Ecker et al., 2024a; ; Bolaños-Puchet and Reimann, 2024; Pokorny et al., 2025; Egas Santander et al., 2025; Tharayil et al., 2025; Laquitaine et al., 2024). We believe that the number of use cases for such a general model is vast and is made larger by the increased size of the model. It is important to highlight where the model deviates from biological data. First, the model omits several known anatomical elements, such as glia and the presence of gap junctions between certain neuron-types, but the model’s spatial context provides a natural scaffold for them. Our modeling techniques assume generality of rules and parameters, unless indicated by the data. Therefore, violations indicate that biological rules and data may be more specific. For example, model thalamic inputs innervate any dendrite placed in the targeted layers equally, but this did not reproduce the timing of inhibitory subpopulation responses. Stronger innervation of PV+ neurons than other inhibitory neurons may rectify this. However, it is unclear whether innervation should be anatomically or physiologically stronger, as the best available data mixes structural and physiological strength (Sermet et al., 2019). In our previous work (Roussel et al., 2023), we linked mouse inhibitory me-models to transcriptomic types (t-types) in a whole mouse cortex transcriptomic dataset (Gouwens et al., 2019). This can provide a direct correspondence in future large-scale mouse models. As we model only a single electrical type for pyramidal cells, there is no one-to-one correspondence between our me-models and the 10 different pyramidal cell types identified there. We are not currently aware of any method which can recreate the electrical features of different types of pyramidal cells using only generic ion channel models. To achieve the firing pattern behavior of more specific electrical types, usually ion channel kinetics are tweaked, and this would violate the compartmentalization of parameters. In future we hope to build morpho-electric-transcriptomic type (met-type) models by selecting gene-specific ion channel models (Ranjan et al., 2019; Ranjan et al., 2024) based on the met-type’s gene expression. Data specific to different neuron sections (i.e., soma, AIS, apical/basal dendrites) of different met-types, such as gene expression, distribution of ion channels, and voltage recordings under standard single cell protocols, would be particularly useful. Finally, when required, we generalized data from different animals and brain regions to build our rat nbS1 model. This is the accepted state-of-the-art in computational neuroscience, for example, all 19 data-driven models of rodent microcircuitry listed in Figure 2 of the recent review of Ramaswamy, 2024 conduct some sort of mixing, including the advanced mouse V1 model of the Allen Institute (Billeh et al., 2020). Whilst a first truly single species model would be a great advance, it is not required to be the basis of in silico research. Our model (and previous ones) can be used to study the many cortical mechanisms that are common to the closely related organisms rat and mouse. This paper and its companion paper serve to present a methodology for modeling micro- and mesoscale anatomy and physiology, which can be applied for other cortical regions and species. With the rapid increase in openly available data, efforts are already in progress to build models of mouse brain regions with reduced reliance on data mixing thanks to much larger quantities of available atlas-based data. This also includes data for the validation of emergent network-level activity. Here we chose to compare network-level activity to data mostly from the barrel cortex, as well as a single study from primary visual cortex. Whilst a lot of the data used to build the model was from the barrel cortex, the barrel cortex also represents a very well characterized model of cortical processing for simple and controlled sensory stimuli. The initial comparison of population-wise responses in response to accurate thalamic input for single whisker deflections was essential to demonstrating that the model was closer to in vivo, and we were unaware of similar data for non-barrel somatosensory regions. Moreover, our optogenetic and lesion study demonstrated the capacity to compare and extend studies of canonical cortical processing in the whisker system. However, we cannot always predict the effect caused by mixing data sources. An interesting example of this is described in the companion paper (Reimann et al., 2024), stating that certain rules of inhibitory connectivity found in a mouse electron-microscopic (EM) connectome (Schneider-Mizell et al., 2025) cannot be explained by an non-selective pruning axo-dendritic appositions. We demonstrated the utility of the model to explore the effect of any connectivity rule, using the example of a mouse EM-compatible ‘SM-connectome’. After rewiring and refitting synaptic, and conductance injection parameters, we found that the purely perisomatic targeting by PV+ neurons and the increased inhibitory targeting by VIP+ neurons shifted how much individual layers were driven by intrinsic and extrinsic excitation. This could not be explored with current techniques in vivo and highlights a strength of in silico experimentation. On the other hand, this data source contradicts some of our rat nbS1 parameters, for example, the number of synapses per connection from light microscopy (Markram et al., 2015) do not match with the ones obtained from mouse V1 EM. Thus we used the rat data sources for building our baseline model, but make both models available in a simulatable format (see ‘Data availability’).