Throughout all experiments, animals were treated in accord with the policies and procedures set forth in The National Institutes of Health Guide for the Care and Use of Laboratory Animals and approved regulations of Cornell University's Institutional Animal Care and Use Committee (IACUC). The experiments described in this paper were conducted under IACUC protocol number 84-55-00.

Fifty-four Syrian hamster pups (Mesocricetus auratus) of both sexes from timed pregnancies in the laboratory colony were used in this study. Animals were fed ad libitum and maintained on a 10L:14D photoperiod.

According to convention, the 24-hour period following birth is designated postnatal day 0 (P0). Only hamsters born within 24 hours of the expected 15.5 day gestation period were used for this study. Biocytin was injected into pup cortex at ages P0, P2, P4, P6 and P8, with transport time optimized at 24 hours. Intracortical transport was principally anterograde with very few cortical cell bodies retrogradely labeled outside the immediate injection area. However, both anterograde and retrograde transport were observed to the thalamus, although at ages earlier than P4, transport to the thalamus was principally retrograde. This thalamic label was used to identify thalamic nuclei with connections to the cortical injection site.

Pups were anesthetized by hypothermia and maintained on an ice blanket in molded head and body restraints. The skull was exposed and a hole made overlying the cortical region of interest. A solution of 5% biocytin was injected through a backfilled micropipette (inner diameter 15–20 µm) using a Picospritzer (General Valve Co.; Fairfield, NJ), with pressure and duration adjusted to deliver >0.1–0.5 microliters of solution. Injections were positioned only in cortical regions that could be clearly viewed, avoiding areas of high vascularization. Because rodent intracortical connectivity originates from both infra- and supragranular layers, injections were centered at a depth adjusted for the different ages to span the full thickness of the cortex while avoiding the underlying white matter. Following injections, the scalp was sutured; pups were rewarmed and returned to the mother. After 24 hours pups were overdosed with sodium pentobarbital and perfused transcardially with 0.9% saline followed by 4% paraformaldehyde and 0.1% gluteraldehyde in 0.1 M phosphate buffer (PB, pH 7.4). Brains were cryoprotected in 30% sucrose at 4° centigrade until processing.

All brains were frozen and sectioned coronally at approximately 60 µm. Sections were treated according to a protocol adapted from Ding and Elberger [26], followed by conventional diaminobenzidine (DAB) processing. Briefly, sections were rinsed in phosphate buffered saline (PBS; pH 7.2,) quenched in 1% H₂O₂, and immersed in 1% Triton-X100 (TX) in PBS. Tissue was incubated overnight in an avidin-biotin solution (1∶100; Vectastain Elite Standard Kit) containing 1% TX. Sections were reacted with 0.004% tetramethylbenzidine (TMB); mounted on chromium-gelatin coated slides, dehydrated, cleared, and coverslipped with Krystalon (Fisher Scientific) with one series counterstained with cresyl violet. Dehydration of the unstained series was kept to less than one-minute immersion in each of three graded alcohols to minimize shrinkage. Following data collection, this series was also lightly stained with cresyl violet to further verify neural divisions and landmarks.

Reconstructions were made using a LeitzDiaplan Microscope and a Neurolucida imaging system with a mechanical stage (Microbrightfield, Inc., Colchester, VT). Measurements were obtained from each traced serial section in each of the fifty-four brains, always including sections containing landmarks comparatively stable across development such as the furthest ventral and caudal levels of the white matter, thalamic complex, and caudate nucleus. To avoid artificially elongating in the medial to lateral plane when converting from coronal sections to dorsal views, in each traced section a midpoint contour was measured using a line drawn intermediate between the superficial white matter and the top of cortical layer I. The measurement began medially at the “point of flexure” (dorsalmedial crest separating the two cerebral hemispheres) and extended laterally to the rhinal fissure. A dorsal cortical surface view was then constructed by plotting each midline measurement to scale using Canvas 6.0 (Deneba Systems, Inc.). In effect, this method generates a flattened or “unrolled” surface from curved serial coronal sections as if viewed from above (the dorsal surface; see Figure 4D).

Multiple sources of information were integrated to position general areal boundaries in the developing cortex and locate injection sites. First, atlases of the adult hamster brain [27] and developing and adult rat brain [28]–[30] were used to establish the overall orientation and form of the map. It is well established that the topology of thalamocortical and corticothalamic projections is conserved from first innervation to adulthood, though maturational gradients in the cortex and the deformation of the cortex by overall growth alters the relative size of cortical regions [31]–[34]. The anterograde and/or retrograde transport of biocytin from injection sites to the thalamus, listed in Table 1, was used to fix the positions of the primary visual, auditory and somatosensory cortex and shift areal boundaries with respect to the adult cortex as required, for each postnatal age. Published and unpublished hamster developmental studies from this laboratory were also used to help align regions of the developing cortex with respect to other adjoining telencephalic regions, such as the striatum and hippocampus, whose topological positions remain fixed from initial generation to adulthood [33], [35], [36].

Thirty-six developing brains with well-labeled axons were completely analyzed microscopically and 24 representative brains were traced using a Neurolucida (25×) for more detailed morphological and statistical analysis, including ages (injected-recovered) P0–1 (n = 7), P2–3 (n = 4), P4–5 (n = 3), P6–7 (n = 4), P8–9 (n = 3). Axons were identified by their coloring, thin uniform appearance, characteristic branching patterns, and, on many occasions, the presence of growth cones. Every visible intracortical axon in each traced section was drawn. Sections to be traced (typically over one half) were determined by the presence or absence of labeled axons, although as noted above, sections containing the furthest ventral and caudal levels of the white matter, thalamic complex, and caudate were always traced to obtain registration measurements for dorsal views.

The tracings of coronal sections were then used to generate dorsal view reconstructions of the furthest distal points where labeled axons were found, as well as axon density plots of projections arising from injection sites. First, radial lines were drawn perpendicular to the middle layers of the gray matter and spanning the entire depth of the white and gray matter, spaced every 200 µm beginning at the point of flexure and ending at the rhinal fissure, with the last measurement the interval between the final 200 µm line and the rhinal fissure. Tangential substrates were then outlined using pseudo phase-contrast on unstained tissue and adjusted using counterstained sections (accounting for shrinkage, which was consistently less than 10%). The tangential substrate boundaries included the subjacent border of the cortex (layer VI), the subjacent border of the infracortical fasciculus (a cell-sparse area above the subplate neurons, also called “channel 2” in [29], and the subjacent border of the subplate neurons, which corresponds to the superficial border of the white matter. Axons crossing each radial line were counted and tangential substrate subtotals obtained. Each of these counts was associated with a location (i,j) on the 60 µm×200 µm grid, indicating that this is the count at the ith radial line in the jth section. The substrate subtotals were labeled , for those axons counted in cortical layers, and , for those in the white matter. We use to denote the collapsed count of all axons encountered at a radial line, i.e. for each count location .

The tangential compartments are not uniformly identifiable in hamster cortex: in far anterior and posterior coronal sections, white matter fibers, subplate neurons and the infracortical fasciculus merge. Moving anterior to posterior in pup brains, white matter fibers are first noted at the level of orbital cortex where the rhinal fissure no longer clearly separates cortex and olfactory bulbs, followed approximately 0.5 mm posterior by a distinct layer of subplate neurons and approximately 0.5 mm further posterior by the band of fibers comprising the infracortical fasciculus. In far lateral regions of cortex, the subplate neurons seem to merge with neurons of the claustrum; in both far lateral and posterior cortex, the fasciculus is quite thick relative to anterior sections (see also [37]). The position of isolated single axons or very sparse projections was always registered independent of the radial lines used for systematic sampling. Care was taken to ensure that counts did not include axons traveling to subcortical or callosal areas; however, although these axons travel different routes, some small uncertainty is unavoidable. Dorsal view reconstructions were produced for projections traveling from an injection site in each of the different substrates (cortical layers I–VI, infracortical fasciculus, subplate and white matter), as well as a “collapsed” view of the first three combined so as to represent the conventionally identified cortical gray matter.

The surface area covered by underlying axons within the borders of the cortex bounded by the point of flexure and rhinal fissure was determined for 20 traced pup brains using NIH Image. Total area of axon coverage (in mm²) was analyzed for cortex, subplate, infracortical fasciculus, and for a category collapsed across these three, as well as for the white matter. These totals were expressed as a percentage of the total dorsal area or anterior/posterior (A/P) or medial/lateral (M/L) length in each individual brain. Schematized cortical areas were not used in statistical analysis; areas were determined for each individual brain at each age (see also Table 2). ANOVA was performed in Statview 5.0 to determine if axon extent in M/L or A/P planes varied based on injection site location or age, followed by Scheffe's Post-hoc test when appropriate.

The data for each animal was recorded as a set of axon counts taken at points on a 2-D grid whose axes aligned with the medial/lateral (ML) and anterior/posterior (AP) axes of the flattened cortical hemisphere. Indexing each grid point (i,j) and calling the corresponding count , these counts document the number of axons originating at the injection site which were detected at location (i,j) on the 200 µm by 60 µm grid of sample points. As described above, for each animal, three sets of count data were analyzed: axons found in the cortical layers only , axons in the white matter only , and a collapsed set , where .

In order to arrive at the desired description in terms of probability distribution functions of axon terminal sites, the following steps were carried out. Any missing counts from the interior of each dataset were interpolated. Each grid was re-centered such that the injection site (detected as the site having the maximum axon count ) had index . We assumed that the injection site is a point source of axons. To simplify further calculations, a 2-dimensional first order interpolating function was fitted to each grid (using Mathematica). With the interpolating function , it was possible to treat the count data as continuous over the 2-D domain with .

Two functions, calculated using , were used to characterize each dataset. The first, , accounts for the angular distribution of axons (see Figure 1A). The data exhibited a prevalence of growth along the direction of the ML axis in preference to the AP axis. To quantify this anisotropy, we fitted a double peaked function ,defined on . Here with parameter is a probability distribution on the circle (see Figure 1B). The distribution is flat at and approaches two delta spikes as . For each animal, we calculated the values of and , or the anisotropy and tilt as we call them respectively, which minimized the least squares error between and (see Figure S2 for fits to data). Thus is a probability distribution for the relative volume of outgrowth in each direction, parametrized in each case by anisotropy and tilt .

Arriving at the radial distribution function, characterizing the length distribution of the axons, requires taking into account the cumulative nature of the count data: is the number of axons one would expect to find passing through any point , not the number of axons terminating there. Given that axons may trace more or less circuitous routes between their origin and terminal arbor, and may also branch en route, it will not be possible to exactly recover the density distribution of endpoints from the count data. To arrive at an approximation to the true distribution, we assume (i) that no branching occurs prior to arrival at the terminal site and (ii) that axons travel from the origin along straight trajectories. These assumptions are generally consistent with the data collected in this study and in earlier work [33]. Given that our network model is constructed using such straight-line axon trajectories, disregarding the (unknown) particulars of axons routes will not affect our simulated networks.

The quantity , under assumptions (i) and (ii) above, can be interpreted as the probability that an axon has length greater than or equal to r (see Figure 1A). The probability that it terminated at a length less than r is simply and so, in principle, is an approximation to the cumulative distribution function (CDF) for the length of the axons. In practice however, portions of the empirical “CDFs” fail to meet the monotonicity property required of distribution functions. The difficulty arises near the origin, at small radii r, where the resolution of the experiment means that few count sites are contributing to the calculated value of . Furthermore, what counts are present may be noisy due to the high concentration of stained axons close to the injection site. For this reason, we chose to disregard the non-monotonic portions of the empirical functions near the origin and fit analytical CDF's to the remaining monotonic data. However, it is that portion of the data, near the origin, which would otherwise provide the normalizing constant for our distributions (i.e. ). Hence, the fitting procedure must also provide an estimate for that normalizing constant. On inspection, the usable portion of the empirical distributions seemed to be well fit by gamma distributions, which have previously be used to model axon length distributions [38]. So as to reduce the number of fitting parameters, we constrained the fitting procedure to use gamma distributions with a shape parameter equal to 2. Thus the radial distribution of each outgrowth pattern was characterized by the mean of its fitted gamma distribution, . See Figure S1 for fits to data.

We study a network model whose nodes are localized populations of neurons, linked by edges which model representative axons. The network is constructed as follows (see also Figure 2). We model the hamster's cortical hemisphere as comprising 2,500 idealized cortical units. Each unit represents the neural population under a 100 µm×100 µm square of cortex. Our model cortex is a grid of 50×50 = 2,500 of such units, thus mimicking a cortical sheet measuring 5 mm×5 mm. These units are the nodes in our network. With each node i we associate a location , chosen uniformly at random within the corresponding unit's square footprint.

Emanating from each node are a fixed number, , of directed network edges, representing efferent axons. Edges are assigned to each node i by using the following procedure: (i) select a length r and a direction from under the empirically derived radial and angular probability distributions, and , respectively; (ii) find the node j nearest to the point at distance r from node i in the direction ; (iii) add a directed edge, pointing from node i to node j, to our network. There are two exceptions to the above procedure. First, if the randomly chosen length r and direction determine a point which lies outside the grid, we chose a new random r and , repeating until the chosen r and determine a point inside the grid. Second, self edges are disallowed – should one occur, we discard it and choose a new random endpoint such that the axon does not terminate at its node of origin. The procedure for generating edges is repeated until the node i has been assigned its full complement of edges. Edges are assigned to all nodes in the network in this manner. Once all edges have been assigned, in the adjacency matrix A, the entry A_ij records the number of directed edges from node i to node j.

We note that although the details of meandering axonal paths were ignored as we deduced and , our method of generating networks does not require such details to reproduce the distribution of axon terminal sites. Our networks were created and visualized using Mathematica 7.0 [39].

We calculated several measures to analyze the networks generated by our spatial model. The out-degree of a node on a directed network counts the number of edges originating at that node. By construction, we have for all nodes in our networks. The in-degree counts all edges terminating at a node: . The number of edges comprising shortest network path between any pair of nodes, i and j, also known as the distance from i to j, is recoded as . The average shortest path, , gives a characteristic length for the network paths. The clustering of a network gives the probability that if the edges i→j and j→k are both present then so too is the edge k→i. The clustering, C, is given by , where, for example, i→j→k is a path of length 2 and a “triangle” refers to a case where the three edges (i→j, j→k and k→i) are present [40]. A small-world network is characterized by having both the short path lengths typical of a random network and the high clustering typical of a more regularly wired (e.g. lattice) network. Using random networks as a baseline, the small-world index, S, makes the classification of networks as being small-world quantitative [41]. The small-world index of a network with N nodes, M edges, clustering C and average shortest path length l, is given by , where and are, respectively, the clustering and average shortest path length of a random network also having N nodes and M edges. For a random network, , but S takes larger values for small-world networks.

Global network efficiency quantifies the efficiency of communication between all pairs of nodes on the network, under the assumption that information flows along the shortest paths available [42]. Considering just one pair of nodes first, if an edge joins the two nodes, the path between them has length 1 and so communication is maximally efficient for that pair: we say that path has an efficiency of 1. If the shortest path between a pair of nodes (i and j) has length , then we say its efficiency is . The average value of that pairwise efficiency, taken over all pairs of nodes in the network, is the global efficiency: . Only a completely connected network (where all possible edges are present) has a global efficiency of 1. All other networks have an efficiency .

The betweenness centrality of a node is a measure of how important that node is for efficient communication on the network [40]. Considering the set of all shortest paths on the network, we see that some “central” nodes may feature in a greater number of shortest paths than do other less central nodes. The betweenness centrality of a node i is the fraction of all shortest paths on the network which pass through i. Specifically, if is the number of shortest paths from s to t, and is the number of such paths containing the node i, the betweenness centrality of i is given by .

It can be useful to think of the nodes on a network as being members of different communities. To investigate the modular nature of our networks, we will assign nodes to non-overlapping communities whose membership is defined by location. If the communities are chosen well, one should observe a greater prevalence of intra-community edges over inter-community edges than would be found in a comparable random network (i.e. a randomly wired network with the same number of nodes and edges). The modularity Q quantifies the extent to which we have such a prevalence [43]. Assigning each node to a community, and denoting the node i's community c_i, we can measure the modularity of the network with respect to that community assignment:

where M is the total number of edges in the network and is equal to 1 if and is zero otherwise. We measure the modularity of our networks with respect to two different, spatially defined, community assignments. First, we partition the nodes into 4 rectangular communities, roughly equal in size, aligned with the medial-lateral axis. Namely, the communities contain the nodes in rows 1 through 13, 14 through 25, 26 through 37 and 38 through 50 of our 50×50 grid, respectively). The communities are thus rectangles whose width spans the medial-lateral dimension of our gird and whose height is about one fourth that of the grid. Second, we use a similar partitioning, but in this case assigning nodes by column number to one of four rectangular communities aligned with the anterior-posterior axis.

The NetworkX package [44] for Python was used to carry out the graph analysis of our networks.

Exploiting the spatially embedded nature of our network, we investigate how anisotropy may affect the volume requirement of axons via an altered number of axon encounters (see Figure 3). Modeling all the axons as straight lines in the plane, we count the number of axon encounters. Only encounters along the body of an axon are considered, those occurring at a node are disregarded. Axon encounters were enumerated using an algorithm implemented in c++.

The greatest numbers of intracortically confined projections are local, extending in a radial fashion for short distances in the gray matter directly adjacent to the injection site. Longer-range projections traveling away from injection sites take multiple paths, coursing through the conventionally identified gray matter, the infracortical fasciculus, and among the subplate neurons, as well as in the white matter itself (see Figure 4B). In these small brains, axons reaching the most distant point from the injection site might equally be found traveling through the cortex or through subcortical white matter — the route taken by an axon, within grey matter or white matter did not dictate the distance traversed. Although this study traces the distribution of a population of labeled axons, the trajectories of individual axons were noted when they could be followed. Some of the individual axons that could be traced from injection site to growth cone or terminal arbor traveled exclusively in one substrate, others switched pathways, for example, from the infracortical fasciculus to the white matter, avoiding the subplate (Figure 4 C); or alternated travel between two substrates.

The mean of the available neural area in each substrate covered from the injections (expressed as a percentage of total dorsal cortical area) is as follows: cortex mean: 61.21%, SE: 3.24, white matter mean: 38.09%, SE: 3.16; subplate mean: 32.34%, SE: 2.59; infracortical fasciculus mean: 30.53%, SE: 2.97 (Figures 5,6). Thus, these small tracer injections resulted in extensive axon spread in the cortex, and insofar as it was possible to rank injection sites by size, the only difference associated with injection size was projection density, not extent.

Even at the earlier age, connections from each site span almost the entire medial/lateral (M/L) distance of the cortex; this coverage persists as total cortical area more than doubles between ages P0–1 and P8–9 (Table 2). This widespread axon extension, evident across all locations, is illustrated in the dorsal reconstructions depicted in Figure 5. Because extension within the subplate and infracortical fasciculus was always less than extension within the cortex itself, in Figure 5 these three divisions are collapsed into one “gray matter” compartment. Qualitative inspection and statistical analyses indicate no significant differences in axon coverage from any site (anterior, middle or posterior) at any postnatal age, when total area of coverage was expressed as a percentage of total cortical area, regardless of the substrate of axon extension.

Tracer placements that happened to bridge more than one cortical area, as judged by retrograde labeling of both primary visual and somatosensory thalamic nuclei, might be expected to produce larger ranges of axon travel if each cortical area specifies a unique list of termination addresses, as contrasted with a model of initial axon outgrowth independent of cortical area identity. Though the number of cases we could use to address this question is small, examination of the area of cortex labeled by tracer injections in the several cases that resulted in retrograde label to both somatosensory and visual thalamic nuclei (45.0%, n = 3, across ages) compared to injections that labeled either one or the other class (70.2%, n = 7, across ages), however, showed the opposite, though non-significant trend.

Despite the general similarity of widespread coverage patterns from P1–P9, some local patterns were evident. As illustrated by the gradient outlines in Figure 5, long-range connections from anterior cortex extend asymmetrically posterior and laterally towards middle cortex, not colonizing far posterior cortex. Axons labeled by injections in middle cortical areas appear reciprocally focused on anterior cortex with few posterior projections. The majority of developing intracortical projections from posterior cortex is confined to the posterior cortex itself, although from all sites a number of axons typically extend to presumed medial limbic regions.

Axon travel in the conventionally recognized pathway, the white matter, is summarized in Figure 6. Travel in the white matter is more confined than travel in the overlying substrates, and more anisotropic, as shown by a comparison of the axon coverage in Figure 6 to the outlines of axon coverage in the collapsed views illustrated in Figure 5.

Because intracortical axons travel in large numbers through the cortex as well as the conventionally identified white matter, these schematized dorsal representations of axon populations do not distinguish axons traveling intracortically from those which exit the cortex, travel in white matter, and re-enter cortex to terminate. The generally uniform picture of axon extension that the dorsal view reconstructions suggest is perhaps at odds with the presence of a large number of abrupt trajectory shifts in axon tracks which might suggest alteration axon extension by detection of an areal boundary (e.g., Figure 4C, left). Trajectory shifts in white-matter axons might indicate regions of cell-to-substrate recognition that specify unique termination zones. For the two distributions which most strongly suggested the development of a particular termination focus, one early and one late in development, we charted every axon that turned from a horizontal trajectory to travel vertically (Figure 7). The locations of abrupt trajectory shifts populated the entire area of extension, however, producing a distribution that was simply a reduced form of the overall distribution. This pattern is more suggestive of random sampling of the substrate by the axon population prior to target selection [45], [46]. Overall these relatively uniform distributions resemble those seen in studies which reconstructed single-axon termination patterns within single areas in larger brains [25].

Given the empirical observations above, we characterize the typical outgrowth pattern as having the following features: (i) having areal coverage larger than half of the cortical hemisphere; (ii) comprising axons travelling in both the white and gray matter which traverse comparable distances; (iii) having a footprint with greater extent along the ML axis than along the AP axis, with travel in the white matter being comparatively more constrained in this regard; (iv) being largely independent of the location of its source. With these features in mind, we developed the following framework to arrive at a quantitative description.

Probability densities for the angular and radial components, denoted and , respectively, were fitted to the observed outgrowth distributions. The resulting values of the characteristic length , anisotropy and tilt for each dataset are shown in Figure 8, and can be summarized as follows.

The radial distribution shows a marked departure from uniformity (see Figure 8B). The measured values of anisotropy for collapsed distribution have mean 0.69 and standard deviation 0.07 on our 0 to 1 scale (see Figure S2 for fits to data). There is a clear tendency for white matter distributions to be more anisotropic () than distributions measured only in the cortex (). The tilt is such that the preferred direction of travel is almost collinear with the ML axis, as illustrated in Figure 8C ().

The length distribution of axons in the white matter, gray matter and collapsed distributions were well fit by gamma distributions (see Figure S1 for fits to data), with mean lengths as shown in Figure 8A. There is a tendency for distributions measured in older animals to be longer and for white matter axons to travel greater distances. Taking the ratio for each animal, we find that is on average 119±33% as long as - indicating a trend of axons travelling in white matter reaching further than those traversing cortex only. This greater length is not incompatible with the lesser areal coverage of white matter axon distributions; travel in this compartment was also noted to be more anisotropic and therefore has a more slender footprint.

We sought to create a model of the early cortical network which was faithful both to the qualitative characteristics and the measured distributions of axon outgrowth. We also wanted the nodes of this network to be biologically meaningful but representing neurons as individual nodes would have made our model computationally unwieldy. We instead take our nodes to be cortical populations or units, comprising all the neuronal cell bodies and local processes within a small volume of the cortical sheet, much like the “cortical output units” of Innocenti and Vercelli [47]. We define these units as 100 µm×100 µm squares on a 2-dimensioal cortical sheet, so that they have roughly the dimensions of a dendritic arbor. Hence our “local” length-scale is set at ∼100 µm. The number of nodes is thus conducive to constructing and analyzing our simulated networks on a desktop computer.

The simulated networks are found to have the small world property but are not scale free. The small world index is measured to be 4.57±0.17, indicating our networks possess short path lengths comparable with those of random networks while also having far higher clustering. This may have been anticipated given the form of the axonal distributions employed, which have a preponderance of local connections and relatively few long links. The networks are not in the class of so-called scale-free networks, with neither their in-degree nor out-degree having the required power-law distributions. The out-degree is, by construction, the same for every node and is equal to the bespoke number of efferent axons per node, . The in-degrees are in close agreement with those for a random network, only differing slightly from the Poisson distribution one would encounter in that case. We can see why such an in-degree distribution arises: if we ignore boundary effects, then nodes attract afferent axons with a probability proportional to the area of their respective Voronoi cells, the polygons in the plane enclosing all points to which a given node is the closest node. Given our method of distributing nodes in the plane, the preponderance of Voronoi cells not adjacent to borders will have area comparable with that of our 100 µm×100 µm grid squares. Such nodes therefore attract afferent edges with approximately equal probability, thus giving rise to the observed narrow, near-Poisson distribution of in-degrees.

Anisotropy in axonal distributions may lead to a more volume-efficient wiring scheme but would seem, prima facie, to entail negative repercussions for the resulting network structure. Given that increased anisotropy in laying down axons leads to a smaller ensemble of possible networks (to see this consider the limiting case where axons are restricted to travel along only one axis), we were curious as to what advantages it might bestow. While investigating the effect of varying anisotropy α on the topology of simulated networks, we also sought to determine whether increased α might lead to a more compact packing of axons in space. We begin with the observation that an axon's length, and thus its volume, is increased if it must deviate to avoid another axon (see Figure 3). Therefore, encounters with other axons tend to increase the total volume requirement.

Our spatial network model predicts that there exists a narrow range of values for the anisotropy parameter α which reduces the extra volume requirement of crossing axons without significantly impacting the ease of communication within the network, as measured by the efficiency E (see Figure 9A). Increasing the anisotropy beyond this range is predicted to reduce network efficiency. We find that the empirically measured values of anisotropy fall within this range for both white matter and gray matter travel.

Increased anisotropy may lead to networks which are more vulnerable in the case of central nodes failing. The presence of such nodes is indicated by a more skewed distribution of node betweenness centrality (see Figure 9B; also Figure S3 for histograms of betweenness centrality distributions). The betweenness centrality b_i of the node i quantifies the relative contribution of individual network elements to the collection of shortest paths on the network (see Figure 10). Conversely, it can be seen as a measure of how detrimental the removal of that element might be to the functioning of the network as a whole. We observe that increased anisotropy in the spatial distribution of links leads to an increased positive skew in the distribution of node betweenness centrality in simulated networks (Figure 9B). The right tail of the distribution becomes increasingly “heavy”, signifying the emergence of a small number of nodes with betweenness centrality much higher than the average value. As studies in other complex networks have shown [48], the existence of such highly central nodes renders the network as a whole more vulnerable to functional disruption in the event of their failure. It is also the case, however, that communication on such a network is more robust to the failure of a typical node, i.e. one which is not among the few highly central nodes.

The possibility of nodes becoming overloaded may be another reason to disfavor highly anisotropic wiring schemes. Assuming that information is propagated along the shortest paths between nodes, betweenness centrality can be interpreted as a measure of how much traffic a node handles. Given the finite neural populations comprising our nodes, they will have a limited capacity to process or propagate information. Hence, with highly anisotropic wiring, nodes which are very central may be at risk of becoming overloaded, thereby affecting the reliability of communication on the network.

Looking at the modularity of our networks as anisotropy is increased, we observe that partitioning the network into communities aligned with the medial-lateral axis of our grid is increasingly favored over choosing communities aligned perpendicular to that axis (see Figure 11). This may favor a layout in which more strongly connected and functionally related cortical areas have that same axial alignment. We investigated modularity with respect to partitioning the nodes into 4 spatially defined communities: first, communities which were (approximately) equally-sized rectangular blocks of nodes, aligned with the medial-lateral axis; second, a similar arrangement but with the rectangular blocks now parallel to the anterior-posterior axis. In both cases, the preponderance of short links, which leads to spatial clustering, ensured that our networks were more modular than the comparable network wired uniformly at random. However, the modularity was seen to increase significantly with anisotropy for the medial-lateral partitioning while decreasing for the anterior-posterior partitioning.