Self-organization and the selection of pinwheel density in visual cortical development

Self-organization of neural circuitry is an appealing framework for understanding cortical development, yet its applicability remains unconfirmed. Models for the self-organization of neural circuits have been proposed, but experimentally testable predictions of these models have been less clear. The visual cortex contains a large number of topological point defects, called pinwheels, which are detectable in experiments and therefore in principle well suited for testing predictions of self-organization empirically. Here, we analytically calculate the density of pinwheels predicted by a pattern formation model of visual cortical development. An important factor controlling the density of pinwheels in this model appears to be the presence of non-local long-range interactions, a property which distinguishes cortical circuits from many nonliving systems in which self-organization has been studied. We show that in the limit where the range of these interactions is infinite, the average pinwheel density converges to $\pi$. Moreover, an average pinwheel density close to this value is robustly selected even for intermediate interaction ranges, a regime arguably covering interaction-ranges in a wide range of different species. In conclusion, our paper provides the first direct theoretical demonstration and analysis of pinwheel density selection in models of cortical self-organization and suggests to quantitatively probe this type of prediction in future high-precision experiments.


Introduction
Neuronal circuits in the mammalian cerebral cortex are among the most complex systems in nature. The biological mechanisms that contribute to their formation in early brain development remain poorly understood. However, it is unlikely that the precise architecture of mature cortical circuits can be attributed to genetic prespecification, since the number of genes in the genome is insufficient [1]. Instead, dynamical self-organization presumably plays a major role in shaping the architecture of neuronal circuits in the cerebral cortex. Dynamical self-organization is most thoroughly described in non-living physical systems driven outside of thermodynamic equilibrium by external forcing. Whereas the emergence of structure is externally driven, the structures formed are primarily determined through interactions within the system itself [2,3]. Neural circuits encompass various positive and negative feedback loops, and these could well form the basis for cortical self-organization. However, evidence supporting this presumption is derived from theoretical considerations [4,5,6] rather than empirical observation. Models for the self-organization of neural circuits have been proposed, but experimentally testable predictions of these models have been lacking.
The system of orientation columns in the visual cortex is a paradigmatic system for studying cortical development and the role of self-organization in this process. Most neurons in the visual cortex respond selectively to a particular orientation of an elongated visual stimulus. Whereas in columns perpendicular to the cortical surface, neurons prefer similar stimulus orientations, the preferred orientation varies mostly smoothly [7] and repetitively across the cortical surface giving rise to a complex two-dimensional pattern called the map of orientation preference (Fig. 1a, b). Throughout the cortical map, there are point-like orientation singularities [8,9,10] called pinwheel centers [11] at which all stimulus orientations are represented in circular fashion. Numerous studies are consistent with the hypothesis that orientation maps develop through activity-dependent self-organization. They form in dark reared animals [12], under substantial manipulation of visual input [13,14], and even in auditory cortex when rewired to be driven by visual inputs [15]. Moreover, an analogy between cortical development and pattern formation appears plausible. Like in other systems where pattern formation has been observed [2,3], the orientation map arises probably from an initially non-selective state, it exhibits a typical periodicity and a spatial extension at least an order of magnitude larger than the basic periodicity length.
The conditions under which orientation maps can arise through self-organization have been thoroughly investigated theoretically [16,17,18]. A recent and highly promising approach stressing the analogy to pattern forming systems showed that pinwheels can be stabilized by activitydependent long-range interactions [19]. A phenomenological order parameter field model based on the Swift-Hohenberg equation [20] was proposed in which orientation maps arise from a supercritical bifurcation of Turing-type. In this class of models, the stabilizing nonlinearity includes only key features of visual cortical organization and is constraint by biologically plausible symmetry assumptions. The model exhibits multiple structurally distinct quasiperiodic attractors resembling orientation maps in the visual cortex.
The qualitative similarity of solutions of this model to orientation maps in the visual cortex appears promising. However, it is unclear at present weather the model accounts for aspects of cortical organization also quantitatively and whether this resemblance is not just superficial. Comparing directly orientation maps in the model and the visual cortex would be difficult, because of the large number of possible map structures. Instead, it will be virtually unavoidable to take a statistical approach to this questions. As discrete entities, pinwheels can be detected, counted and localized with high accuracy. Here, we calculate in the model the average density of pinwheels in the limit of infinite interaction range. We show that this density is representative for a large regime of intermediate interaction ranges covering the estimated ranges in various mammalian species. Therefore, the quantity pinwheel density appears particularly well suited for testing for signatures of long-range dominated self-organization in experiment. In the following, we briefly describe the system of orientation preferences in the visual cortex and present the long-range self-organization model for its activity-dependent development.

Orientation preference maps
In the visual cortex, as in most areas of the cerebral cortex, information is processed in a 2-dimensional (2D) array of functional modules, called cortical columns [21,22]. Individual columns are groups of neurons extending vertically throughout the entire cortical thickness that share many functional properties. Orientation columns in the visual cortex are composed of neurons preferen- Figure 1: Patterns of orientation columns and long-range horizontal connections in the primary visual cortex of tree shrew visualized using optical imaging of intrinsic signals (modified from [23]). a, Activity patterns resulting from stimulation with vertically and obliquely oriented gratings, respectively. White bars depict the orientation of the visual stimulus. Activated columns are labeled dark gray. The used stimuli activate only columns in the primary visual area V1. The patterns thus end at the boundary between areas V1 and V2. b, The pattern of orientation preferences calculated from such activity patterns. The orientation preferences of the columns are color coded as indicated by the bars. A part of the pattern of orientation preferences is shown at higher magnification. Two pinwheel centers of opposite topological charge are marked by arrows. c, Long-range horizontal connections extend over several millimeters parallel to the cortical surface (tree shrew, superimposed on the orientation preference map). White symbols indicate locations of cells that were filled by a tracer (biocytin); labeled axons are indicated by black symbols. tially responding to visual contours of a particular stimulus orientation [7]. In a plane parallel to the cortical surface, neuronal selectivities vary systematically, so that columns of similar functional properties form highly organized 2D patterns, known as functional cortical maps (Fig. 1a, b).
Experimentally, the pattern of orientation preferences can be visualized using optical imaging methods [8,9]. Optical imaging of intrinsic signals is based on the fact that the optical properties differ in active vs. less active parts of the cortex [24]. This is utilized to record patterns of activity from light reflectance. In a typical experiment, the activity patterns E k (x) produced by stimulation with a grating of orientation θ k are recorded. Here x = (x, y) represents the location of a column in the cortex. Using the activity patterns E k (x), a field of complex numbers z(x) can be constructed that completely describes the pattern of orientation columns: The pattern of orientation preferences ϑ(x) is then obtained from z(x) as follows: Typical examples of such activity patterns E k (x) and the patterns of orientation preferences derived from them are shown in Fig. 1a, b. Numerous studies confirmed that the orientation preference of columns is an almost everywhere continuous function of their position in the cortex. The domains formed by neighboring columns with similar orientation preference are called isoorientation domains [25].

Pinwheels
At many locations the iso-orientation domains are arranged radially around a common center [9,10]. Around these pinwheel [11] centers, stimulus orientations are represented in circular fashion. Such an arrangement had been previously hypothesized on the basis of electrophysiological experiments [26,27] and theoretical considerations [28]. The regions exhibiting this kind of radial arrangement were termed pinwheels [11]. The centers of pinwheels are point discontinuities of the field ϑ(x) where the mean orientation preference of nearby columns changes abruptly. They can be characterized by a topological charge which indicates in particular whether the orientation preference increases clockwise around the center of the pinwheel or counterclockwise. Here, C j is a closed curve around a single pinwheel center at x i . Since ϑ is a cyclic variable within the interval [0, π) and up to isolated points is a smooth function of x, q i can only have the values where n is an integer number [29]. If its absolute value |q i | is 1/2, each orientation is represented once in the vicinity of a pinwheel center. Pinwheel centers with a topological charge of ±1/2 are simple zeros of z(x). In experiments only pinwheels that had the lowest possible topological charge q i = ±1/2 are observed. This means there are only two types of pinwheels: those whose orientation preference increases clockwise and those whose orientation preference increases counterclockwise. This organization has been confirmed in a large number of species and is therefore believed to be a general feature of visual cortical orientation maps [30,31,32,33,23,34].

Hypercolumn and pinwheel density
The pattern of preferred orientations is roughly repetitive [7,21]. The column spacing Λ, i.e. the spacing between adjacent iso-orientation domains preferring the same stimulus orientation, is typically in the range of ∼ 1mm. The column spacing Λ determines the size of the cortical hypercolumn, which is considered to be the basic processing unit of the visual cortex [7,22,35]. We define the size of a hyper column by Λ 2 . The pinwheel density is defined as the number of pinwheels per unit area Λ 2 . Thus, by this definition, the pinwheel density is independent of the spacing of columns and dimension-less.

Intra-cortical connectivity
Visual cortical neurons are embedded in densely connected networks [36]. Besides a strong connectivity vertical to the cortical sheet between neurons from different layers within a column, neu-rons also form extensive connections horizontal to the cortical surface linking different orientation columns. These connections extend for several millimeters parallel to the cortical surface and are therefore called long-range horizontal connections. As shown in Fig. 1b for the example of the tree shrew, these connections are clustered primarily connecting domains of similar orientation preference. They have been observed in various mammals [37,23,38,12] and repeatedly hypothesized to play an important to role in visual processing tasks such as contour integration.

Activity-dependent development
In normal development, orientation columns first form at about the time of eye opening [39,13,12] which for the ferret is approximately at post natal day (PD) 31. This is just a few days after neurons first respond to visual stimuli. A subset of neurons show orientation preference from that time on, but the adult pattern is not attained until seven weeks after birth [40]. Roughly clustered horizontal connections are present by around PD 27 [41]. Many lines of evidence suggest that the formation of orientation columns is a dynamical process dependent on neuronal activity and sensitive to visual experience [17,42]. This is suggested not only by the time line of normal development, but also receives support from various experiments manipulating the sensory input to the cortex. Most intriguingly, when visual inputs are rewired to drive what would normally become primary auditory cortex, orientation selective neurons and a pattern of orientation columns even forms in this brain region that would normally not at all be involved in the processing of visual information [15,42]. This observation suggests that the capability to form a system of orientation columns is intrinsic to the learning dynamics of the cerebral cortex given appropriate inputs. Moreover, the comparison of development under conditions of modified visual experience demonstrates that adequate visual experience is essential for the complete maturation of orientation columns and that impaired visual experience, as with experimentally closed eye-lids [13,12] or by rearing kittens in a striped environment consisting of a single orientation [14], can suppress or impair the formation of orientation columns (but see also [43]). Viewed from a dynamical systems perspective, the activity-dependent remodeling of the cortical network is a process of dynamical pattern formation. In this picture, spontaneous symmetry breaking in the developmental dynamics of the cortical network underlies the emergence of cortical selectivities such as orientation preference [16]. The subsequent convergence of the cortical circuitry towards a mature pattern of selectivities can be viewed as the development towards an attractor of the developmental dynamics [18]. This is consistent with the interpretation of cortical development as an optimization process. In the following, we will briefly describe a model [19] that is based on this view.

Modeling cortical self-organization
Self-organization has been observed to robustly produce large scale structures in various complex systems. Often, the class of patterns emerging depends on fundamental system properties such as symmetries rather than on system specific details. Pattern formation can therefore often be described by abstract models incorporating these properties only. For systems undergoing a Turing type instability, canonical model equations are of the Swift-Hohenberg [20,2] type where the linear part is and z(x, t) is a complex scalar field. If the bifurcation parameter r < 0, the homogeneous state z(x) = 0 is stable. For r > 0, a pattern with wavelength close to Λ = 2π/k c emerges. The lowest order nonlinearities N 2 and N 3 are quadratic and cubic in z, respectively. The form of these nonlinearities determines the class of the emerging pattern, i.e. whether hexagons, rotating stripes, spiral waves, or another type of pattern emerges.
Following this paradigm, we adopted a model of the form Eq. (5) with nonlinearities derived from key features of the visual cortex (following [19,44]). As in experimental recordings, orientation columns are represented by a complex field [28,18] where ϑ is the orientation preference and |z| a measure of selectivity at location x = (x, y) in the map. The factor 2 in the exponent accounts for the π-periodicity of stimulus parameter orientation. Constructing the nonlinearity of the model relies on the following assumptions. The model includes the effects of long-range intracortical connections between columns with similar orientation preference (Fig. 1c). Unlike many non-living systems in which interactions are purely local, long-range interactions are an important and distinctive feature of neuronal circuits in the cortex and particularly in the primary visual cortex. Further, based on the spatial homogeneity of circuits across cortex, it is assumed that the dynamics is symmetric with respect to translations, and rotations of the cortical sheet. This means that patterns that can be converted into one another by translation or rotation of the cortical layers belong to equivalent solutions of the model, Eq. (5) , by construction. It is further assumed that the dynamics is symmetric with respect to shifts in orientation Thus, two patterns are also equivalent solutions of the model, if their layout of orientation domains is identical, but the preferred orientations differ everywhere by the same constant angle. Solutions shall contain representations of all stimulus orientations. For simplicity, couplings to other visual cortical representations such as ocular dominance or retinotopy are neglected. Considering only leading order terms up to cubic nonlinearities a nonlinearity fulfilling these requirements is given by with long-range interactions mediated through convolutions with a Gaussian with range σ. The second parameter 0 ≤ g ≤ 2 controls local and nonlocal interactions. The first term is the only strictly local term consistent with the required symmetries, the second nonlocal term represents the simplest non-local term that is symmetric with respect to these and with respect to permutations Here the non-linear operator is written in a trilinear form as introduced in [19,44]. This additional symmetry implies that all two-orientation solutions, for instance real valued solutions, of the model, are unstable, which in turn guarantees that all stimulus orientations are represented.
One should note that in this model the spatial range of the nonlinearity σ is a control parameter independent of the wavelength Λ. The patterns selected for different ratios σ/Λ are displayed in Fig. 2b. It can also be derived as an approximation to a model for the combined development of orientation preference and long-range horizontal connections [44]. The model as defined in Eq. (5-13) is variational. It is consistent with synaptic models based on Hebbian plasticity, e.g. [17,45,46]. It is the only model known to the authors that exhibits stable aperiodic solutions dominated by a single spatial scale Λ. These solutions resemble orientation maps observed in the visual cortex. In the following, we discuss the structure of these solutions and the phase diagram of the model.

Weakly nonlinear analysis
Solutions of the model close to the bifurcation point r = 0 are known in closed form, derived by means of a perturbation method called weakly nonlinear stability analysis [47,2]. When the dynamics is close to a finite wavelength instability, the essential Fourier components of the emerging pattern are located on the critical circle. Solutions are planform patterns composed of a finite number of Fourier components with wavevectors on the critical circle, |k j | = k c . By symmetry, the dynamics of amplitudes A i of a planform are governed by amplitude equa-tionsȦ where j − denotes the index of the mode antiparallel to mode j. The form of Eq. (15) is universal for models of a complex field z satisfying symmetry assumptions (8)(9)(10). All model dependencies are included in the coupling coefficients g ij and f ij and may be obtained from F [z] by multiscale expansion [47,2]. Denoting the angle between the wave vectors k i and k j by α and δ ij the Kronecker delta, the coefficients read [19,44]  where are called angle-dependent interaction functions. Stationary solutions of Eq. (15) are given by families of planforms [19,44] z of order n with wavevectors distributed equidistantly on the upper half of the critical circle and binary values l j = ±1 determining whether the mode with wave vector k j or with wavevector −k j is active. These planforms cannot realize a real valued function and are called essentially complex planforms (Fig. 2). For these planforms the third term in Eq. (15) vanishes and the effective amplitude equations for the active modes reduce to a system of Landau equationṡ with stationary solutions Eq. (18) with amplitudes of equal modulus and an arbitrary phase φ i independent of the mode configuration l j . If the dynamics is stabilized by long-range nonlocal interactions (g < 1, σ > Λ), large n planforms are the only stable solutions.
In this long-range regime, the order of planforms grows as approximately linear with the interaction range. For a given order n, different planforms are degenerated in energy. This is a consequence of the permutation symmetry Eq. (13). This symmetry also implies that the relevant stable solutions are essentially complex planforms which in turn guarantees that all stimulus orientations are represented.

Calculation of pinwheel density 2.1 Large range limit of interactions: Planform anisotropy
First, we calculate the average pinwheel density ρ l for an ensemble of planforms, Eq. (18), with a fixed set of wavevector directions l = (l 0 , l 1 , . . . , l n−1 ) but arbitrary phases φ j in the limit n → ∞.
Here, and until noted otherwise, shall denote average over phases φ j . In this limit, z and local linear functionals of z have Gaussian statistics such that the density of pinwheels is determined by the second order statistics of the field. Second, we evaluate the expectation value of ρ l over all sets of l. For large σ, to good approximation g ii ≈ 1 and g ij ≈ g and hence |A i | ≈ 1/ √ ng. Planforms (18) simplify to where for later convenience the constant √ 2g was absorbed into z(x). Pinwheels are the zeros of the field z(x). The number of pinwheels in a given area A is obtained by where δ(x) denotes Dirac's delta function and is the Jacobian of the field split here for later convenience into its real and imaginary part. Averaging Eq. (24) over the ensemble of phases φ j reads implying that is the expectation value of the pinwheel density for a fixed set of l. This expectation value only depends on local quantities, namely on the field, Eq. (23), and its spatial derivatives such that knowledge of the joint probability density p (z, ∇z) is sufficient to evaluate Eq. (28).
Owing to the central limit theorem, this probability density becomes Gaussian in the large n limit. Eq. (28) is then determined by the first and second order statistics of z and ∇z. Furthermore, since their statistics is the same at each location x, it is sufficient to evaluate Eq. (28) for z(0), ∇z(0). The spatial dependency is thus omitted in the following. The average in Eq. (28) is given by an integral over the joint probability density with covariance matrix C which shall be analyzed in the following. First, the diagonal elements of C are evaluated. Using (18) and (26), the auto-correlations of the field are for the real part and for the imaginary part. For the spatial derivatives one obtains where k xj is the x-component of k j and, likewise, The equality of these correlations follows from inserting Eq. (19) into (34) yielding and also such that all auto-correlations become when choosing without loss of generality the column spacing to be Λ = 2π/k c ≡ 1.
Most off-diagonal elements of the covariance matrix C vanish. All non-vanishing contributions are related to the planform anisotropy defined by which depends on the set l. The covariance between the field and its derivative reads with χ ≡ 1 n n−1 j=0 l j k j ≥ 0. The modulus χ = | χ| is small for an isotropic distribution of wavevectors l j k j . To estimate its upper bound χ max , consider the most anisotropic case with all l j k j situated in the right plane (l j = 1 for j ≤ n/2, l j = −1 for j > n/2). For large n this upper bound is such that the modulus ξ = ξ of the anisotropy is bounded by 0 ≤ ξ ≤ 1. In the following, without loss of generality, we assume ξ = ξ(1, 0) implying that all correlations involving one derivative in y-direction vanish. Correlations involving ∂ x are obtained by writing Eq. (40) and in the form and comparing both imaginary parts showing that does not vanish for anisotropic planforms. Expression (44) are the only non-vanishing non-diagonal elements of the matrix C. Indeed, follows from comparing both real parts in Eq. (43). Furthermore, correlations between the real and imaginary part and between their derivatives, e.g.
vanish since they contain terms of the form sin φ j cos φ j ′ = sin φ j cos φ j = 0. Finally, because vanishes for arbitrary n, also correlations between derivatives in different directions, do not contribute to the density of pinwheels. Altogether, the covariance matrix for the vector v = (R, I, ∂ x R, ∂ x I, ∂ y R, ∂ y I) reads The integral we want to solve is For the exponent we find Performing the integral over R and I, we get where w = (∂ x R, ∂ x I, ∂ y R, ∂ y I). Substituting and the Jacobian reads such that we obtain after integrating over angles. Evaluating the determinant as and performing the final integration steps, we obtain the pinwheel density of an ensemble of planforms with a fixed set of wavevector directions l. This result shows that the pinwheel density ρ l only depends on the anisotropy ξ of the considered planform. Since the anisotropy ranges within 0 ≤ ξ ≤ 1, the pinwheel density is confined by 1.36 ρ l ≤ π.

Distribution of planform anisotropies
What is the distribution of anisotropies ξ in the large n limit? To address this question consider the ensemble defined by the different sets l = (l 0 , . . . , l n−1 ). From now on, shall denote the expectation value over this ensemble. In the following we assume that the distribution of the vector anisotropies ξ is isotropic and Gaussian in the large n limit. The isotropy follows from the rotation symmetry of the model equations. It implies that The assumption of Gaussian statistics of ξ is justified by the fact that by its definition (39), the vector anisotropies ξ results from various wavevectors l j k j with pairwise independent directions l j . The distribution thus reads with variance ν ξ given by where l j l j ′ = δ jj ′ . The probability density for ξ follows from the distribution of the vector anisotropy ξ by where in the first equation the prefactor accounts for the change to polar coordinates.

Pinwheel density in the large n limit
Distribution (64) and Eq. (60) yield a distribution of pinwheel density ρ l by means of the coordinate transform ξ → ρ l . With and ξ = (π 2 − ρ 2 l )/8 (Eq. (60)) the pinwheel density distribution reads The expectation value of the pinwheel density follows as where Φ i is the imaginary error function. In the limit n → ∞, the second term vanishes implying Moreover, the pinwheel density is δ-distributed at ρ = π in the large n limit as its variance ν converges towards 0. Having which goes to π 2 for n → ∞, one finds suggesting that for large n not only the average pinwheel densities ρ but the density ρ of almost every realization is close to π. Thus, in the limit of infinite interaction range, the value of the average pinwheel density is π in almost every realization.

Intermediate range of interaction
Next, we ask to which extent these results do remain valid for finite interaction ranges and weather they are fairly robust against variations in range. Real cortices are of course finite and the effective interaction range for a given species or animal may depend on the range of intracortical long-range horizontal connections (see Fig. 1 Figure 3: Pinwheel densities for finite interaction range n ≈ 2πσ/Λ. a, Average pinwheel density ρ in planforms of order n (diamonds). Planforms with random sets of wavevector directions l j and phases φ j were synthesized on 2048 × 2048 grids with aspect ratio Γ = 32. For each planform the pinwheel density was determined in a quadratic Γ = 8 subregion. Planforms were drawn until a achieving SEM ∆ = 0.03 for n < 15 and ∆ = 0.01 for n ≥ 15. The dashed line represents ρ = π, valid in the limit n → ∞. Note that the average pinwheel density ρ is close to π even at intermediate orders of n. b, Standard deviation (SD) s of densities ρ. Note that s decreases with order n.
of size 8Λ × 8Λ. Realizations are collected until sufficient precision of average pinwheel densities ρ is reached, measured by the standard error measure (SEM) ∆ = s/ √ N , were s is the standard deviation (SD) of densities ρ and N the number realizations. The afforded precision of ∆ < 0.03 for 5 ≤ n ≤ 14 and ∆ < 0.01 for 15 ≤ n ≤ 20 requires between 100 and 1000 realizations per order n. At a given order n, we determine the average pinwheel density ρ and the SD s from the ensemble of realizations. Note that by this method, both quantities can be evaluated with arbitrary precision by using sufficient large number of realizations. Fig. 3a shows the average pinwheel density ρ of planforms of various orders ranging between 3 ≤ n ≤ 25. Among smaller orders n ≤ 7, average densities ρ fluctuate substantially covering the range 2.5 ≤ ρ ≤ 3.5. However, they are much more confined, between 2.9 ≤ ρ ≤ 3.2, for intermediate orders 8 ≤ n ≤ 15. For large order n > 15, the ensemble average ρ appears to converge towards π from below. This is consistent with Eq. (60) showing that already for intermediate orders of n averaging over l leads to pinwheel densities smaller but close to π. That average pinwheel densities ρ are smaller than π is explained by the upper bound of ρ l in Eq. (60) implying an upper bound also for its average. That they are, in fact, not much smaller than π is suggested by the "relativistic" form of Eq. (60). Even moderate anisotropies up to, e.g., ξ l = 0.3 result in pinwheel densities ρ l > 3. Furthermore, as shown in Fig. 3b, the variation s of pinwheel densities ρ in different realizations decreases drastically with n. For large n, s became successively smaller consistent with the limiting value s 2 = µ = 0 derived above. These results show that even for intermediate interaction ranges, both the average pinwheel density ρ of an ensemble of maps and the pinwheel density ρ of almost every single realization are close to π.

Discussion
We derived signatures of cortical self-organization that can be tested experimentally. In a model for the self-organization of the system of orientation preferences in the visual cortex we calculated the density of pinwheels (topological defects) that the model is predicting. We find that pinwheel densities close to π are robustly selected if interactions between remote contour detectors are prevalent. Near criticality (r ≪ 1), in the limit of large intracortical interaction range σ, the average pinwheel density converges to the fixed number ρ = π. For intermediate ranges 1Λ σ 4Λ, average pinwheel densities ρ remain smaller but close to this limit value. For successively larger ranges σ, the pinwheel density approaches π from below. Moreover, the variation of pinwheel densities ρ across realizations decreases with interaction range such that when the interaction range is large, almost every map exhibits a pinwheel density close to π. Thus, for a broad parameter regime of interaction ranges, the model predicts an average pinwheel density close to π even in individual cortical orientation maps.
The results presented here are obtained close to the bifurcation point (r ≪ 1). In this regime, the model is well approximated by the amplitude equations and its behavior is representative for a large class of models. Thus, the results derived in this regime are robust against variation of the model. Investigations further away from criticality at r > 0 would be instructive, but can at present only be carried out by numerical methods. While dependencies on the details of the model are expected to play a role for r > 0, solutions in this regime may still be categorized by those obtained for r ≪ 1.
In the visual cortex, distant neurons can directly and strongly interact by long-range horizontal connections and by feedback connections from one visual area to another. In mature cortical circuits, such connections link neurons that share similar receptive field properties and represent somewhat displaced locations in visual space [37,23]. In fact, in a wide range of mammals, the range of these connections is much larger than the spacing of columns [37,23,12]. Thus, these species are candidates for testing the value of the average pinwheel density predicted by our study.
The approach presented in this study, to our knowledge, for the first time allows to analytically characterized the statistics of pinwheels in populations of asymptotically stable states of visual cortical self-organization. The results obtained are valid for the general class of Turing-type systems exhibiting the four symmetries Eq. (8)(9)(10)13) and dominated by long-ranging interactions when studied close to the instability threshold. Asymptotically stable states have been studied in two competing model classes. In one class of models, pinwheels annihilate dynamically, such that pinwheel density is a time dependent quantity and no particular finite pinwheel density is intrinsically selected [18,48,49,50]. In other models, a substantial number of pinwheels are preserved in the final state, but the pinwheels typically crystallize into repetitive arrangements [48,51,49,52]. These states are easy to characterize, but are not even qualitatively consistent with the aperiodic arrangement of visual cortical orientation columns. The model studied here is distinguished from these competing models by exhibiting in the long-range interaction dominated regime (g < 1, σ > Λ) a large multiplicity of aperiodic stable solutions. Our study presents a transparent analytical approach for studying pinwheel density selection in such models. It will be interesting to generalize this approach to further model classes to comprehensively clarify whether genuinely different models of visual cortical development can be distinguished through their pinwheel statistics, as the above results suggest.
It is natural to ask whether and how the predicted pinwheel density can be tested experimentally. The currently available published data is consistent with a pinwheel density around 3 but appears insufficient to obtain a high precision estimate of this quantity (see e.g. [53,54], for a summary of the earlier literature see [18]). To obtain a precise estimate both pattern wavelength and absolute pinwheel number need to be reliably quantified. Over the past years, image analysis methods have been devised that enable estimating the pattern wavelength of orientation maps with a precision in the range of a few percent ( [35,55]). For instance, using a wavelet method for the estimation of local pattern wavelength Kaschube and coworkers showed that genetically related cats often differ in the mean wavelength of orientation columns by less than 4 % ( [35]). As the count variance of pinwheel number estimates is of order (number of pinwheels) 1/2 , estimating the pinwheel density with an accuracy in the percent range will require data sets encompassing on the order of 10000 genuine pinwheels. In many model animals, this size of data set is equivalent to more 100 brain hemispheres that need to be imaged under consistent experimental conditions. Whereas this may appear a large animal cohort for any individual study, such large size data sets will eventually accumulate in many laboratories that use intrinsic signal optical imaging as a standard technique, enabling to test models of visual cortical development with quantitative precision.