Overview of the physiology and simulation workflow.

1. Anatomical model: Summary of the anatomical nbS1 model described in the companion paper. 2. Neuron physiology: Neurons were modeled as multi-compartment models with ion channel densities optimised using previously established methods and data from somatic and dendritic recordings of membrane potentials in vitro. 3. Synaptic physiology: Models of synapses were built using previously established methods and data from paired recordings in vitro. 4. Compensation for missing synapses: Excitatory synapses originating from outside nbS1 were compensated with noisy somatic conductance injection, parameterized by a novel algorithm. 5. In vivo-like activity: We calibrated an in silico activity regime compatible with in vivo spontaneous and stimulus-evoked activity. 6. In silico experimentation: Five laboratory experiments were recreated. Two were used for calibration and three of them were extended beyond their original scope. 7. Open Source: Simulation software and a seven column subvolume of the model are available on Zenodo (see data availability statement). Data generalisations: Three data generalisation strategies were employed to obtain the required data. Left: Mouse to rat, middle: Adult to juvenile (P14) rat, right: Hindlimb (S1HL) and barrel field (S1BF) subregions to the whole nbS1. Throughout the figure, the corresponding purple icons show where these strategies were used.

Modeling and validation of neuron physiology.

A1: Example excitatory neuron morphologies. See companion paper for exemplar morphologies of all morphological types. A2: Example morphologies placed within the atlas-based volume. Axons shown in grey. B: Optimized conductance densities for two exemplary e-types. B1: cADpyr e-type on L5 TPC:A m-type, B2: bSTUT and cSTUT e-types on L5 NBC m-type. (The L5 NBC m-type is combined with more e-types than the two shown, see panel D.) Morphologies were visualized with NeuroMorphoVis (Abdellah et al., 2018). Neurite diameters are enlarged (x3) for visibility. Soma and dendrites in black, axon in red. C: e-types used in the model. As used in Markram et al. (2015) and similar to the Petilla terminology (Ascoli et al., 2008)). D: me-type composition. Heatmap shows the proportion of e-types for each m-type. Each me-type is assigned to one of the three I subpopulations. Assignments depending not only on m-type are highlighted by a green box. Particularly, for BTC and DBC m-types, the cACint e-type belongs to the Sst+ subpopulation. E: Validation of dendritic physiology of all L5 TTPCs. Panel reproduced from Reva et al. (2023). E1: Validation of back-propagating action potential (bAP) amplitude (i.e., the dependence of bAP amplitude on distance from the soma) for basal (green) and apical (blue) dendrites. Reference data (in red) comes from Stuart and Sakmann (1994); Larkum et al. (2001) (apical) and Nevian et al. (2007) (basal). Lines show exponential fits for the in silico (green and blue) and in vitro (red) data. Color bar indicates dendritic diameter. E2: Validation of EPSP attenuation. Reference data (in red) comes from Berger et al. (2001) (apical) and Nevian et al. (2007) (basal). Lines and colorbar same as in E1.

Previously published modeling techniques and data to parameterize them that were combined in this work to model the physiology of nbS1 and conduct in silico experimentation.

Modeling and validation of synaptic physiology.

Caption on the following page. A1: Exemplary pair of L5 TTPCs (visualized with NeuroMorphoVis (Abdellah et al., 2018)). Presynaptic cell in gray, postsynaptic cell in red, synapses between them in purple. Neurite diameters are enlarged (x3) for visibility and axons were cut to fit into the figure. Pre- and postsynaptic voltage traces on the top right. A2: Exemplary postsynaptic traces with different STP profiles. A3: Assignment of STP profiles to viable pathways. (Pathways were considered viable if there were at least 10 connections in all eight subregions of the model.) B1: Validation of first PSP amplitudes (see also Table S5). Dashed gray line represents perfect correlation between experimental and model values. Error bars show one standard deviation (also for C1, D and F). B2: Predicted PSP amplitudes of all viable pathways in the circuit. Postsynaptic cells were held at -70 mV using an in silico voltage-clamp. Means were calculated over 100 pairs of neurons with 35 repetitions each. C1 and C2: same as B1 and B2, but showing the CV of the first PSP amplitude (corresponding Table is S6). D: Validation of mPSC frequencies (see also Table S7). E: Location of synapses from VPM fibers (purple) and POm fibers (red) on 38 neurons (dark gray) in a 5 µm radius column (visualized with BioExplorer). F: Validation of thalamocortical EPSP amplitudes as in B1. The four pathways used for the validation are marked with a black rectangle on H1 to its right. G: EPSP latencies (time from presynaptic spike to the rise to 5% of peak EPSP amplitude in the postsynaptic trace). H1 Left: mean VPM evoked EPSP amplitudes on postsynaptic cell types (over 50 pairs). Right: Comparison of normalized in silico amplitudes (normalized by L4 E as in Sermet et al., 2019) to in vitro reference data from Sermet et al. (2019). Heatmap shows model minus reference values, thus positive values indicate a higher normalized EPSP amplitude in our model than in the reference experimental dataset. H2: same as H1 but for POm (normalized by L5 E as in Sermet et al., 2019).

Spontaneous activity: calibration, in vivo-like dynamics, linking structure to function.

A: Seven column subvolume. B: OUµ and OUσ of somatic conductance injection (B1) parameterized by population (B2) determine FRs (B3). Target FRs equal to in vivo firing rates multiplied by PFR. Inset: Table summarising meta-parameters. C: Target PFR (x-axis) plotted against the observed PFR (y-axis) for each E (red) and I (blue) population after calibration. The number of sides of the shape indicates the layer (triangle used for L2/3). D: Euclidean distance between target and observed PFR values (over populations) decreases over iterations for two meta-parameter combinations (upper). Values after final iteration shown for all meta-parameter combinations (lower; dashed line shows termination condition). E: (left) Firing rates, (centre) mean conductance injection and (right) OUµ by population for the non-bursting simulations (1 line per simulation; E and I values separated). F: Mean missing number of synapses vs. the mean conductance injection, averaged over all combinations, for each population. Markers as in C. G1: Max normalised histograms for two meta-parameter combinations (5 ms bin size, 1σ Gaussian smoothing). G2: Effect of meta-parameters on correlation between E and I histograms (5 ms bin size, 1σ Gaussian smoothing; central column). G3: Fourier analysis of spontaneous activity for the 59 non-bursting simulations. H: Activity in flat space over seven consecutive 100 ms windows. Activity smoothed (Gaussian kernel, σ = 1pixel). I: Left: Correlation r-value between the histograms of layer-wise E and I populations for the 59 non-bursting simulations. Right: RC/U by population for the 59 non-bursting simulations. J: Top: The mean node participation of E and I populations in each layer for dimension one (left) and six (right) vs. RC/U for the parameter combination [Ca2+]o = 1.05 mM, ROU = 0.4 and PFR = 0.3. Markers as in C. Line shows linear fit. Bottom: Correlation r-values between mean node participation and RC/U when neurons with highest node participation are not included in the calculation of mean node participation (via a sliding threshold).

Layer-wise population responses to single whisker deflection and active whisker touch stimuli closely match in vivo millisecond dynamics and response amplitudes.

A1: Stimuli are simulated by activating a proportion (FP) of fibers projecting to the central hexagon of the seven column subvolumne. A2: Each fiber is assigned a spike time from a PSTH of VPM activity recorded in vivo for either a single whisker deflection or active whisker touch paradigm. B: Spiking activity and histograms (baseline and max normalised) for each layer-wise E and I population for the parameter combination [Ca2+]o = 1.05 mM, ROU = 0.4, PFR = 0.3 and FP = 10% for a 2.5s section of the 10 whisker deflection test protocol. C: Trial-averaged PSTHs (baseline and max normalised) in response to the single whisker deflection stimulus for all of the parameter combinations which passed the initial criteria, and the in vivo references. D: Spatiotemporal evolution of the trial-averaged stimulus response in flat space for the same parameter combination in B. E: Left: Latencies of different layer-wise E and I populations for all criteria passing combinations in response to single whisker deflection stimuli. Right: Latencies for inhibitory interneuron subpopulations in response to active whisker touch stimuli. Correponding in vivo references are shown. F: Heatmap showing the effect of the meta-parameters on the difference of RE for the entire central column from the corresponding in vivo reference. Parameter combinations which failed the criteria tests are masked. G: Left: Ratio RE between maximum evoked response and spontaneous baseline activity by population. One line per simulation. Simulation line colour shows the difference of RE for the entire central column from the corresponding in vivo reference (color values shown in G). Thick black line shows median of in silico values. White shows in vivo reference. Right: Normalised ratios of each simulation by dividing by RE for the entire central column.

Reproducing and extending detailed experiments exploring the canonical model and encoding of contrast, rate-coded and synchronous information.

First experiment: Reproducing the effect of optogenetic inactivation of L4 pyramidal cells on L2/3 stimulus responses and predicting the role of other pathways through circuit lesions. A: Schematics of whisker kinetics and VPM fiber rates during 500 ms long whisker hold stimulus. Fraction of VPM fibers coding for each kinetic feature are taken from Petersen et al. (2008). B: Mimicking the effect of activation of the Halo inhibitory opsin in silico. Injected somatic hyperpolarizing current mimicking opsin activation (top), and the resulting somatic voltage trace from a combination of injected conductance, current, and synaptic PSPs from the network (bottom). C: Comparison of average traces from selected L2/3 PCs (see Methods) in control (black) and optogenetically inhibited (green) conditions. D: Same as C, but instead of mimicking the optogenetic inhibition of L4 PCs, only the connections to L2/3 PCs are “cut” (compare inset with the one in C). The right part depicts connections systematically cut from PCs in all layers, while the left shows L4 only for a better visual comparison with the conditions of Varani et al. (2022) in C. Second experiment: Contrast tuning modulation by optogenetic activation of PV+ interneurons, and encoding of synchronous and rate-coded by interneuron subtypes. E: Spatial distribution of firing rates of VPM fibers projecting to a single simulated column at different points in time, corresponding to linear drifting gratings with a temporal frequency of ftemp = 2 Hz and a spatial frequency of fspat = 0.001 cycles/µm. F: Firing rate signal of a single VPM fiber corresponding to a random sequence of drifting grating patterns with different contrast levels C, as indicated by the legend. All optogenetic stimuli targeting PV+ (or Sst+) interneurons were completely overlapping with the grating stimulus intervals, as indicated by the blue bars. G1: Spiking activity and PSTH of an exemplary PV+ interneuron over 10 repetitions (trials) in response to a grating stimulus at maximum contrast (C = 1.0) without (control) and with medium strength (150%) optogenetic stimululation. H1: Contrast tuning modulation of an exemplary PV+ interneuron by different levels of PV+ optogenetic stimulation, ranging from 0% (control) to 300%. The individual markers denote actual data points (mean SD over 10 repetitions, normalized to baseline response at maximum contrast) while curves illustrate sigmoidal fits to these data points. G2/H2: Same as C1/D1, but for exemplary PCs. H2 inset: Illustration of sigmoidal parameters m, n, Rmax, and c50.

Third experiment: Encoding of synchronous and rate-coded by interneuron subtypes. I: A binary signal is encoded either through changes in the rate or synchronicity of optogenetic pulses applied to a set of PCs. We define a direct correspondence between single optogenetic pulses and single spikes. Spiking activity is measured in PV+, Sst+ and 5HT3aR+ cells, and the mutual information between the binned activity of individual cells and the original binary signal is measured. J: Visualization of the rate (upper) and synchronous (lower) coded stimulus experiments stimulating 1000 PCs, showing the binary signal (top), input neuron spike trains for 40 neurons (middle), and responses of the three L2/3 neuron types (bottom). K: Upper: Results for the rate coded stimulus experiment. Left: Mutual information between spiking activity and the binary signal (one point for each cell that spiked at least once). Only activity in the 50 ms following the change of the binary signal is considered. Cells with mutual information not significantly different from that of a shuffled control are coloured grey. Centre: Same as left but considering all time bins. Right: The effect of the number of stimulus neurons on the mean single cell mutual information for each subpopulation. Lower: The same as upper but for the synchronous stimulus type.

The effect of inhibitory targeting trends observed in electron microscopy on spontaneous and stimulus-evoked activity.

A1: Panel illustrating interneuron synaptic targeting features and caption taken from Schneider-Mizell et al. (2023) (under CC-BY 4.0 license). i) For determining inhibitory cell subclasses, connectivity properties were used such as an axon (green) making synapses (green dots) the perisomatic region of a target pyramidal cell (purple). ii) Dendritic compartment definitions for excitatory neurons. iii) Cartoon of a multisynaptic connection (left) and the synapses within the multisynaptic connection considered ”clumped” along the presynaptic axon (right). A2: Exemplar neurons before and after rewiring showing synapse locations from different inhibitory m-type groupings: Distal targeting (DistTC: martinotti, bipolar, double bouquet cells), Perisomatic targeting (PeriTC: large basket, nest basket cells), Sparse targeting (SparTC: descending axon, neurogliaform, horitzontal axon cells), Inhibitory targeting (InhTC: bitufted, small basket cells). A3: Heatmaps showing the mean fraction of efferent synapses from the inhibitory m-type groupings onto inhibitory neurons, somas, proximal dendrites, apical dendrites, distal dendrites, or that are part of clumped or multisynaptic connections. Heatmaps show mean fractions (left to right) for the original, Schneider-Mizell et al. (2023) characterization of the MICrONS dataset ((MICrONS-Consortium et al., 2021)), SM-connectome, and the change between the original and SM-connectome. B1: Proportional change in mean conductance injection between original and SM-connectome by population. B2: Distribution of ratios between connected and disconnected firing rates (RC/U) by population over non-bursting meta-parameter combinations for original and SM-connectome. B3: Comparison of correlation (r-values) between E and I histograms (5 ms bin size, 1σ Gaussian smoothing; central column) between original and SM-connectome (each point represents a specific meta-parameter combination). C: Scatter plot showing the change in the mean number of afferent inhibitory synapses onto each E (red) and I (blue) population after calibration versus the proportional change in mean conductance injection. The number of sides of the shape indicates the layer (triangle used for L2/3). D1: Proportional change in the the ratio between the peak of the evoked response and the pre-stimulus baseline activity level (RE) between the original and SM-connectome. Plot shows distribution over meta-parameter combinations which passed the the initial assessment of in vivo similarity to the in vivo data based on latencies and decays (see Methods). D2: Change in the latencies of populations peak responses between original and SM-connectome (distribution over same meta-parameter combinations shown in D1).

Full nbS1 simulations: independent functional units, millisecond-precise stimulus-responses, selective propagation of stimulus-evoked activity through mid-range connectivity, linking structure to function.

Caption on following page. A1: Mean FRs across the model. Columns 0 and 39 (subsequently used) are highlighted. A2: Spiking correlations of E/I populations in 240 hexagonal subvolumes (diameter: 460 µm). B: Spiking activity and max normalised histograms for layer-wise E and I populations for the parameter combination in A, and the 0th (upper) and 39th (lower). C: Target vs. observed PFR values for the four non-bursting simulations. The number of sides of each marker corresponds to the layer (triangle represents L2/3). Red: E, blue: I.. D: Correlations between E and I populations by layer for the four non-bursting simulations for the two hexagons (top and bottom). Color of line from light to dark represents PFR. Middle: Distributions of spiking correlations of layer-specific populations in columns (520 µm) of the nbS1 model (black dots; mean values: grey). Compared to the central column when the seven column subvolume is simulated in isolation (green). E: Trial-averaged activity in flatspace following single whisker deflections. Lower: Subset of time windows shown with only the top 2% of bins for each time window after baseline activity level subtracted. F: Trial-averaged PSTHs (baseline and max normalised) for column 0. For FP = 15% and FP = 25% respectively, the evoked responses passed and failed the latency and decay-based criteria tests for similarity to in vivo. Response for FP = 20% passed the tests for all populations except L5 E which showed a more sustained response. G: Trial-averaged response of column 39. H: Correlations of mid-range inputs. Upper: Distribution of spiking correlations of inputs into 240 subvolumes. Calculated at reduced spatial resolution, based on connection counts and correlations between the subvolumes (see Methods). Subvolumes are sorted by the E/I correlation of neurons within them. Lower: Mean of the correlations of inputs into 240 subvolumes vs internal spiking correlation for all subvolumes during spontaneous activity. Black line: linear fit. Calculated based on connection counts and correlations between 50 µm hexagonal subvolumes (see Methods). Data for evoked activity is shown in Figure S12. I: Upper: Spiking correlations of subvolumes against their distances in spontaneous (red) and evoked (orange: Fp = 0.15%, blue: Fp = 0.2%, green: Fp = 0.25) activity. For increased spatial resolution, smaller (58 µm) subvolumes were used, hence correlation values are not comparable to values in B. Lower: Distribution of soma distances of neuron pairs connected by local (orange) or mid-range (blue) connectivity. Green: distances of all pairs, independent of connectivity. J: Upper: Classification of subvolumes based on low or high internal correlation and membership in mid-range rich club Data for evoked activity shown in Figure S12. Lower: E/I correlation of subvolumes for members of the rich club (red) and non-members (blue). K: FRs of subvolumes following a stimulus with FP = 0.15% (top) and FP = 0.2% (bottom). Mean (lines) and SEM (shaded area) over 10 repetitions shown for subvolumes directly innervated by the VPM stimulus (black), strongly innervated by directly innervated subvolumes (with over 106 synapses, red), or by medium strength (over 105 synapses, yellow) or weak (over 104 synapses) indirect innervation. Dashed lines indicate locations of peaks

Key resources table

Predictions

Excitatory synaptic pathways.

Average class parameters are marked in bold and are used predictively (in lack of reference in vitro data) for the remaining pathways belonging to the same class. Physical dimensions are as follows: peak conductance g^: nS, depression and facilitation time constants D, F, and the EPSC τdecay: ms, the Hill coefficient of the nonlinear [Ca2+]o dependent scaling of release probability UHill: mM, the release probability USE, the average number of vesicles in the release-ready pool NRRP, and the NMDA/AMPA ratio g^ratio are dimensionless.

Inhibitory synaptic pathways.

Average class parameters are marked in bold and are used predictively (in lack of reference in vitro data) for the remaining pathways belonging to the same class. Physical dimensions are as follows: peak conductance g^: nS, depression and facilitation time constants D, F, and the IPSC τdecay: ms, the Hill coefficient of the nonlinear [Ca2+]o dependent scaling of release probability UHill: mM, the release probability USE, the average number of vesicles in the release-ready pool NRRP, and the GABAB/GABAA ratio g^ratio are dimensionless.

Thalamocortical synaptic pathways.

Values taken from the internal connectivity (Table S2) are marked in bold. Physical dimensions are the same as in Table S2.

Validation of PSP amplitudes (see Figure 3B1)

Validation of first PSP amplitudes’ CVs (see Figure 3B2)

Validation of mPSC frequency (see Figure 3E2)

Spontaneous activity calibration: estimating χ and ϕ.

A: Estimation of χ and ϕ for L5 E (illustration). ΨCa2+,ROU : CFR 1→ OUµ is estimated by dividing it into the two functions χ and ϕ, which are determined separately. Left: χ : (ROU, UFR) OUµ. The color of each point shows the mean unconnected FR in a single column for L5E for different combinations of OUµ (x-axis) and OUσ (y-axis). For a given value of ROU interpolation is used over the datapoints to estimate the predicted FRs along a line of gradient ROU starting at the origin. As FRs increase monotonically with increasing OUµ and OUσ, the interpolation can be used to estimate the OUµ and OUσ combination that will produce a given target unconnected FR. Right: ϕCa2+,ROU : UFR 1→ CFR. In each iteration, 10 simulations of the seven hexagon subvolume are run for a range of target PFR values. For each population (i.e. L5 E) an exponential is fit to unconnected vs. connected firing rates, and used to estimate unconnected firing rates which will achieve target connected firing rates on the next iteration. B: Visualisation of the metaparameter space. The algorithm is run for a single combination of Ca2+ and ROU (blue), and is then generalised for other combinations (red). C: Data for estimation of χ for each population. Population FRs from unconnected simulations of a single column. Simulations vary by OUµ and OUσ. Colour shows FR. Black circles show region where FRs are greater than 0 and less than the population’s in vivo reference FR. The dependence of unconnected FR on OUµ and OUσ is heterogeneous across populations. D: Final estimation of ϕ1.1,0.4 for each neuron group. Exponential fit of unconnected vs connected FRs for each population after 5th iteration. The final fits are heterogeneous across populations.

In silico experimental methods.

A: Synaptic pathway lesions. After selecting a pre- (green) and post-synaptic (blue) population of neurons, connections between them in the specified direction are removed (red arrows). B: Optogenetic manipulations. B1: Decay with depth of the strength of a stimulus applied to mimic the effect of optogenetic manipulations with 470 nm (blue) and 595 nm (yellow) light. Indicated for exemplary intensities (values of I0). B2: Targeted neuronal populations in two in silico optogenetic experiments. Density of affected cells (combination of expression level and light depth) indicated by the corresponding color of light. Unaffected populations in grey. Left: Excitation of inhibitory populations. Right: Inhibition of PCs in L4. C: Stimulation with thalamic inputs. C1: Each thalamic fiber in the model is assigned a position in a flat coordinate system. Based on the coordinates, type of fiber (VPM vs. POm), or randomly, time series of sensory stimulation are assigned to the fibers. C2: Sensory stimuli are transformed by a transfer function capturing pre-thalamic processing. The transfer function can simply be identity. The transformed signal is then used with a spiking model to generate stochastic spike trains that serve as inputs into a simulation.

Firing rate distributions.

A: Firing rate distributions of different subpopulations for one of the parameter combinations. B: Non-spiking cells are excluded for improved visualisation of the distributions.

Whisker deflection responses - supplementary.

A: PSTHs in response to the single whisker deflection stimulus for all parameter combinations. PSTHs of simulations which passed the initial criteria are coloured green, whilst those which did not are coloured grey. B: Layer-wise PV+ and Sst+ subpopulation PSTHs for simulations which passed the criteria. C: Analysis of evoked responses of L5 E under two simulation cases. The first case uses the mean of the VPM to L4 E and L6 E synaptic conductance for the VPM to L5 E conductance, whilst the second uses the maximum of the two values (VPM to L6 E). C1: For the two simulation cases, L5 E PSTHs for the criteria passing (green) and failing (grey) simulations. For the mean case, some of the simulations that do not burst or do not have a secondary rebound, fail the criteria because of low signal to noise ratio. C2: Ratios between evoked and spontaneous activity of L5 E population for criteria passing simulations. There are fewer criteria passing simulations for the mean case, and those that do are much lower than the in vivo reference (40.2). D1: Proportion of spiking cells by population in the central column over 100ms following stimulus onset. Each line represents a different criteria passing simulation. D2: Heatmap showing the effect of the meta-parameters on the proportion of spiking cells over the entire central column. E: Heatmap showing the effect of the meta-parameters on the proportion of spiking cells which spike only once over 100ms following stimulus onset.

Active whisker touch responses - supplementary:

Same as Fig. 5 but using an active whisker touch VPM stimulus and in vivo reference for the active whisker touch paradigm (see Methods).

Additional validations: reproducing the experiments of Varani et al. (2022).

A: Raster plots of the microcircuit’s activity and population firing rates below. A1: Control conditions, A2: 95% of L4 PCs inhibited (by direct somatic current injection depicted in B above, see Methods). B: Voltage traces of all L2/3 (top) and all L4 (bottom) PCs. Panels show spiking traces (top), and subthreshold traces (bottom). B1 and B2 depict the same conditions as C above.

Additional validations: reproducing the experiments of Varani et al. (2022).

A1: Scatter plots showing the fitted sigmoidal parameter values in control condition vs. PV+ optogenetic stimulation at lowest (50%) and highest (300%) level respectively for all 259 robustly tuned PV+ interneurons. Diamond markers indicate mean values over neurons, including all optogenetic stimulation levels from 50% to 300%. Colors as in panel D1 legend. B2: Same as B1, but for all 228 robustly tuned PCs.

Calibration of the drifting grating stimulus based on peak firing rates.

The drifting grating stimulus was calibrated by running a parameter scan of 45 simulations to determine optimal values for FP, Rbk, and Rpeak. The maximum (top) and average (middle) of the first and second peak firing rates r1 and r2 (connected points) over the whole population of PCs, together with the average of the normalized peak difference (Michelson contrast) r^diff (bottom) are summarized. Different colors indicate different contrast levels C ∈ [0.06, 0.12, 0.5, 1.0]. The rate threshold rth = 30 Hz at maximum contrast, the calibration target values that were taken into account for selecting the optimal parameter set, as well as the optimal parameter set (FP = 1.0, Rbk = 0.2 Hz, and Rpeak = 10.0 Hz), are highlighted.

Calibration of the drifting grating stimulus based on sigmodal parameters.

The drifting grating stimulus was calibrated by running a parameter scan of 45 simulations to determine optimal values for FP, Rbk, and Rpeak. The values of the four parameters m, n, Rmax, and c50 for fitting a sigmoidal function to the average firing rates of the whole population of PCs within the full stimulation interval are summarized. Dashed lines indicate fitted parameter values when considering only the first or second half of the response interval. Calibration target values that were taken into account for selecting the optimal paramter set, as well as the optimal parameter set (FP = 1.0, Rbk = 0.2 Hz, and Rpeak = 10.0 Hz), are highlighted.

Modelling the effects of optogenetic stimulation of interneurons.

A1-A3: Divisive scaling (Div), subtractive shifting (Sub), and saturation additive (SA) model fits (dashed lines) to responses of an exemplary PV+ interneuron (solid lines) at different levels of PV+ optogenetic stimulation, ranging from 0% (baseline) to 300%. B1: Goodness of Div/Sub/SA model fits (µ SEM of r2 score) for different depolarization levels over all 259 robustly tuned PV+ interneurons, showing that direct photostimulation effects on interneurons can be best described by the saturation additive model. C1: Polar plot of S (saturation) and A (additive) parameter values of the SA model fits to all 259 robustly tuned PV+ interneurons at different levels of PV+ optogenetic stimulation. Diamond markers indicate mean values (in polar coordinates) over neurons, colored arc lines the circular SD of the polar angles. B2: Same as B1, but for conductance-based model fits describing indirect photostimulation effects on all 228 robustly tuned PCs, assuming a saturating additive model description of interneurons. C2: Same as C1, but for the conductance-based model fits to all robustly tuned PCs.

Re-calibration of SM-Connectome.

A: The SM-connectome took 6 iterations to re-calibrate for the full 60 meta-parameter combinations. The mean conductance injections found for the original connectome for the 60 meta-parameter combinations were used for the first iteration of the re-optimisation and lead to decreases in firing rates for all inhibitory populations except L1 I (A1). In A1 and A2 blue and yellow points refer to Ca2+ 1.05 and 1.1 respectively. After 6 iterations the firing rates converged to the target firing rates (A2-A5) and there was minimal difference in population firing rates over the 60 meta-parameter combinations between the original and SM-connectome (A4). B Correlation r-values between the E and I populations, as for the original connectome in the main text. C: Unconnected vs connected MFRs for the different populations after the 6th iterations. Unlike for the original connetome, firing rates decrease following connection of the network.

Full nbS1 - supplementary:

A: Evoked responses of a single hexagonal subvolume over 10 single whisker deflections. B: Trial-averaged activity of the full circuit for FP = 15% following the single whisker deflection stimulus (top). Middle and bottom: trial-averaged activity of the full circuit for FP = 15% following the single whisker deflection stimulus for FP = 15% and FP = 20%, showing only the top 2% of bins for each time window after the baseline activity level has been subtracted. C: Distribution of spiking correlations of inputs into 240 subvolumes. Calculated at reduced spatial resolution, based on connection counts and correlations between the subvolumes (see Supp. Methods). Subvolumes are sorted by the E/I correlation of neurons within them. For evoked activity with FP = 20%. D: Mean of the correlations of inputs into 240 subvolumes (as in C) against internal spiking correlation for all subvolumes during spontaneous activity. Black line: linear fit. Calculated based on connection counts and correlations between 50 µm hexagonal subvolumes (see Supp. Methods). E: Classification of subvolumes based on low or high internal correlation and membership in the mid-range rich club. Data from evoked activity with FP = 20%. F: Indegree from the subvolume innervated by the VPM stimulus against correlation with that subvolume, indicated for all other subvolumes. Color indicates the Δt with maximum correlation.