Multi-scale account of the network structure of macaque visual cortex

Cortical network structure has been extensively characterized at the level of local circuits and in terms of long-range connectivity, but seldom in a manner that integrates both of these scales. Furthermore, while the connectivity of cortex is known to be related to its architecture, this knowledge has not been used to derive a comprehensive cortical connectivity map. In this study, we integrate data on cortical architecture and axonal tracing data into a consistent multi-scale framework of the structure of one hemisphere of macaque vision-related cortex. The connectivity model predicts the connection probability between any two neurons based on their types and locations within areas and layers. Our analysis reveals regularities of cortical structure. We confirm that cortical thickness decays with cell density. A gradual reduction in neuron density together with the relative constancy of the volume density of synapses across cortical areas yields denser connectivity in visual areas more remote from sensory inputs and of lower structural differentiation. Further, we find a systematic relation between laminar patterns on source and target sides of cortical projections, extending previous findings from combined anterograde and retrograde tracing experiments. Going beyond the classical schemes, we statistically assign synapses to target neurons based on anatomical reconstructions, which suggests that layer 4 neurons receive substantial feedback input. Our derived connectivity exhibits a community structure that corresponds more closely with known functional groupings than previous connectivity maps and identifies layer-specific directional differences in cortico-cortical pathways. The resulting network can form the basis for studies relating structure to neural dynamics in mammalian cortex at multiple scales. Electronic supplementary material The online version of this article (10.1007/s00429-017-1554-4) contains supplementary material, which is available to authorized users.

Supplemental material to "Multi-scale account of the network structure of macaque visual cortex" Maximilian Schmidt · Rembrandt Bakker · Claus C. Hilgetag · Markus Diesmann · Sacha J. van Albada the date of receipt and acceptance should be inserted later Black curve in all panels, fit using a beta-binomial model with probit sigmoidal function (eq (7) in manuscript) (a 0 = −0.152, a 1 = −1.534, φ = 0.214). (A): Red dots, real data from Markov et al (2014b). (B): Red dots, random data generated from a betabinomial distribution with fitted parameters a 0 = −0.152, a 1 = −1.534, φ = 0.214. Blue boxes, random data generated from a beta-binomial distribution with parameters a 0 = −0.152, a 1 = −1.534, φ = 0.001. (C): Red dots, real data from Markov et al (2014b). Gray dashed curve, fit using a beta-binomial model with logit sigmoidal function. : r = 0.14, p = 0.69). The thickness data is the same as in Fig. 4. Cortical depth data obtained from F99 surface statistics available through the Caret Software (Van Essen, 2012). Values for each area are averaged across cortical surface and both hemispheres. The data is obtained using the F99 Sulcal depth tool on http://cocomac.g-node.org and can be directly accessed via these two links: http://cocomac.g-node.org/cocomac2/services/f99_sulcal_depth.php? atlas=FV91&shape=Depth-Right&mode=avg&output=tsv&run=1 and http://cocomac.g-node.org/cocomac2/ services/f99_sulcal_depth.php?atlas=FV91&shape=Depth-Left&mode=avg&output=tsv&run=1.     Table S2: Scheme for translating architectural types, overall neuron densities and cortical thicknesses given in Table 4 of Hilgetag et al (2016a) to the modeled areas in the parcellation scheme of Felleman and Van Essen (1991). Entries marked with a star are used to translate the overall neuron density and cortical thickness which are not available in the finer of the two parcellations used by Hilgetag et al (2016a).

Derivation of local average connection probability
The average connection probability is obtained by integrating over all possible positions of the two neurons, with C 0 the connection probability at zero distance and polar coordinates r 1 , r 2 , θ 1 , θ 2 . This can be reduced to a simpler form (Sheng, 1985),

Local connectivity
The indegrees of the microcircuit model (Potjans and Diesmann, 2014) K ′ ij (R) are adapted to the area-specific laminar compositions of the multi-area model with an area-specific factor c A (R), where i, j denote single populations in the 1 mm 2 patch of the cortical area. The total number of synapses local to the patch (type I) is the sum over the projections between all populations of the area: We thus obtain c A (R) by determining N syn,I . To this end, we use retrograde tracing data from Markov et al (2011) consisting of fractions of labeled neurons (F LN ) per area as a result of injections into one area at a time.
The fraction intrinsic to the injected area, F LN i , is approximately equal for all 9 areas where this fraction was determined, with a mean of 0.79. For areas modeled with reduced size, this fraction is smaller because, in that case, synapses of both type I and II contribute to the value of 0.79 (Fig. 4E). We approximate the increasing contribution of type I synapses with the modeled area size as the increase in indegrees averaged over population pairs, where in the last step we use eq. (7). Using N syn,I (R full )/N syn,tot (R full ) = F LN i , we obtain where N syn,tot (R) = ρsynπR 2 D with D the total thickness of the given area. The conversion factor can thus be obtained with We substitute this into eq. (7) for the modeled areas where R = R 0 and obtain the population-specific indegrees for type I synapses:

Processing of CoCoMac data
We use a new release of CoCoMac, in which mappings from brain regions in other nomenclatures were scrutinized to ensure a consistent transfer of connections into the FV91 name space. The CoCoMac database provides information on laminar patterns on the source side from retrograde tracing studies as well as on the target side from anterograde tracing studies. The data was extracted by using the following link, which specifies all search options: http:// cocomac.g-node.org/cocomac2/services/connectivity_matrix.php?dbdate=20141022&AP=AxonalProjections_ FV91&constraint=&origins=&terminals=&square=1&merge=max&laminar=both&format=json&cite=1 Furthermore, we obtained the numbers of confirmative studies for each area-level connection with the following link: http://cocomac.g-node.org/cocomac2/services/connectivity_matrix.php?dbdate=20141022&AP=AxonalProjections_ FV91&constraint=&origins=&terminals=&square=1&merge=count&laminar=off&format=json&cite=1 To process these data, we applied the following steps: -A connection is assumed to exist if there is at least one confirmative study reporting it.
-A connection from layer 2/3 is modeled if CoCoMac indicates a connection from either or both of layers 2 and 3. -In the database, some layers carry an 'X' indicating a connection of unknown strength. We interpret these as '2' (corresponding to medium connection strength). -We take connection strengths in CoCoMac to represent numbers of synapses in orders of magnitude, i.e., the relative number of synapses N ν syn in layer ν of area A with connection strength α(ν) is computed as N ν syn = 10 α(v) / v ′ ∈A 10 α(ν ′ ) .

Mapping of synapse to cell-body locations
Detailed calculation in section Materials and Methods. The numbers are listed in Table S11. We map the 17 different cell types of Binzegger et al (2004) to the 8 cortical populations considered in our network using Table  S10.   Table S11: Conditional probabilities P(i|scc ∈ v) for the target neuron to belong to population i if a corticocortical synapse scc is located in layer v, computed with eq. (11) applied to the data set of Binzegger et al (2004). Empty cells signal zero probabilities.

Derivation of the cortico-cortical connectivity
The number of cortico-cortical synapses from excitatory population j in area B to population i in area A is calculated as Here, S = 2/3E and I = {5E,6E} respectively denote the supragranular and infragranular excitatory populations, v i is the layer containing population i, and α(v) is the layer-specific labeling density estimate from CoCoMac. Z i is a factor which takes into account that cortico-cortical feedback connections preferentially target excitatory rather than inhibitory neurons (Johnson and Burkhalter, 1996;Anderson et al, 2011). For each feedback connection in the model, we redistribute the synapses across the excitatory and inhibitory target populations and determine Z i such that 93 % of synapses in each cortico-cortical projection are established on excitatory neurons.

Validation of network characteristics
To validate our findings, we tested the main characteristics of the network model against alterations of three main heuristics used in its construction:

Exponential decay of connection densities
To test the robustness of the area-level connectivity, we excluded a random subset of the F LN data (Markov et al, 2014a) and repeated the fit of the exponential decay of connection densities (Fig. 4C) for these reduced datasets ( Figure S4A). We then used the altered fitted parameters to estimate both unknown and excluded data and proceeded with the network construction. To test if this had a significant effect on the network, we repeated the clustering analysis of Fig. 7 for 50 different trials ( Figure S4B). The community structure is generally stable if 5% or 10% of the data are ignored. Ventral and dorsal stream areas are assigned to separate large clusters. Frontal areas are consistently clustered together. The two lower visual areas V1 and V2 are clustered together except in two trials, where the data on their connectivity is ignored so that it is estimated using the EDR rule, resulting in a lower connection density than suggested by the data. Areas STPp, STPa, PITd, and MSTd change clusters frequently. However, overall, the community structure proves to be robust against leaving out small portions of the underlying experimental data. If half of the experimental data is excluded, the resulting fits for the exponential decay of connection densities show clear deviations. Still, ventral and dorsal stream areas as well as frontal areas are well separated into clusters, but the cluster assignment of the remaining areas varies across trials. As an additional validation, we used the full F LN data to fit the exponential decay of connection densities and then introduced a random offset uniformly distributed within ±∆ for the logarithm of each estimated F LN to reflect the spread of experimental data around the fitted exponential decay of connection densities. We chose ∆ = 1.55 to be the median of the absolute deviation of the logarithm of experimental data and exponential fit. This leads to a degree of variation between trials that is comparable to the case of pruning 10% of the data. Areas PITd and MSTd frequently and FST, STPa, STPp, VP, and VOT to a lesser extent change clusters.
In conclusion, the community structure of the network does not rely on single data points, but the data leave open the possibility of different cluster assignments especially of areas PITd and MSTd.   shows the regression using the full dataset (as shown in Fig. 4C). b Visualization of community structure in the resulting network connectivity using the fitted parameters for the exponential decay of connection densities on randomly drawn subsets of the data. Clusters are determined by applying the map equation method (Rosvall et al, 2009) on the outdegrees, as in Fig. 7. The colors and cluster names are chosen as in Fig. 7. For comparison, the first column of each matrix shows the clustering of the network based on the full dataset as shown in Fig. 7 Threshold for hierarchical categorization based on SLN In Fig. 5B Figure S5 shows that the chosen threshold has little influence on the laminar connectivity patterns in the network. Both the target patterns and the resulting population-specific connectivity patterns show only marginal differences. The laminar distribution of shortest paths reveals small differences depending on the chosen threshold: For a broader range of SLN corresponding to lateral connections (middle column), feedforward paths solely originate in layer 2/3. For a narrower range of SLN corresponding to lateral connections (right column), more feedback connections are formed in layer 2/3. The qualitative picture remains unchanged: There are clear differences between feedforward and feedback connections, and lateral connections are more similar to feedforward connections than to feedback connections in their start-end patterns. Figure S5: Laminar connectivity patterns and shortest paths are robust against altering SLN threshold. Top panels: Same display as in Fig. 5B. Laminar target patterns of synapse locations in relation to the fraction of supragranular labeled neurons (SLN ) of the source pattern. Target patterns are taken from the CoCoMac database (Felleman and Van Essen, 1991;Barnes and Pandya, 1992;Suzuki and Amaral, 1994;Morel and Bullier, 1990;Perkel et al, 1986;Seltzer and Pandya, 1994) and SLN data from Markov et al (2014b) mapped to the FV91 scheme. Middle panels: Same display as in Fig. 5D. Laminar patterns of cortico-cortical connections in the feedback, lateral, and feedforward direction, measured as the indegree of the population pairs divided by the sum of indegrees over all pairs, and then averaged across area pairs with the respective connection type  Figure S5: Laminar connectivity patterns and shortest paths are robust against altering SLN threshold.
Sigmoidal relation of SLN with log cell density ratios To test the robustness of our path analysis (Fig. 8) against the sigmoidal relation of SLN versus logarithmic cell density ratios, we used a similar procedure as for the exponential decay of connection densities: We ignored random subsets of the data (5%, 10% and 50%) from Markov et al (2014b) for the sigmoidal fit and then estimated missing data as well as the ignored data based on this modified fit. The fit of the sigmoid is only slightly influenced by excluding subsets of the data ( Figure S6A). The resulting patterns of shortest paths also show little variation across trials. For 50% of data pruning, there is visible variation in the results, but the qualitative picture still holds. This shows that the heuristic of the SLN sigmoid delivers a robust basis for laminar connection patterns in the network model. Figure S6: Laminar connectivity patterns and shortest paths are robust against altering SLN threshold. We randomly choose 5% (top part), 10% (middle part), and 50% (bottom part) of the SLN data of Markov et al (2014b) and repeat this for 50 trials. Left part Same display as in Fig. 5A. Fraction of source neurons in supragranular layers (SLN ) vs. logarithmized ratio of the overall neuron densities of the two areas. SLN from Markov et al (2014b), neuron densities from Hilgetag et al (2016b). Black curve, fit using a betabinomial model eq. (1) on a randomly drawn subset of data points (blue dots) and ignoring the remaining data points (red dots). Overlapping gray lines show linear regression lines for different trials and the dashed line (overlapping with black line) shows the regression using the full dataset (as shown in Fig. 5C). Right part: Same display as Fig. 8. Top panels Population-specific patterns of shortest paths between directly connected pairs of areas categorized according to their hierarchical relation as defined by the SLN . Arrow thickness indicates the occurrence of the particular pattern. Middle panels Population-specific patterns of shortest paths between all pairs of areas categorized according to the difference between their architectural types. Arrow thickness indicates the relative occurrence of the particular pattern. In top and middle panels, gray arrows show the mean value across trials and the underlying red arrows show the mean plus one standard deviation (partly not visible because of the small variation across trials). Lower panels Occurrence of population patterns in areas that appear in the intermediate stage in the shortest path between two areas. Error bars indicate standard deviation across trials (partly not visible because of small values)  High to low type Horizontal Low to high type Figure S6: Laminar patterns of shortest paths are robust against using only a subset of SLN data.

Path analysis using only connections with experimental SLN
To test whether the predicted laminar pattern of shortest paths depend crucially on the SLN predicted with the sigmoid, we repeated the path analysis shown in Fig. 8 but took only connections into account for which SLN data were provided by the dataset of Markov et al (2014b). This amounts to approx. 23% of the connections. Figure  S7 shows that the laminar patterns of shortest paths largely match those determined from data plus estimates. The most obvious differences occur in low-to-high-type connections with a smaller proportion of paths starting in supragranular layers, and in horizontal connections with a larger proportion of paths starting in infragranular layers. Furthermore, layer 6 becomes more prominent compared to layer 5 in the origin of horizontal, feedback, and low-to-high-type paths. This does not reflect a shortcoming of the sigmoidal fit, but occurs because for this subset of connections there are relatively more source patterns stored in CoCoMac where L6 sends more connections than L5 than if all connections are taken into account. The overall characteristics, however, are still preserved, showing that the laminar patterns predicted by the model are not artificially introduced by the estimation of SLN through the sigmoidal relationship (Fig. 5A) Figure S7: Laminar patterns of shortest paths are not introduced by estimated SLN . Same display as in Fig. 8 but taking only connections with available experimental data on fractions of supragranular labeled neurons (SLN ). a Population-specific patterns of shortest paths between directly connected pairs of areas categorized according to their hierarchical relation as defined by the SLN . Arrow thickness indicates the relative occurrence of the particular pattern. b Population-specific patterns of shortest paths between all pairs of areas categorized according to the difference between their architectural types. Arrow thickness indicates the occurrence of the particular pattern.