Mammals show a wide range of brain sizes, reflecting adaptation to diverse habitats. Comparing interareal cortical networks across brains of different sizes and mammalian orders provides robust information on evolutionarily preserved features and species-specific processing modalities. However, these networks are spatially embedded, directed, and weighted, making comparisons challenging. Using tract tracing data from macaque and mouse, we show the existence of a general organizational principle based on an exponential distance rule (EDR) and cortical geometry, enabling network comparisons within the same model framework. These comparisons reveal the existence of network invariants between mouse and macaque, exemplified in graph motif profiles and connection similarity indices, but also significant differences, such as fractionally smaller and much weaker long-distance connections in the macaque than in mouse. The latter lends credence to the prediction that long-distance cortico-cortical connections could be very weak in the much-expanded human cortex, implying an increased susceptibility to disconnection syndromes such as Alzheimer disease and schizophrenia. Finally, our data from tracer experiments involving only gray matter connections in the primary visual areas of both species show that an EDR holds at local scales as well (within 1. 5 mm), supporting the hypothesis that it is a universally valid property across all scales and, possibly, across the mammalian class.
Understanding brain networks is arguably one of the major challenges of the 21st century [1]. The mammalian cortex is an extraordinary computational device, and analysis of its network properties with 107–1010 neurons and 1011–1015 synaptic connections is still largely unresolved. In the brain, activity of a single neuron encodes relatively little information; instead, that is achieved via population coding, through spatially distributed temporal activity patterns of cell assemblies. This contrasts with packet-switching information technology (IT) networks, which encode information directly into the packets and the network merely ensures routing between any two nodes. Since the spatiotemporal activity of cell populations is strongly determined by their connectivity and physical layout, cortical network structure and its spatial embedding play a significant role in the brain’s processing algorithm, in sharp contrast with IT networks. A purely bottom-up approach to deriving global brain function from local circuitry is currently intractable [2]. In contrast, a meso-scale approach is more feasible, focusing on the network of interactions between the elements of a mosaic of distinct areas representing the loci of function-specific computation (visual, auditory, somatosensory, motor, etc.). As the mammalian brain is shaped by evolution, morphological and areal network level inter-species comparisons will help identify those features that are conserved across species from those that are species-specific. This will lead to a better understanding of network structural properties and provide valuable clues to the evolution of brain function [3]. However, progress in this direction has been hindered due to the absence of (i) the necessary data to address the physical properties of the network between areas and (ii) adequate theoretical network comparison methods. Published connectivity maps using consistent interareal tract tracing studies, first in the macaque [4] and more recently in the mouse [5,6], allow consideration of the network as a directed, spatially embedded and weighted graph (weights representing neuronal connection densities projecting between areas). The absence of full homology between the nodes (areas) and edges (projections) of the networks of the two species makes it difficult to determine commonalities and similarities between them. However, if generic, global organizational principles exist (constraining the adaptation and growth of cortical connections in similar ways), then we expect to see similarities at the statistical level between the network features in the two species. Here we show that the cortical networks in the macaque and the mouse in fact do exhibit a common organizational principle despite their very different evolutionary trajectories and large differences in brain size. Supplemented by partial tract tracing data in the microcebus (the mouse lemur) we suggest that this principle and the associated network model is a universal determinant of the interareal network across mammals, allowing tentative predictions for the human brain. Expansion of the cerebral cortex is accompanied by an increase in the proportion of white matter relative to brain size [7–10]. However, this increase is not rapid enough to maintain a constant neuronal connection density (defined as the fraction of neuron-to-neuron connections compared to all possible ones). Thus, an increase in brain size is expected to result in a reduction in the long-distance connectedness of cortical areas [11–14]. The reduction of the fraction of connections with cortical expansion and the minimization of the metabolic costs are important design features of the cortex [4,15–26]. One can hypothesize that this wire minimization constitutes a critical constraint for the optimal placement of areas in the cortex, serving to increase communication efficiency in larger brains [11,27–29], and is supported by recent evidence suggesting reduction of long-distance connectivity with increases in brain size [28]. Recent retrograde tract tracing data in macaque [30] provides supporting evidence precisely of such a wiring constraint, in the form of an exponential decay of the wiring probability p (d) with projection distance d: p (d) ~e−λd, with a decay length (~1/λ) that is short relative to hemispheric dimensions (in the macaque λ ≅ 0. 19 mm−1, corresponding to a decay length of 1λ≅5. 2 mm). A simple way to think of the decay length 1λ is that every increase by 1λ in projection length leads to a decrease in the number of projections by a factor of 1e≅0. 37 (i. e. , 37%). Note that using the base of the natural logarithm is convenient, as in this case 1λ is equal to the average projection length, providing a simple, intuitive interpretation. We refer to this decay property of connection density with distance as the Exponential Distance Rule (EDR). Retrograde labeling using fluorescent tracers (see Materials and Methods section) is an accurate labeling method that reveals all incoming connections j→i to an injected (target) area i by labeling the cell bodies of the neurons in source area j whose axons make connections in area i. Importantly, there is no transneuronal labeling, so the retrograde labeling method used yields only one-step incoming connections to the injected nodes of the network. Note that the EDR is purely a property of the distribution of the physical lengths of individual axons, without regard to any network topological structure. The EDR states that there are many fewer long-range axons than short ones and quantifies this: the number of axons of length d that we find in the cortex is proportional to e−λd. In general, to experimentally establish the EDR, we do not need to work with brain areas as nodes of a network; we only need to be able to count neurons and measure the corresponding axon lengths. In this sense, the EDR is a more basic and general property than the description of cortical connectivity as a network at some coarse-grained (e. g. , mesoscale) level. Once the level of description is defined (e. g. , areal), the network properties are, however, consequences of the distribution of the axonal lengths connecting the vertices. Since connectomes are embedded in physical space, the EDR property effectively constrains the topological structures that connectomes can form across different levels, ranging from the single neuron to the areal level [31]. In addition to the discovery of the EDR in the macaque, the consistency and completeness of this tract tracing data [32] has led to a deeper insight into the interareal network properties of the macaque cortex [30,33]: it revealed a much denser (ρ = 0. 66) interareal cortical graph than previously reported (network density is defined as ρ = MN (N−1), where N is the number of areas and M is the number of connected ordered area pairs, see glossary). High density graphs have low specificity at the binary level (areas connected or not), so that what distinguishes one area from another is the particular combination of areas it is connected to, combined with the weights of the connections, i. e. , their connectivity profile or fingerprint [33–36]. Because the range of weights spans many orders of magnitude (five in the macaque), the specificity of individual connectivity profiles is actually very high [5,30,37].
To what extent does the EDR, as a connectivity constraint, determine the properties of the interareal network? To address this issue, one needs (i) a family of EDR-based network models and (ii) a method of comparison between the model-generated networks and the experimental data network. The exponential decay rule ~e−λd in the macaque was obtained from collating all the labeled neurons (over 6. 4 million) following tract tracing experiments in different areas and constructing an interareal distance matrix, the latter estimated as the distances between the area barycenters through the white matter (WM), along the shortest paths. Here axonal p (d) should be interpreted as an average property (see Fig 1A), the probability that an axonal bundle projects to a distance d, independently of the specific functional nature of the areas. At this level of description, the strength of the connection between areas, expressed as the fraction of labeled neurons (FLN), depends uniquely on their geometrical separation. Thus, the network is viewed as a spatial, directed, and weighted graph dependent on the matrix D = {dij} of interareal distances dij. We emphasize here that the EDR arises from the estimated probability distribution of axon lengths. Although the strength-distance relation is consistent with the EDR, the probability distribution of axons lengths provides a more compelling demonstration of the property and leads naturally to the parametric EDR model described below. The probability density function, q (d), of the distances in the matrix D is typically a unimodal distribution (Fig 1B), which, when combined with the exponential decay p (d), leads to a log-normal distribution of edge weights, confirmed by the empirical FLN data [4,39]. The EDR distribution with the corresponding distance matrix D in a given brain naturally defines a parametric family of random graphs, called EDR random graphs (Fig 1C), parameterized by the decay rate λ. For these model graphs we make the choice p (d) = λe−λd, where now λ is the (only) model parameter. To distinguish the decay rate parameters in these models from the experimentally measured ones, we denote the latter as λexp, e. g. , for macaque λexpmac = 0. 19 mm−1. We also employ, as a null model, the constant distance rule (CDR) family of random graphs, where there is no dependence of connection probability on distance, corresponding to the λ→0 limit, i. e. , to the choice p (d) = const. The EDR family of random graphs is defined via a simple algorithm [4] in the spirit of the Maximum Entropy Principle, i. e. , it is based only on the given information (p (d) and D), while all else is uniformly random. The algorithm proceeds as follows: First, we randomly draw a connection length d from the distribution p (d). Second, we choose uniformly at random an area pair whose separation distance in the matrix D falls in the same distance bin as d, according to some binning criterion (bin sizes used in this study were typically 5 mm for the macaque and 0. 4 mm for the mouse) and finally, insert a randomly oriented connection between them. Multiple connections between the same area pair in the same direction generate the weights for the directed edges with a log-normal distribution. These steps are then iterated until the graph density in the model reaches the observed value in the experimental network. We denote the data network obtained from the experiments by Gexp (e. g. , for the macaque we use Gexpmac, and for the mouse Gexpmac). Our goal is to compare the properties of the EDR model networks with the properties of Gexp. Since the model networks are only based on distance-dependent connection probabilities, one cannot expect perfect agreement (edge-by-edge) with the biological connectivity graph Gexp, however, if the distance rule is a strong determinant of the interareal network, the model graphs should be statistically similar to Gexp. The comparisons are performed via parameter matching of network properties [4]: for a given network property P, the interareal distance matrix D and parameter λ is used to generate a large ensemble of EDR graphs GEDR (λ). By varying λ we determine the value λP via minimizing the deviation |P (Gexp) − 〈P (λ) 〉|, with the average 〈∙〉 taken over at least 103 EDR graph realizations from GEDR (λ). Thus the model parameter is determined so that the average of P in the model is as close as possible with the value of P observed in the data network. We then compare the fitted value λP with λexp, the decay rate obtained directly from the experiments. If the two are close, then the EDR is a strong determinant for the measure P of the cortical network. Thus, the extent a particular measure in the EDR model and in the data network agree, i. e. , |P (Gexp) − 〈P (λP) 〉| with respect to the same comparison with the CDR model, i. e. , with |P (Gexp) − 〈P (λ = 0) 〉|, expresses the degree to which the EDR influences that particular measure in the cortical network. This analysis is repeated with several local and global network measures. The more measures for which there is an agreement between λP and λexp, the stronger the effect of the EDR in shaping the interareal network. This method also has the added advantage of identifying those network properties that are not well described by the EDR, and thus, based on the nature of these measures, providing us with clues for additional network mechanisms. In the macaque, the EDR model predicts very well many local, global and weighted network properties of the interareal network (see [4] for details), and thus it is a strong determinant for the large-scale network organization of the macaque cortex. It also captures its pronounced core-periphery organization (i. e. , a densely connected set of areas—core, with feedback and feedforward links to/from a more loosely connected set of peripheral areas), with the core strongly dominated by associative areas [4,40]. The EDR network model of cortical connectivity represents a radical departure from previous, purely topological models of cortical networks, which do not take into account their physical, i. e. , weighted and spatially embedded nature, and this has now been well documented in the recent literature [41,42]. The spatial clustering and geometrical positioning of the nodes in the EDR model in the macaque is observed to strongly echo the functional layout of the cortex as revealed by numerous physiological and anatomical studies [36,43]. To determine whether a similar description is valid for the mouse (Mus musculus) cortex, we first conducted retrograde tracer experiments in the mouse neocortex (S2 and S3 Figs), in order to determine the projection length distribution p (d), which, indeed, shows a clear exponential decay (Fig 2A). The decay rate, λexpmus, was determined from an exponential fit as λexpmus = 0. 78 mm−1, with a 95% confidence interval of (0. 72,0. 83) (see Fig 2A, and inset). This exponential decay is to be compared with the same distribution for the macaque from Fig 2B in Ref [4] (see Table 1 for the λ parameter estimates). The distance matrix Dmus was determined from flattened cortex measurements. The corresponding distance distribution q (d) is unimodal, as shown in Fig 2B, which is to be compared with the same distribution for the macaque from Fig 2C in reference [4] (the consequences of the differences in these distributions are discussed in more detail in the section “Functional Layout in Terms of Spatial Clustering of Cortical Areas”). We then applied the network analysis described in the section above to the largest available edge-complete graph (the status of connectivity between all pairs of nodes is known) of 33 areas in one hemisphere of the mouse neocortex [5,6], denoted by Gexpmus from here on (S1 Fig). This mouse dataset contains 719 directed pathways and has an interareal network density of ρmus = 0. 68, similar to that reported in the macaque. Fig 3 shows the proximity of λPmus obtained using the parameter matching method to the decay rate λexpmus for several network measures including the number of area pairs connected uni- (M1) or bidirectionally (M2), 3-motifs, clique distributions and the second largest eigenvalue of the symmetrized form AAT of the adjacency matrix A. These measures have been selected in part because they probe graph properties from local to global scales, and are of varying complexity. Additionally, these measures (see glossary for definitions), and in particular the deviations from their values in random graphs carry functional significance: unidirectionally connected areas depict an asymmetric role in information processing (driver versus driven nodes), the 3-motifs have extensively been studied as building blocks of functional organization in complex networks [38], cliques identify maximally connected network regions usually representing activity-specific strongly correlated communities or clusters, and the second largest eigenvalue is related to the rate of spreading processes (e. g. , epidemics or information) on the network [4]. The different comparisons and fits based on these measures are highly consistent and indicate λ to be in the range 0. 78–0. 93 mm−1 (purple vertical band in Fig 3A–3D). The broader range in mouse of |λP−λexpmus| compared to macaque might be due in part to the fact that the mouse connectivity matrix was generated by anterograde tracing [5,6]. Other statistical network properties, such as degree distributions, are likewise well captured by the EDR network model with λ = 0. 78 mm−1, (see S4 Fig). A clique (see glossary in S1 Text) is a complete subgraph of a network, i. e. , it carries the maximum number of possible edges between its nodes. In dense graphs (thus with many cliques) the size (number of nodes) distribution of the cliques provides insight into the network’s heterogeneity [4]. The largest cliques in dense graphs can be used to define the cortical network core [4,40]. As in macaque [4,40], the clique distribution analysis in the mouse (Fig 4A) reveals a distinct core-periphery structure. The mouse connectome, Gexpmus, includes a dense core of 12 nodes organized into the two largest cliques each of size 11, plus a periphery of 21 nodes. There are a total of Mcc = 131 links within the core, Mcp = 190 links from the core to the periphery, Mpc = 170 from periphery to core, and Mpp = 228 links within the periphery. Densities for the mouse are the following: core 99% (versus 92% in macaque), periphery 54% (versus 49%) and the links between the core and periphery, 71% (versus 54%). The likelihood of a core having 12 nodes in a random graph on 33 nodes with the same density ρ = 0. 681 as in Gexpmus is vanishingly small: (3312) (1321) p131 (1−p) 1 = 2. 07×10−12 (versus 10−17 in the macaque). To see how well the EDR model reproduces the clique distribution, we define a scalar deviation measure σcl (λ) between the clique-size distributions in the data and the EDR model as the root mean square (RMS) of the clique-count log-ratios. The best agreement between the two distributions is achieved at λclmus = 0. 93 mm−1 (Fig 3C) and the clique distributions in the model and data are rather close at this value (Fig 4A). Anatomically, the cortical core in the mouse shows significant differences with that previously reported in macaque, the most striking being that the mouse core includes portions of primary somatosensory cortex (SSp-ll and SSp-tr) and primary motor cortex (MOp) (Fig 4B). While additional injections may well expand the core membership in macaque, primary areas in the macaque core are extremely unlikely, given the rarity of connections linking primary areas [30]. This contrasts with the mouse where the inter-primary area subgraph has a density of over 80% [5,6]. In agreement with the presence of primary areas in the mouse core, the two-dimensional map of the flattened cortex (Fig 4D) shows that the mouse cortical core might be spatially more widespread across brain regions compared to that of the macaque, where the core appears concentrated in frontal and parietal areas [4]. Note that in both mouse and macaque, the core areas have overall, higher in-degrees than non-core areas (Fig 4B for the mouse). The wider spatial spread of the mouse compared to the macaque core may reflect the relative expansion in primates of higher-level association cortex with respect to the primary areas [3]. These differences in the cortical core of the mouse and macaque need to be considered in light of the proposal that in primates at least, the core is related to cognitive architectures such as the global workspace, thought to be involved in consciousness [40,44]. Network motifs refer to the different possible connectivity patterns of a small, fixed number of nodes. For example, in the 33-node mouse cortical network there are (333) = 5456 triplets of nodes, each of which has one of the 16 connectivity patterns shown in Fig 5A. Three-node motifs have been proposed as the building blocks of network circuits and their pattern of variation in frequency to reflect functional properties of the networks [38]. For instance, motif 10 (oriented 3-cycle) is significantly under-represented in the cortex, while motif 3 is significantly over-represented (lone, bidirectional link) when compared to a random network, in both species (see Fig 5D). As in macaque [4], the mouse EDR predicts the observed motif frequency distributions in this species significantly better than does the CDR (Fig 5B). Despite the marked quantitative differences in motif distributions between mouse and macaque (Fig 5C), there could, however, be qualitative similarities. Testing this requires comparing the observed motif profiles to that of a randomized null model [38] consisting of an ensemble of random networks having the same degree sequence as the data. Graphs were uniformly sampled from this ensemble by repeatedly rewiring edges [45]. Fig 5D shows how the motif counts of the empirical connectomes differ from such a randomized null model. As similar patterns are observed in both species, these findings suggest that they are part of the same class of large-scale networks with similar architectural and functional constraints. Repeating this analysis for networks generated by the EDR model (see Fig 6) we find a remarkable similarity in the motif profiles not just between the networks of the real data network and the EDR model but also between the model networks of the two species (Fig 6). This confirms the existence of a common network architectural invariant in these two species. This is unexpected, insofar as the decay rates and the distance matrices are very different between the two cortices. Since the motif profiles are binary measures, these findings indicate graph structural similarity between the two brains. In order to test whether there are significant similarities in the large-scale connectomes of the two species beyond the constraints imposed by the EDR, we used the EDR model as a null model [46]. S5 Fig shows that motif counts continue to look similar between the mouse and macaque, although their similarity is now less pronounced. To further probe the wiring similarity between the mouse and macaque connectomes, we next study the connectivity similarity profile measure. Elsewhere we have demonstrated that a quantitative measure of the similarity of the connectivity profiles of target cortical areas decreases in a regular fashion with increasing distance between them, i. e. , the closer two target areas are, the more their source areas overlap [4,33]. We have also shown that changes in similarity reflect the functional layout of the cortex [33], and thus it is natural to compare the behavior of this measure between the mouse and macaque. A similarity index [4] can be defined for both incoming (in-link similarity) and outgoing (out-link similarity) connections. In order to compare macaque and mouse similarity indices, we focus here on the incoming connections, as those are the ones fully specified for all the injected areas in the macaque dataset. Next, we analyze the similarity between the connectivity profiles, for all possible target area pairs. The in-link similarity index for any two (target) areas is a measure describing the extent to which both targets receive/or do not receive in-links from the same source areas, compared to a fully randomized state of the network (see Materials and Methods section for details). Fig 7 shows the distribution of in-link similarity indices as function of the distance between all area pairs in the mouse (Fig 7A) and the macaque (Fig 7B). In both species, in-link similarity decreases with increasing distance between the area pairs, i. e. , areas that are further apart on the cortical sheet have increasingly dissimilar in-link connectivity profiles on average, while the opposite is true for areas that are closer to one another. The colored regions in both Fig 7A and 7B are the probability densities of in-link similarity indices generated by the corresponding EDR models, with red corresponding to higher, and blue to lower probabilities; in both cases the EDR model captures the average behavior rather well. In order to compare distance-dependent quantities between brains of very different sizes, all distances are rescaled by the average interareal distance in each species (〈d〉mus = 4. 54 mm and 〈d〉mac = 26. 35 mm). Interestingly, as the largest distances are dmaxmus = 10. 1 mm and dmaxmac = 58. 2 mm respectively, this fits both brains onto the same adimensional template, as dmaxmusdmus = 2. 22 mm and dmaxmacdmac = 2. 21 mm. Fig 7C shows the corresponding distribution of adimensional distances q (d/〈d〉). When plotting the in-link similarity indices against the rescaled distances (Fig 7D), we find a remarkable overlap between the clouds of points in the two species. This is rather surprising given the fact that they have very different decay rates λ. They also have rather different interareal distance matrices as the macaque cortex is folded, resulting in it having a more peaked distance distribution than the mouse (Fig 7C). Fig 8 shows the sensitivity of the in-link similarity indices using the EDR models in both the mouse (panels 8A–8D) and the macaque (8E–8H). For a given distance matrix, the point clouds are observed to rotate as a function of λ in both species, and hence there is no a priori reason for the overlap in Fig 7D. This overlap, however, is an indication of the existence of a network architectural invariant, present in both species, also reflected in the motif profiles discussed earlier. Further explanation for the significant, overall overlap between the similarity distributions for the two species is provided in the Discussion section. With the help of the common adimensional template defined above we now discuss species-specific characteristics in our comparison of cortical networks. The EDR decay p (d) can simply be recast in terms of adimensional distances, by writing (d) ~e−λd = e−γ dd, where γ = λ〈d〉 is the adimensional (or normalized) decay rate. Accordingly, γmus = 0. 78 × 4. 54 = 3. 54 and γmac = 0. 19 × 26. 35 = 5, showing that on the common template, the mouse has a shallower connectivity decay than the macaque. The distribution of distances in the mouse is broader compared to the macaque (Fig 7C), which when coupled with the shallower connectivity decay contributes to the mouse cortex experiencing a less constraining effect of the EDR than does the macaque. This difference in the EDR between the two species explains some of the differences in the functional layout of the cortex in mouse and macaque. In Fig 9A and 9B we show the same similarity indices for all area pairs as before but also indicate which area pairs are connected (black circles) and which are not (white circles) and provide smooth estimates (colored regions) of connection probability as a function of similarity and adimensional distance. Comparing Fig 9A and 9B we see that in macaque, spatially clustered, presumably functionally related neighboring areas are heavily interconnected and share similar connectivity profiles, while more distant areas show weaker probability of connectivity and similarity index. This relationship between probability of connectivity, spatial separation and similarity is, however, weaker in the mouse. In both species, connection probability changes as a function of distance. Fig 9C and 9D show how this relationship differs in the two species. Consistent with a steeper EDR in the macaque, neighboring areas exhibit 100% connectedness, and the probability of connections (density) decreases smoothly and consistently with distance to around 10% density at maximum distances [47]. This contrasts with the mouse, in which neighboring areas do not quite reach densities of 100% and widely separated areas have densities in the region of 50% to 80% (Fig 9C). Hence, these results show that compared to macaque, in the mouse, widely separated areas are more likely to be interconnected. These differences in the probability of being connected as a function of distance between the two species appear highly significant (smooth curves in Fig 9C and 9D). Numerous studies point to the cost of long-distance connections as an inherent design challenge associated with differences in brain size [48]. One way to define total wire length is: Λ = Σi, jAijDij, where A denotes the binary adjacency matrix and D is the interareal distance matrix. Yoked permutations of the rows and columns of the adjacency matrix reassign the distances to each pair of areas while maintaining the connectivity unchanged. As in macaque [4], the total wire length of the mouse inter-areal network is significantly shorter than a random permutation of the areas (S6 Fig). Simulated annealing methods [4] showed that optimization of area placement can lead to a 12% reduction in total wire length in the mouse, significantly higher than the 5% reduction obtained in macaque [4]. Next, we address the strength of connections with the expectation that long-range connection strengths (expressed as FLNs) would decrease in the larger brain. Due to the EDR, the FLN clearly decreases with distance. Distinguishing interareal association and canonical connections allows an improved understanding of the effect of distance on connection weight (for definition of associative and canonical connections see [3]) (Fig 9E and 9F). This suggests that the decline in FLN is steeper in canonical cortex compared to association cortex, so that the long-distance association cortex connections are one to two orders of magnitude stronger than the connections between canonical cortex areas with the same separation (see [3]). However, the results suggest that the decline in weight with distance is steeper in the macaque compared to the mouse. Together these findings show that compared to the larger macaque cortex, in the smaller mouse brain long-distance binary interareal connections are marginally more numerous. By contrast there is a highly significant increase in the weight of the long-distance connections in the mouse, and this species difference is more pronounced in the projections of association than in the canonical connections of the primary areas.
The present meso-scale network investigation of the neocortex, with appropriate network comparisons, provides detailed information on a common organizational principle that explains numerous network features in two widely separated species, with distinct evolutionary histories. Based on phylogenic considerations, and the fact that evolution is essentially a tinkerer [49], one expects to find evolutionarily preserved features embedded in these networks, i. e. , architectural invariants. Evolutionarily preserved features, however, often are expected to manifest themselves as organizational principles tied to biophysical constraints. The success of the mammalian class includes adaptation to diverse habitats and lifestyles, which is in part attributed to the behavioral flexibility ensured by the neocortex [50]. The modulation of corticogenesis [51] has led to extant mammals exhibiting a five-orders of magnitude range of brain size [52], going from small-brained mammals that include miniaturization of ancestral forms to the expansion and additional arealization that characterize primates, especially humans. The present results suggest that the EDR plays a key role across the mammalian order to optimize the layout of the inter-areal cortical network allowing larger-brained animals to maintain communication efficiencies combined with increased neuron numbers. Our results indicate that the EDR and the associated network model provide a unifying framework to capture common network properties but also some of the differences across the mammalian branch and thus allow network comparisons between species. The EDR decay rate λ and cortical geometry (interareal distances) significantly impact on the structural heterogeneity of the cortical network with important consequences for the general functional layout and core-periphery structure, that we speculate, could be involved in higher cognitive processes [40]. The limitation of the EDR model stems from the fact that the EDR describes an overall, or average property. At this level, without additional determinants, it should not be used as a generative model of individual connections as we have emphasized elsewhere [4,40]. If we plot the decay of the probability of connections for several target areas, as shown in Fig 10A for the mouse, we find significant variability. The black line in Fig 10A, the average decay, is the same as that in Fig 2A. The fluctuations for a given target, however, are not noise, but rather they are part of a signal. This we illustrate in macaque: Fig 10B shows the consistency of fluctuations following repeat injections in area V1, in five different individuals. There are numerous factors that one might need to take into account to better understand this variability. For example one may need to consider the observed systematic variation in neuron numbers across the cortex [53,54], the anisotropy of axon outgrowth distributions [55] and possibly diverse developmental factors [56]. Overall, however, these considerations emphasize that the EDR network serves as a framework, upon which other details are imposed. Note that in order to assess the ability of the EDR model (or any connectome model) to reproduce properties of empirical network data, it is crucial that the data is as edge-complete as possible, i. e. , that the connectivity between any two nodes is known. Otherwise the lack of fit between model and data cannot be used to discard the model [57]. This holds for two reasons: (a) the EDR network model produces complete connectivity information between its nodes, it cannot generate “untested” connections, by default, and (b) many network measures can be sensitive to the absence or presence of an even a small fraction of connections in the network. It is also important to emphasize the roles of cortical geometry [58] and that of areal segmentation in shaping the network properties of the connectome. Since the connection probability depends on distance, network properties are influenced by the relative proximity of areas. In turn, the strength of connections between functionally defined areas correlate with the amount of signaling activity between them and therefore with their functional roles within the information processing hierarchy in the brain. Ad-hoc segmentations, however, will generate ad-hoc distance matrices for the EDR model, and accordingly, the model networks would no longer be interpretable from a functional circuitry point of view, and in this sense, it is important to use optimally defined functional parcellation of the cortex. Our comparative analysis of motifs and connectivity similarity indices demonstrates the existence of network architectural invariants, which in turn imply that the EDR parameter λ and areal positioning (geometry) are not independent parameters: while both change during evolution, the changes are orchestrated in such a way as to ensure that certain network/circuitry properties are preserved. As argued in the introduction, the network, i. e. , the graph connectivity (form) must play a significant role in the information processing algorithm itself (function), and thus these network invariants are a reflection of common processing dynamics in the cortex. Our use of a normalized or adimensional distances facilitates comparisons across brains of different sizes. Fig 10C shows directly the fingerprint of such universal principles in neocortical organization: it shows the connection probability decay on the adimensional template brain from a common target area (area V1) in macaque (data from reference [30]), mouse as well as microcebus. At short to medium distances where the vast majority of neurons are located, decays are identical, but are observed to change in a species dependent fashion for the long-range connections. Microcebus belongs to a group that contains the smallest existing primates, with a brain under 2 cm in length. Although the microcebus data is only for V1, it remarkably fits to the same adimensional template, with a decay rate λ between that of mouse and macaque, suggesting that the quantitative differences that distinguish the species are due to both brain size, and primate-rodent differences. The EDR could be the expression of the consequence of a universal information processing principle implemented in the cortex across several scales, specifically to include single neurons in the local circuit, which present over 80% of the cortical connectivity [31,39]. Hence, the two major ingredients for the EDR are found in the local circuitry, the log normal distribution of synaptic weights [59] and an exponential decay of connection distances as reported here. Further, the experimental evidence presented in Fig 11, shows that p (d) follows a nearly identical, exponential decay out to within 1. 5 mm for both mouse and macaque, with λexplocal≅4. 54±0. 08 mm−1. These are gray matter, non-myelinated connections, and are observed to have a very different decay rate than white matter connections. Thus, at least in area V1, the decay of connectivity with distance seems to behave in a very similar fashion in both mouse and macaque, and therefore surprisingly the decay rate in the gray matter does not appear to be related to brain size. Using the reported data in [59] for the rat visual cortex obtained from quadruple whole cell recordings, the local decay rate in the rat can be determined to be λexplocal≅4. 96 mm−1, a value consistent with the one found above in mouse and macaque, above. Table 1 summarizes the EDR related parameters in the mouse and macaque, for both white matter and gray matter connections. The universal character of the EDR is further supported by mathematical arguments. The exponential distribution (EDR) is memoryless (Markov property), i. e. , in our case, the probability that an axon of some length grows by an additional amount is independent of its current length (within cutoff limits). In this way, it has the property that f (d+l) = f (d) f (l), where f (l) = ∫l∞p (l′) dl′ = e−λl is the probability of an axon growing to a length beyond l. The exponential distribution is the only continuous distribution with this property [60]; for all other distributions, growth depends on the current length, i. e. on past growth history. This also implies that the EDR is the maximum entropy probability distribution for axonal lengths with given expectation value (= 1/λ), see [61]. These properties are evolutionarily advantageous, conferring maximum adaptability during cortical expansion. Moreover, as more neurons are added, the probability distribution of the shortest connection among an arbitrary number of other connections also obeys an exponential distribution [62], making the EDR an invariant property locally as well, supported by the experimental data quoted above. The present findings could have important consequences for understanding the human brain. The recognized limitations of current tractographic analysis of diffusion MRI data [63,64], means that direct observation of long-distance connections in the human brain is not presently feasible. Given the specificity of long-range cortico-cortical connectivity [47], this technical limitation has important consequences for understanding the human connectome, and we believe that comparative connectomics as developed in the present study will be a necessary step for determining universal principles of cortical connectivity. Hence, an in-depth understanding of the influence of changes in brain size will play an important role in better understanding the human brain. Since the EDR leads to a decrease in the strength of long-range connections in macaque compared to mouse, we may hypothesize that increase in brain size leads to increased reductions of weight in long-range projections for the whole mammalian branch. In the human brain the small number of fibers in such long distance connections will pose an acute problem for detection for some time. This could constitute an important limitation. For example, one could speculate that the low weight of human long-range connections may contribute to an increased susceptibility to disconnection syndromes, such as have been proposed for Alzheimer disease and schizophrenia [65–67].
Experiments were performed in male and female PV-Cre [68] (Jax: 008069), x Ai9 reporter mice (Jax: 007905), harboring the loxP-flanked STOP cassette, which prevented the transcription of the tdTomato protein driven by the chicken β-actin (CAG) promoter [69]. The crossing produced Cre-mediated recombination, which resulted in the expression of the red fluorescent protein in the subset of parvalbumin (PV) -positive GABAergic neurons. All experimental procedures were approved by the institutional Animal Care and Use Committee at Washington University and conformed to the National Institutes of Health guidelines. Injections were made in Microcebus murinus in area 10 and area V1. Surgical and experimental procedures were in accordance with European requirements 2010/63/UE and approved by the ethics committee CELYNE (ref 00439. 02). For tracer injections, mice were anesthetized with of a mixture of Ketamine (86 mg · kg−1) and Xylazine (13 mg · kg−1, i. p) and secured in a head holder. The body temperature was maintained at 37°C. Intracortical connections within the left hemisphere were retrogradely labeled by inserting a glass pipette (20 μm tip diameter) into the brain and injecting Diamidino Yellow (50 nl, 2% in H2O; EMS-Chemie, Gross-Umstadt, Germany) by pressure (Picospritzer, Parker-Hannafin). Injections were performed stereotaxically 0. 35 mm below the pial surface, using a coordinate system whose origin was the intersection between the midline and a perpendicular line drawn from the anterior border of the transverse sinus at the posterior pole of the occipital cortex. The injections were made in the following areas: V1, RL, AL, LM, P, RSD, ACAd, MOs, SSp-bfd, SSs, Au. Area AM was injected twice (S3 Fig). The parcellation (names and locations of the areas) is based on Wang et al. [70] (S2 Fig) but differs from those used in Fig 4D and associated analyses. Four days after the tracer injection, mice were deeply anesthetized with an overdose of Ketamine/Xylazine and perfused through the heart with phosphate buffered saline, followed by 1% paraformaldehyde (PFA) in 0. 1 M phosphate buffer (PB, pH 7. 4). Immediately after, the cortex was dissected from the rest of the brain, completely unfolded, flat-mounted and post fixed overnight in 4% PFA at 4°C. Next, the tissue was cryoprotected in 30% sucrose and cut at 40 μm on a freezing microtome in the tangential plane. To survey the injection site and the distribution of labeled neurons across cortical areas, sections were wet-mounted in PB and imaged in St. Louis under a dissection scope equipped for UV- and red-fluorescence illumination. For plotting DY labeled neurons, the sections were permanently mounted onto glass slides and stored at 4°C. The distribution of DY-labeled neurons was analyzed in Bron. Plots of DY neurons were made at 20× under a fluorescence microscope equipped for UV illumination (excitation: 387–398 nm, emission: 435–475 nm), using the Mercator software package running on ExploraNova technology. Labeled neurons were contained in 12–16 sections per hemisphere. Digital charts of the coordinates of DY labeled neurons across each section were stored in the computer. Next, the regional pattern in the density of PVtdT expression was imaged under fluorescence optics. Finally, the sections were stained for Nissl substance, imaged under bright field illumination and superimposed onto the digital maps of DY and PVtdT fluorescence. In Bron, all the images were acquired using MorphoStrider software (ExploraNova). The digital charts were saved in PDF files and were scaled in Adobe Illustrator. The charts and the corresponding images were brought to a common scale, allowing reconstruction of the sections. Sections were stacked in order, and then aligned. The landmark for the alignment of the sections was the injection site, followed by rotation around this point, allowing a 3-D reconstruction of the flattened brain. The injected area was delimited, as were the borders of the neocortex. Automated processing was performed using in-house software, written in Python. For each case, the positions of labeled neurons inside neocortex, but outside the injected area (i. e. , extrinsic neurons) were extracted in digital format, The fraction of labeled neurons per area (FLN) was estimated as the number of labelled neurons extrinsic to the injected area expressed as a fraction of the total number of labeled neurons in the cortical hemisphere [39]. Unlike the template matching procedure used in previous studies [5,6] we parcellated each cortex individually based on multiple markers expressed across different tangential sections. In a stepwise procedure we first used density differences in the expression of PV-tdT- labeled cell bodies and processes to delineate borders of single areas such as V1, S1, S2, Au, PD, UF, PV, GU, ORBI, MM. RSD, MOp, MOs, and ENTm (S2 Fig). Next, the PVtdT-expression pattern was used to outline regions which from previous studies are known to include multiple areas. The list includes (1) LM, LI, P, POR, 36p; (2) AL, LLA, RL, A, AM, PM; (3) TEp, TEa, ECT, PERI; (4) AIp, AId, AIv; and (5) ACAd, ACAv, PL, ILA, ORBm, FRP (S2 Fig). Each of these regions was further partitioned into areas based on the topographical distribution of DY labeled neurons, the size and location relative to readily identifiable areas, the rhinal sulcus, the crest of the medial wall [70–74] and the cytoarchitecture revealed by Nissl staining [5,75,76]. The segmentation shown in Fig 4D was carried out in Adobe Illustrator, combining Allen Brain Atlas boundary criteria (visualized with Brain Explorer 2) with photos of PVtdT and Nissl staining, for one section of the flattened mouse brain. The contours of the cortical areas are non-self-intersecting closed polygons; therefore, computing their centroids is straightforward. The distances between areas were considered as the distances between their respective centroids. The Allen Brain Institute (ABI) dataset was collected on their website, offered as link in their original research publication [5]. The University of Southern California (USC) matrix was, on the other hand, extracted directly from their original article [6]. The ABI mouse atlas possesses 40 isocortical areas according to their Supplemental Table 1. Out of these 40 areas, 2 did not correspond to any line or column in the data as structured in S1 Fig (i. e. the connectivity matrix). An additional 4 areas were not considered as primary target of an injection by the authors, leading to our decision to exclude them from our analysis. We then extracted from the ABI data a 34 × 34 weighted and directed connectivity matrix. The larger USC dataset has a finer grained parcellation than that of ABI, although based on the same fundamental scheme. We contracted the USC final matrix down to a level of 42 × 42 by merging areas together in both rows and columns so as to obtain a squared matrix. At this point, the ABI and USC matrices had 33 areas in common, which corresponds to 97% and 79% of their full respective matrices. A similar parcellation scheme was extracted out of the two datasets, allowing complete, connection-by-connection comparison between the two matrices, see S1 Fig for the final connectivity matrix. The database in the mouse has been generated following tracer injections in all cortical areas. The macaque data, however, was obtained from 29 injections using a 91-area atlas. Because in macaque we are using an edge-complete subgraph, the statistical features are predicted to reflect those of the, as yet unavailable, fully connected graph. However, the presently available dataset cannot give complete information on detailed areal relationships, such as for example the full membership of the cortical core. Fig 7C shows the histograms of connection distances for mouse and macaque after normalization by the mean distance for each species. By construction, both distributions have mean equal to 1 and can be reasonably well described by truncated normal distributions. When fitting the distributions by maximum likelihood using functions from the truncnorm package [77] in R [78], the variance of the macaque normalized distances appears smaller than for the mouse data with a ratio of 0. 608. Is this significant? First, we examine this question with an F-test on the ratio of variances. The test is two-sided because we do not specify a priori which variance is greater. This is a more conservative approach. The F-statistic is the variance ratio with degrees of freedom (405,527) giving a highly significant p = 1. 58 × 10−7. The test assumes normality, however. To verify the conclusion, then, we performed a permutation test that does not make the normality assumption [79]. In short, we permute the macaque and mouse labels a large number of times and recomputed the variance ratio for each new permutation. Under the null hypothesis that both distance distributions are the same, we expect a large number of variance ratio estimates on the permuted datasets that are more extreme than the variance ratio computed on the data. The p-value is computed from the proportion of ratio estimates more extreme or equal to the obtained value. For ratios, the definition of more extreme is based on the values that are less than the estimate and greater than its reciprocal. We include the ratio estimate from the dataset in the distribution of permutation estimates. S7 Fig shows the value of the ratio of variances for 100,000 permutations of the two datasets. The vertical line indicates the value obtained from the data, which is lower than all of the other values of the permutation distribution, indicating that the obtained ratio is highly unlikely under the hypothesis that both distributions are the same. The p-value is indicated in the graph. The p-value is smaller than 10−5, which is the resolution of the test for 100,000 permutations. Thus, the width of the distribution of distances for macaque is significantly narrower than that for mouse. To analyze the density of connectivity with distance, we estimate the probability of a connection with distance. This can be done with a logistic regression. By performing the analysis on the binary connectivity (that is, presence/absence of a connection) at each distance, no binning is involved. Standard logistic regression implemented via a Generalized Linear Model with a binomial family [80] specifies that the expected value of the connection probability is related to a linear predictor through a link function that is often taken to be the log of the odds ratio or logit function. The model fit would be g (E (Y = 1) ) = β0, i+β1, iDistance where Y is a binary variable indicating whether a connection is present between two areas, g is the link function, here log (p/ (1 − p) ) with p the expected value or probability of a connection, and β0, i and β1, i are intercept and slope, respectively, of the linear predictor, with i varying with the species. There is no a priori reason to suppose, however, that the sigmoid function of distance that this model implies will provide an adequate description of the change in probability with distance. To allow for a more flexible description of this relation, we fit the data with a Generalized Additive Model (GAM) using a binomial family [81]. The GAM model is given by g (E (Y = 1) ) = Imousefmouse (Distance) +Imacaquefmacaque (Distance) where fi are smooth functions of the covariates constructed from sums of spline curves with increasing complexity and Ii are indicator variables taking on the value of 1 for i = mouse (or, respectively, macaque) and 0 otherwise. The complexity (or wiggliness) of the fitted model is controlled by including a penalty in the fitting criterion based on the integrated square of the second derivatives of the f’s. The choice of degree of penalization (or smoothness) is controlled by minimizing a criterion related to prediction error (i. e. , fitting some of the data and calculating the error on the remaining portion) called the un-biased risk estimator (UBRE) that is closely related to Aikake’s Information Criterion (AIC). Like AIC, UBRE favors a model that maximizes the predictability of future rather than the actual data and serves to minimize the tendency to overfit the data. The fits were performed with functions from the mgcv package [81] in R [78]. The estimates of the smooth curves for macaque and mouse are plotted in Fig 9C and 9D for macaque and mouse, respectively, with twice the estimated standard errors of the fits. To estimate the significance of the differences between the two estimates, we also fit the simpler nested model in which a single smooth curve described connectivity dependence with distance for both species. A likelihood ratio test of the nested models gave a χ2 (2. 26) = 54. 04 with p = 3. 05 × 10−12, strongly supporting that the differences in the curves are significant. Note that the generalization of degrees of freedom in the case of GAM fits are not necessarily integer valued. We extended this analysis to consider the connection probability as a smooth function of both distance and similarity. The GAM framework is used again but now to model surfaces of two variables, here giving the log-odds ratio of the connection probabilities as function of similarity and normalized distance. The model is given as g (E (Y = 1) ) = Imousefmouse (Distance, Similarity) +Imacaquefmacaque (Distance, Similarity), where fi are now smooth 2D functions of the covariates constructed from sums of spline surfaces with increasing complexity and Ii, as before, are indicator variables taking the value of 1 for i = mouse (or, respectively, macaque) and 0 otherwise. Contour plots of the estimates of the connection probability as a function of normalized distance and similarity are shown in Fig 9A and 9B. The color gradient indicates connection probability, passing from high (yellow, near 1) to low (green, near 0) connection probability. The curves indicate estimates of contours of constant connection probability (notated on the curves as probability values). To evaluate the significance of the species difference displayed in Fig 9A and 9B, we also fit the simpler nested model in which a single smooth surface described connectivity dependence with distance and similarity for both species. A likelihood ratio test of the nested models gave a χ2 (8. 4) = 61. 6 with p = 3. 5 × 10−10, strongly supporting that the difference in the surfaces are significant. Note that as above the generalization of degrees of freedom in the case of GAM fits are not necessarily integer valued. The method used to compute binary similarity indices with macaque data has been described previously [4]. Our published macaque database is made of 29 injected areas for a 91 parcellation scheme, thus giving a 91 × 29 connectivity matrix. In this context, only the in-degrees of injected areas are completely known, the out-degrees of source areas remaining incomplete. Therefore, if one wants to compare macaque and mouse using a degree-based binary similarity measure, one has to restrict oneself to in-degrees, in order to use complete data. For this reason, we detail here only the in-degree based similarity measurement calculations. The union of ABI and USC databases used here provides information about all 33 areas in terms of the in and out-going connections between 33 areas. We compared the similarity of the input pattern of pairs of areas by evaluating the number of sources areas from which both receive projections or neither do (i. e. , similarity implies both projections exist or are absent; dissimilarity implies one is absent and the other is present). We define a normalized in-link similarity measure, Sxyin as follows: For any pair of areas (x, y), let nxyin denote the number of projecting areas from which either both x and y or neither x nor y receive an incoming link. Because nxyin≤33, we compute the ratio nxyin/33 for every area pair (x, y). Clearly, this number will depend on the in-degrees of x and y, denoted by kxin and kyin (0≤kx (y) in≤33). We define the in-link similarity as: Sxyin = nxyin33−pxyin, where pxyin is the expected value of the ratio (nxyin/33) if the incoming connections of x and y were distributed uniformly at random across the 33 source areas. Thus: pxyin = (kxin33) (kyin33) + (1−kxin33) (1−kyin33), where the first term is the probability that both x and y receive a link from a given source, and the second term is the probability that neither of them receive a link from a given source.
It was recently shown that the network of connections between different areas of the macaque cortex has strong structural specificity in terms of the strength of connections as a function of the distance between areas. This has led to a model of cortex connectivity that predicts many observed architectural features, including the existence of a strong core-periphery organization. When viewed across species, increases in brain size are accompanied by a relative decrease in connectivity, and thus an important question is whether there are architectural commonalities in the cortical networks within the mammalian branch. Here, based on tract tracing data from the folded macaque brain and the smooth mouse brain, we introduce a common model framework that allows network comparisons between species. We show that despite important differences in size, the cortices of both species share several network invariants, suggesting that the mammalian cortex exhibits universal architectural principals. This framework also captures differences between the two brains, including the fact that, unlike the macaque, the mouse core includes primary areas and that there is a relative decrease in the frequency of long-distance connections in the large macaque cortex compared to mouse. This approach allows network architectural extrapolations to the human cortex.
infographics
medicine and health sciences
neural networks
nervous system
vertebrates
mice
neuroscience
animals
mammals
primates
network analysis
old world monkeys
computer and information sciences
monkeys
network motifs
animal cells
macaque
cellular neuroscience
rodents
cell biology
data visualization
anatomy
central nervous system
neurons
graphs
biology and life sciences
cellular types
amniotes
organisms
Spatial Embedding and Wiring Cost Constrain the Functional Layout of the Cortical Network of Rodents and Primates