References to publications of input data and methods employed for individual modeling steps. An asterisk next to a reference indicates that substantial adaptations or refinements of the data or methods have been performed that will be explained in this manuscript. In the other cases, a basic summary will be provided and an exhaustive description in the Methods.

Overview of the model building and analysis workflow.

Step 1: Building was based on a volumetric atlas of the modeled regions: (1) S1J; (2) S1FL; (3) S1Tr; (4) S1HL; (5) S1Sh; (6) S1DZ; (7) S1DZO; (8) S1ULp. Additional atlases of biological cell densities and local orientation towards the surface were built. Step 2: Neuron morphologies were reconstructed and classified into 60 morphological types (m-types). They were placed in the volume according to the densities and orientations from step 1. Step 3: The anatomy of intrinsic synaptic connectivity was derived as the union of one algorithm for local connectivity and one for mid-range connectivity. Step 4: Extrinsic inputs from two thalamic sources were placed on modeled dendrites according to published methods. Step 5: Taken together, these steps allowed us to predict the topology of connectivity at scale with (sub-)cellular resolution. Step 6: The anatomical model served as the basis of a physiological model, ready to be simulated. This is presented in an accompanying manuscript. Step 7: The model, simulation and analysis tools have been made publicly available. Left: During modeling, three types of generalization had to be made to fulfill input data requirements: from mouse to rat, from adult to juvenile, and from one cortical region to another. Generalizations used are indicated in each step.

Anatomy of the model

A: Exemplar 3D reconstructions of the 18 excitatory m-types. Rendering and visualization was done in NeuroMorphoVis (Abdellah et al., 2018). Dendritic diameters are scaled (x3) for better resolution. B: Modeled cortical layers, with exemplars of excitatory morphological types placed in the model. C: The placement of each morphological type recreates the biological laminar anatomy of dendrites and axons, which then serves as the basis of local connectivity. Inset shows modeled brain regions (solid colors) in the context of the non-modeled regions (transparent).

Volumetric anatomy of the model

A: M-type composition per layer: For each layer, stacked histograms of the relative fractions of each m-type, comparing the model (right bar for each layer) to the input data (left bar). For simplicity, the layer designation is stripped from each m-type. Left: For excitatory types; right: Inhibitory types. B: Stacked histograms of the fraction of space filled by neurites at various depths in the model. The y-axis indicates the distance in µm from the bottom of layer 6. B1: For neurites of neurons in different layers indicated in different colors. Grey: estimated lower bound for the volume of axons supporting the mid-range connectivity. B2: For neurites of different m-types. Colors as in A, but inhibitory types are grouped together. B3: For different types of neurites. C: Comparing fractions for axons and dendrites to the literature (Santuy et al., 2018). X-marks indicate overall means, other marks means in individual layers; teal: reference, black: model.

Intrinsic connectivity as union of local and mid-range connectivity.

A: Local connectivity was derived as in (Markram et al., 2015; Reimann et al., 2015), i.e. based on neuron placement and morphologies. Right: Axo-dendritic appositions were considered as potential synapses. They were then filtered to yield among other biological constraints biologically realistic bouton densities. B: Range of resulting connection probabilities within 100µm. Indicated are the minimum (inner circle) and maximum (outer circle) found over 20 repetitions of sampling 1000 pairs in a simulated multi-patch procedure (Fig S9A). See Methods for mapping from morphological types to inhibitory subclasses (PV/SST/INH). C: An exemplary layer 5 PC axon reconstruction from Janelia MouseLight (mouselight.janelia.org, Gerfen et al. (2018)) that includes long-range collaterals, shown in the context of mouse somatosensory regions. The blue circle highlights local branches all around the soma. Red highlights depict more targeted collaterals into neighboring regions. D1: Predicted pathway strengths as indicated in Figure S2C3. D2: Pathway strength emerging from the application of the apposition-based connectivity algorithm described in Reimann et al. (2015). E: Values of the diagonal of (D1) compared to the diagonal of (D2). F: Schematic of the strategy for connectivity derivation: Within a region only apposition-based connectivity is used; across regions the union of apposition-based and mid-range connectivity. G: Connection strength constraints for the mid-range connectivity derived by subtracting D2 from D1 and setting elements < 0 to 0 (colors as in D). H: Resulting total density from both types of modeled synaptic connections in individual regions compared to the data in (D1).

Recreating an electron-microscopic specificity analysis

In the top part we analyzed synapses, axons and dendrites in an approximately cubic volume in L4, emulating the techniques of Motta et al., 2019. In the bottom part we analyzed a 100 × 100 µm volume of the MICrONS dataset defined by Schneider-Mizell et al., 2023. A: Fraction of synapses placed on somata (SOM), proximal dendrites (PD), smooth dendrites (SD), apical dendrites (AD) and axon initial segments (AIS) for axon fragments in the sampled volume. Black bars indicate mean values over axons, arrows indicate binomial probabilities fit to observations of axons forming at least a single synapse; pink: for excitatory axons, black: inhibitory. Fractions missing from 100% are onto other compartment types. B: Overall count on synapses on different postsynaptic compartment classes inside the studies volumes, comparing the data of Motta et al., 2019 to the model. Colors as indicated in the legend. C: Distributions of the fractions of synapses onto different compartment types over excitatory (left) and inhibitory (right) axons. Expected from the binomial model in A (grey) against the observations in our anatomical model (black outlines). D: Fraction of axons with significantly increased synapse counts onto different compartment types compared to the binomial control. Indicated for two values of the false detection rate criterion (q=0.05, 0.3, Storey and Tibshirani, 2003). Comparing the reference data to our anatomical model. E: Top numbers of neurons in four connectivityderived classes in the 100 × 100 µm volume of the MICrONS dataset, defined by Schneider-Mizell et al., 2023. Bottom: Numbers in a comparable volume of the model. For assignment of m-types into the four classes see the main text. F: Derivation of alternative connectivity with targeting specificity, using the example of the “PeriTC” class. The axo-dendritic appositions of an axon (top) are classified as matching the targeting (green box; here: appositions with somata or proximal dendrites) or not. Non-matching appositions are removed with probability pnt (middle). This is followed by a non-specific removal of connections (formed by single or multiple synapses) until the biological density of synapses on the axon is met (bottom).

G: Bouton densities resulting from the process in F. Left to right: For PeriTC neurons, targeting somata and proximal dendrites; for InhTC neurons, targeting inhibitory neurons; for SparTC neurons, targeting the first synapse of a connections; for SparTC neurons, without targeting preference. Respective optimized pnt values reported in the main text. H: Resulting targeting of postsynaptic compartments, fraction of synapses in multisynaptic connections and clumped synapses (see Schneider-Mizell et al., 2023). Red: Default model; blue: optimized alternative model with specific targeting; grey: optimized alternative model without specific targeting; orange: data of Schneider-Mizell et al., 2023 Boxes outline the characteristic feature of each class. I: Inhibitory targeting specificities combining the best performing models.

Anatomy of thalamic innervation.

A: Depth profiles of synapse densities (dashed lines) and mean number of thalamic inputs per neuron (solid lines) for core- (red) and matrix-type (blue) thalamocortical projections. Shaded area indicates the standard error of mean. For a region-specific validation of synapse densities (dashed lines) against experimental data, see Figure S13. B: Mean and standard deviation of the number of thalamic inputs for neurons in individual layers or all neurons. Colors as in A. C: Common thalamic innervation (CTI) of an exemplary neuron (black dot) and neurons surrounding it, calculated as the intersection over union of the sets of thalamic fibers innervating each of them. Scale bar: 200 µm. D: CTI of pairs of neurons at various horizontal distances. Dots indicate values for 125 randomly picked pairs; lines indicate a sliding average with a window size of 40 µm. We perform a Gaussian fit to the data, extracting the amplitude at 0 µm (A) and the standard deviation (σ). E: Values of A (unitless) and σ (in flatmap coordinates) for pairs in the individual layers or all pairs. Colors as in A.

A1, A2: Parcellation of the modeled volume into 230 µm radius columns. Exemplary slice of columns highlighted in green. B: Geometrical metrics of column subvolumes in the flat view. Peripheral columns masked out (grey); green outline: highlighted columns in A. Left: Conicality, defined as the slope of a linear fit of depth against column radius. Negative values indicate narrowing towards L6. Right: Column height, i.e. cortical thickness at the location of the column. C1: Modularity of the networks of connections within each column. C2: Column conicality and height colored by modularity; D1: Conicality of columns against their laminar neuronal composition, normalized against the overall composition of the model. Colored lines indicate linear fits. D2: Conicality against the density of connections in subnetworks given by the intersections of columns with individual layers. E: Counts of afferents formed onto neurons in individual columns from neurons in the entire model. E1: Normalized in-degrees originating from neurons in individual layers, plotted against conicality. E2: In-degree of neurons in individual layers normalized by the overall in-degree in the model into each layer plotted against column height. F: r-values of linear fits against generalized, n-dimensional in-degree as in E2. F1: Of generalized in-degree against height; F2: Of generalized in-degree against conicality.

Global connectivity structure local vs. mid-range.

A1: Simplex counts of the local connectivity network and several types of random controls (see text). Examples of d-simplices for d = 1, 2, and 3. Intuitively, an n-simplex is formed by taking an n − 1-simplex and adding a new node to which all neurons connect (sink). A2: In orange shades, simplex counts of the mid-range connectivity network and several types of random controls. In gray, simplex counts of the combined network and in blue the sum of the local and the mid-range simplex counts. Inset: local, mid-range and the sum of simplex counts on a linear scale. B1: Node participation per layer, normalized to add up to one in each dimension. Top/bottom row: Local/mid-range connectivity network. Left to right: Participation as source, in any position, and as a sink of a simplex. B2: From left to right: Spatial location of the cells in the simplicial n-cores in flat coordinates for all n. Depth coordinates for each connected component in the simplicial cores (cells participating in simplices of maximal dimension) each dot marks the depth of a neuron in that component. Adjacency matrix of the simplicial cores with respect to local (green dots) and mid-range (orange dots) connections. B3: Number of unique neurons in the simplicial cores present in each position of a simplex (i.e., position along the order from source to target. numbers: total counts, bars: counts relative to the total). C1: Right: Distribution of path distances between pairs of neurons in the local simplicial core; green: along only local edges; grey: along all edges. Left: Euclidean distances between pairs at a given path distance. Arrows: See C2. C2: Total Number of neurons in all paths of length 3 between neurons of the local core which are at distance 3 in the combined circuit (black arrow in C1) split by location. Red/purple: at positions 2 and 3 respectively; grey: expected from randomly assigned m-types while maintaining the global distribution. Cross-/stripe-patterned: For pairs at path distances greater/less-or-equal than 3 along local edges only (green arrows in C1).

Overview of the building workflow and tools used. Names of individual software tools are bold in the description (see also Key resources table); descriptions of required inputs in italics. Where not listed, the inputs of a step are outputs of the previous step. For completeness, also physiological modeling steps not described in this work are also depicted, but indicated semi-transparent. Note that the Key resources table only lists tools related to anatomical modeling.

Anatomical, morphological and other aspects affecting connectivity and our predictions for their relevance for efferents of different neuron types. See main text for an explanation of the individual aspects. The signs in the table indicate whether we predict a certain aspect to be not relevant (-), relevant (+), or highly relevant (++) for a given neuron type.

Morphological types (m-types) used in the model

Predictions

Derivation of neuron density depth profiles. Left: A vertical profile of neuron densities, calculated from antibody stains of neuronal nuclear protein (NeuN). From left to right: Neuron densities in each bin are split into individual morphological types through antibody stains for various markers (blue boxes) and an established mapping of markers to types (green), or by applying established neuronal compositions (red). Information in this flow diagram is illustrative, not quantitative.

Data sources for connectivity modeling.

A: Mapping from mouse to equivalent rat subregions. Left: flat view of mouse somatosensory subregions, Right: equivalent flat view for rat. Corresponding subregions are indicated with matching colors. Triangles are drawn on top of each subregion (shown for SSp-ul and S1FL) using manually annotated points that are assumed to correspond to each other. Inset: These triangles define affine transformations between mouse and rat subregions. B: Transformed rat subregions (colored outlines) are shown overlaid to their corresponding mouse subregions (colored areas). The topographical mapping of connections between subregions predicted in Reimann et al. (2019) (black arrows) can thus be generalized to rat connectivity. C: Mean connection densities between rat somatosensory subregions were derived from the Allen Mouse Brain Connectivity Atlas by summing over corresponding mouse somatosensory subregions. Labels on the left side of the connection density matrix describe the mapping applied from mouse to rat subregions. D: The topographical mapping of long-range connectivity was parameterized as in (Reimann et al., 2019): For each pair of subregions, e.g. SS1FL and S1HL, affine transformations TS1FL and TS1HL were defined that transformed the top-down views of the regions into a common coordinate system. Connectivity was then distance-dependent in the common coordinate system. Transformations that were optimized for mouse data in (Reimann et al., 2019) were adapted through the process in A and B for use in rat. This described connectivity where parts of one subregion were predominately connecting to specific parts of another subregion. The plot indicates parts that are preferentially connecting to each other in the same color. Sharp transitions at subregion boundaries are the result of modeling connectivity separately for each pair of subregions in this iteration of the model. E: For the vertical dimension, predicted layer profiles of connections between subregions of (Reimann et al., 2019) are generalized to rat. One example is illustrated. E1, left: Exemplary slice through the SSp-ll region of the Allen Mouse Brain Atlas. Right: Layer profile of synapses for a feedforward connection targeting layer 4. E2: Generalized version we used for rat. Right: Slice through the S1HL region of the atlas. Left: Generalized layer profile, still targeting layer 4.

Input data for the derivation of the locations of thalamic inputs.

A1: Layer profile of bouton densities of VPM axons in the rat barrel field. Light red: Digitized from Meyer et al. (2010). Dark red: Binned into 20 bins per peak with a lower cutoff of 1/mm3 applied. A2: Same for POm axons. B: Total length of 11 reconstructed axons of Janelia MouseLight neurons with somas in VPM in various somatosensory subregions. Dashed line: Median C: Standard deviation of a Gaussian fit of the reconstruction points of the same axons around their centroid in somatosensory areas. Dashed line: Median

Normalized depths used for layer boundaries

Topological comparison as in (Kanari et al., 2018, 2019) of in vivo stained and in-vitro (i.e., from slices) axonal (A) and dendritic (B) reconstructions of rat somatosensory cortex of pyramidal cells from layers 2-6. Top row presents the topological profiles of in vivo stained reconstructions, second row presents the topological profiles of in slices stained reconstructions and bottom row the difference between them (red: in vivo stained, blue: in slices stained). Number of cells per layer are reported in individual pannels.

Comparison of representative reconstructructed morphologies that have been stained in slice (top) or in vivo (bottom). One pyramidal cell morphology per layer is depicted.

Classification and morphometrics of excitatory morphologies.

A: Same as Figure 2A. B: Morphometrics of all morphologies belonging to the 18 excitatory m-types.

Placement features used for modeling. The right column denotes the vertical target interval by specifying a layer and a relative offset within the layer, with 0.0 indication the bottom and 1.0 the top of the layer.

Scoring the placement of a neuron morphology for a voxel.

A: Neurite features, here the apical tuft, were manually given a vertical annotation interval (grey) and assigned a target interval, expressed as a layer interval (red). (Note: parts of the apical dendrite visually shortened.) B: Then the placement of a neuron morphology in a given voxel (blue) is scored as follows: Normalized depth values of the target interval are calculated (red); the normalized depth of the voxel is looked up in an atlas (blue); it is used to calculate the normalized depth of the annotation interval when the morphology is placed in the voxel (grey); a score is calculated as in (Markram et al., 2015) based on the overlap of the intervals.

Number of thalamic fibers providing inputs to the model and each of its subregions.

Exemplary neurons in the model, one per m-type, rendered in the context of a slice spanning all cortical layers (grey borders).

Comparing local connection probabilities to the literature.

A: Local connection probabilities were measured by emulating a multi-patch clamp sampling procedure (see Methods). Resulting relative locations of 1000 sampled pairs of excitatory neurons in layer 5 (L5EXC) are indicated. Black triangle location of the potential presynaptic partner. Grey: Locations of unconnected neurons. Orange: Of connected neurons. Top: In the x-y plane that is roughly parallel to layer boundaries. Bottom: In the x-z plane; the z-axis is roughly orthogonal to layer boundaries. B: Violinplots of distributions of connection probabilities reported in 124 literature sources, gathered by Zhang et al. (2019). Colors of the violins indicate means of the distributions. Numerals indicate the number of literature sources used. Arrowheads indicate the mean values measured in the model. They are black, if the value falls inside the distribution of literature values; purple if they are outside, but within 25% of the variability; red if they outside that. Where only a single literature source was available, 10% variability was assumed. Literature sources for mouse and rat were used. C: As B, but only literature sources where the reported animal age range contained P14, the target age of the model. D: As B, but only literature sources for rat were used.

Validation of modeled connectivity.

A, B: Touch connectivity is validated by comparing emerging bouton densities (A) and mean numbers of synapses per connection (B) to the target values from the data. Each data point depicts a single morphological type (A) or type-specific pathway (B). Black arrows: Bouton density for Chandelier Cells (ChC) in layers 2, 3 and 4. C: Mean densities of synapses from the union of touch connectivity and mid-range connectivity in pathways within and between regions (C2) compared to that target values from the data (C1). D: Structure of the topographical mapping of connections; each part of the model predominately connects to neurons at equally colored locations. (D1) Target mapping from the data. (D2) Analyzed in the union of touch connectivity and mid-range connectivity. E: Layer profiles of synapses placed in mid-range connections between regions (blue bars) compared to the predicted target profiles from (Reimann et al. (2019), pink lines). Depicted is the density of synapses in a depth bin relative to the overall mean density over the entire depth.

Caption

Modeling thalamo-cortical synaptic inputs.

A: Input vertical profiles of thalamocortical synapses from core-type projections (A1) and matrix-type projections (A3). A2: Exemplary model neuron with projection synapses placed on dendrites according to the prescribed densities. B: Modeling of afferent thalamic fiber locations: Depicted is an exemplary slice of the model. Shaded vertical areas indicate the upper edges of layers, with colors consistent with the rest of the manuscript (L1: yellow, L2: orange, L3: red, L4: pink, L5: blue, L6: green). For each fiber a location at the border between layers 4 and 5 is randomly chosen (10 examples shown; white spheres). A point 1500µm towards the bottom of layer 6 is chosen as the starting point of each fiber (black dots). A domain of influence is then assigned around a line starting at that location with the indicated direction (black arrows and red areas). Influence weakens with distance from the line with a Gaussian profile; for details, see D, E. C: Locations of fibers in the flat view. 10% of the full density shown. Inset: Locations of approximately 1mm2 shown in full density. The green circle indicates a single standard deviation of the Gaussian of influence strength used for the core style projections (i.e., 120µm). D: For an exemplary synapse its distances to surrounding fibers is measured and the value of the Gaussian for those distances calculated. It is cut off at two standard deviations. E: The probability that a given fiber is chosen as innervator of a synaptic location (blue bar) is proportional to its value in D (shown only for fibers with nonzero probability).

Validation of thalamo-cortical density profiles.

Synapse density profiles of thalamocortical synapses from core-type (VPM) projections (A) and matrix-type (POm) projections (B), compared with the target profiles from Meyer et al. (2010) (pink lines). Depicted are the mean densities (black lines: ± SD) over voxels within depth bins relative to the respective total cortical thickness in each of the eight sub-regions.

Deviance and diversity of the morphological composition.

The morphological composition of a subvolume is a vector that specifies for reconstructed morphology how often it was placed into the subvolume. A: For each layer and columnar subvolume the deviance of its morphological composition, a measure of its difference to the global composition of the entire layer. It is measured as the negative logarithm of the p-value of a chi-square test of the hypothesis that the composition of the layer/column intersection matches the composition of the layer. B: For each layer and columnar subvolume the diversity of its morphological composition, a measure of the degree to which many different morphologies were placed. It is measured as the binary entropy of the morphological vector. C: For all morphologies of layer 5 neurons their cumulative frequencies. Grey: in layer 5 globally; orange: in the column with the highest layer 5 diversity; blue: in the column with the lowest layer 5 diversity. Cumulative frequencies calculated after sorting by frequency.

Modularity across sub-columns

A: Left: Heat map of modularity values of the networks of connections within each column as parcellated in Figure 7. Right: The same but for ER-controls of each column, thus fixing the number of nodes and connections but modelling random connectivity. B: Modularity is strongly correlated to the volume of the column and the number of neurons in it. The modularity values are much higher than those obtained in columns of the same size with random connectivity. Moreover, for the controls modularity is in fact anti-correlated with volume and number of neurons. C: Modularity values of random subnetworks of the columns with highest modularity (higher than 0.4). For each column, 50% of the neurons were sampled at random 50 times, providing 50 subnetworks for each column. The modularity of these was computed in two ways. Blue restricting the modules of the full column to the sample. Pink: Re-optimizing to find the modules that maximize the modularity value (see Methods). The modularity distributions are virtually the same between the two methods and coincide in their mean with the modularity value of the full column (purple cross).

Basic properties of network connectivity.

A:Top/bottom degree distributions of the local/mid-range networks and their corresponding control models. On the second row additional inserts are provided showing the bi-modality of the out and total degree distributions of the SBM and DBM controls. B: Matrices showing the probability of connection between given pre/post-synpatic m-types. Vertically at the right of each matrix the probability of connection for a fixed pre-synaptic m-type. Horizontally on top of each matrix, the probability of connection for a fixed post-synaptic m-type. C: Left and center, percentage of special nodes in both networks, sources (nodes of in-degree 0) and proper sinks (nodes of out-degree 0 that are not isolated nodes). Right: layer distribution of the special nodes of the mid-range network.

Rich club analysis

A: Depiction of simplices and non-simplices in various dimensions with distinction between the source and target nodes. B:Rich-club analysis, top/bottom rows for the local/mid-range networks. Left column, rich-club curves for the original networks as well as their CM controls. Top left additionally shows the rich-club curves of the DBM controls of the local-networks as well as their corresponding CM controls. Middle column shows normalized rich club curves obtained by dividing the rich-club curves for the original networks by the mean of the rich-club curves of their corresponding controls. Right column shows the simplicial rich club curves for color-coded dimensions ≥ 2.