User:Jack D. Cowan/Proposed/Models of visual hallucinations
Contents |
Form constants and visual imagery
Seeing vivid visual hallucinations is an experience described in almost all human cultures. Painted hallucinatory images are found in prehistoric caves (Clottes and Lewis-Williams, 1998) and scratched on petroglyphs (Patterson, 1992). Hallucinatory images are seen both when falling asleep (Dybowski, 1939), and on waking up (Mavromatis, 1987), following sensory deprivation (Zubek, 1969), after taking ketamine and related anaesthetics (Collier, 1972), after seeing bright flickering light (Purkinje, 1823; Helmholtz, 1925; Smythies, 1960), or on applying deep binocular pressure on the eyeballs (Tyler, 1978), in “near death” experiences (Blackmore, 1992), and most strikingly, shortly after taking hallucinogens containing ingredients such as LSD, cannabis, mescaline, or psilocybin (Siegel and Jarvik, 1975). The images are stable with respect to eye movements and are seen both by blind subjects and in sealed dark rooms. This suggests that they are generated in the brain. One possible location for their origin is provided by fMRI studies of visual imagery which suggest that primary visual cortex (V1), which is the first cortical area to receive visual information from the retina, is activated when human subjects are instructed to inspect the fine details of an imagined visual object (Miyashita, 1995). In 1928 Klüver, a University of Chicago neurologist, organized such images into four groups called form constants (Klüver 1966): (I) tunnels and funnels, (II) spirals, (III) lattices, including honeycombs and triangles, and (IV) cobwebs, all of which contain repeated geometric structures. Figure 1 shows examples of each of the form constants as they appear in the visual field.
The retinocortical map
Assuming that the form constants are generated in V1, the first step in constructing a theory of their origin is to compute their appearance in V1 coordinates. This is easily done since the visual world is mapped onto the cortical surface in a topographic manner, which means that neighboring points in a visual image evoke activity in neighboring regions of visual cortex. Moreover, one finds that the central region of the visual field has a larger representation in V1 than the periphery (Drasdo 1977, Sereno et al 1995), partly due to a non-uniform distribution of retinal ganglion cells. The retinotopic map is defined as the coordinate transformation from points in the visual world to locations on the cortical surface. In order to describe this map, it is first necessary to specify visual and cortical coordinate systems. Since objects located a fixed distance from one eye lie on a sphere, it is natural to introduce spherical coordinates with the "north pole" of the sphere located at the fixation point, the image point that focuses onto the fovea or center of the retina. In this system of coordinates, the latitude angle is called the eccentricity \(\epsilon\) and the longitudinal angle measured from the horizontal meridian is called the azimuth \(\theta\ .\) In most experiments the image is on a flat screen such that, if the curvature of the sphere is ignored, the pair \((\epsilon,\theta)\) approximately coincides with polar coordinates on the screen. In primary visual cortex the visual world is split in half with the region \(-90^o \leq \theta \leq 90^o\) represented on the left side of the brain, and the reflection of this region represented on the right side brain. Note that the eccentricity \(\epsilon\) is usually specified in terms of the angle subtended at a fixed distance from the screen. The structure of the retinotopic map in monkey is shown in figure 2. This suggests that eccentricity is mapped on to the horizontal coordinate x of the cortical sheet, and \(\theta\) is mapped on to its y coordinate An approximate equation for the retinotopic map can then be obtained through specification of a quantity known as the cortical magnification factor. This determines the distance across a flattened sheet of cortex separating the activity evoked by two nearby image points. It can be shown that sufficiently far from the fovea, the coordinate transformation is well approximated by a scaled version of the complex logarithm (Schwartz 1977, Cowan 1977): \[ x=\alpha \log [\epsilon/\epsilon_0],\quad y = \beta \theta \] where \(\alpha,\beta, \epsilon_0\) are constants. Under this complex logarithm, the funnel and spiral form constants shown in Figure 1 are mapped to stripes of neural activity in cortex. This simple observation has led to a mathematical theory of geometric visual hallucinations based on spontaneous pattern formation in large-scale cortical models (Ermentrout and Cowan 1979, Bressloff et al 2001).
Neural Turing mechanism
A possible mechanism for the spontaneous formation of stripes of neural activity under the action of hallucinogens was originally proposed by Ermentrout and Cowan (1979). They studied interacting populations of excitatory and inhibitory neurons distributed within a two-dimensional cortical sheet. Modeling the evolution of the system in terms of a set of Wilson-Cowan equations (Wilson and Cowan 1972, 1973), they showed how spatially periodic activity patterns such as stripes can bifurcate from a homogeneous low-activity state via a Turing-like instability analogous to that found in reaction diffusion systems. The model also supports the formation of other periodic patterns such as hexagons and squares; under the retinocortical map these generate more complex hallucinations in the visual field such as checkerboards, see figures 3 and 4. Similar results can be obtained in a reduced one-population model provided that the synaptic interactions are characterized by a mixture of short-range excitation and long-range inhibition (the so-called “Mexican hat distribution”).
The neural field equation for a one-population model takes the form \[ \frac{\partial u({\mathbf r},t)}{\partial t} = - u({\mathbf r},t)+\int_{\R^2}w(|{\mathbf r}-{\mathbf r}'|)f(u({\mathbf r}',t))d{\mathbf r}' \] Here \(u(\mathbf r,t)\) represents the activity of a population of neurons in a given slab of cortical tissue located at position \({\mathbf r}\) at time \(t\ .\) The second term on the right-hand side represents the total synaptic input into the population at \({\mathbf r}\) with \(w\) the synaptic weight distribution and \(f\) a nonlinear firing rate function. The latter is taken to be the sigmoidal function \[ f(u)=\frac{f_0}{1+\exp(-\eta(u-\kappa))} \] where \(\eta\) determines the gain and \( \kappa \) the threshold. The weight distribution \(w\) is assumed to depend on the Euclidean distance \(|{\mathbf r}-{\mathbf r}'| \) between the presynaptic and postsynaptic populations in the cortical plane. It follows that \(w\) is invariant with respect to the action of the Euclidean group \(E(2)\) – rigid body translations, rotations and reflections in the plane. This symmetry plays a crucial role in determining the types of spontaneous patterns that can be generated by a Turing-like instability.
Suppose that for a given set of parameters, the neural field equation supports a stable, spatially uniform low activity state given by the fixed point solution \(u({\mathbf r},t) = u_0\) for all \({\mathbf r},t\ .\) One possible effect of the action of drugs on certain brain stem nuclei is to increase the excitability of the network through changes in the effective gain or threshold of the neurons, which could then destabilize the fixed point. The local stability of the fixed point can be determined by linearization. That is, set \(u({\mathbf r},t)= u_0+\mbox{e}^{\lambda t}a({\mathbf r})\) and Taylor expand the neural field equation to first order in \(u\ .\) This leads to the eigenvalue equation \[ \lambda a({\mathbf r}) = - a({\mathbf r})+\mu \int_{\R^2} w(|{\mathbf r}-{\mathbf r}'|)a({\mathbf r}')d{\mathbf r}', \] where \(\mu=f'(u_0)\ .\) Since the local weight distribution \(w\) is homogeneous, it follows from continuous translation symmetry that the eigenmodes are in the form of plane waves \(a({\mathbf r})=\mbox{e}^{i {\mathbf r}\cdot{\mathbf k}}\) with wavevector \({\mathbf k} \ .\) The corresponding eigenvalues satisfy the dispersion relation \[ \lambda = \lambda(q) \equiv -1 +\mu {W}(q), \] where \(q=|{\mathbf k}|\) and \({W}\) is the two-dimensional Fourier transform of \(w\ ,\) \[ {W}(q)=\int_{\R^2} w(|{\mathbf r}|)\mbox{e}^{-i {\mathbf k}\cdot{\mathbf r}}d {\mathbf r}. \] It is then relatively straightforward to set up the conditions under which the homogeneous state undergoes a Turing instability as the excitability parameter \(\mu\) increases, namely, that \(W(q)\) be bandpass. This can be achieved with the “Mexican Hat” function shown in Figure 5a, representing short-range excitation and long-range inhibition. It is simple to establish that λ then passes through zero at the critical value \(\mu_c = 1/W(q_c)\) signaling the growth of spatially periodic patterns with wavenumber \(q_c\ ,\) where \(W(q_c) = \max_q\{W(q)\}\ .\) For sufficiently small \(\mu\ ,\) the dispersion curve \(\lambda(q) \) is negative for all q and the uniform low-activity state is stable. As \(\mu\) increases, the dispersion curve crosses zero at the non-zero wavenumber \(q_c\ ,\) which is illustrated graphically in Figure 5b.
Amplitude equations
Close to the bifurcation point the emerging patterns can be represented as linear combinations of plane waves \[ u({\mathbf r}) =\sum_{n=1}^N\left [z_n\mbox{e}^{i {\mathbf k}_n \cdot {\mathbf r}}+ z_n^*\mbox{e}^{-i {\mathbf k}_n \cdot {\mathbf r}}\right ] \] where the sum is over all wave vectors on the critical circle \(|{\mathbf k}|=q_c\) and \((z_n,z_n^*)\) is a complex conjugate pair of amplitudes. Rotation symmetry implies that the space of such modes is infinite-dimensional. That is, all plane-waves with wavevectors on the critical circle are allowed (see Figure 6a). However, translation symmetry means that we can restrict the space of solutions to that of doubly-periodic functions corresponding to regular tilings of the plane. The symmetry group is then reduced to that of certain crystal lattices: square, rhomboid and hexagonal lattices (see Figure 6b). The sum over n is now finite with N = 2 (square, rhomboid) or N = 3 (hexagonal) and, depending on the boundary conditions, various patterns of stripes or spots can be obtained as solutions. Amplitude equations for the coefficients \(z_n\) can then be obtained using perturbation methods and weakly nonlinear theory, although their basic structure can be deduced from the underlying symmetry (Hoyle 2006).
In the case of a square or rhombic lattice, we can take \({\mathbf k}_1 = q_c(1,0)\) and \({\mathbf k}_2 = q_c(\cos \varphi, \sin\varphi)\) such that to cubic order \[ \frac{dz_n}{dt} = z_n\left [\mu - \mu_c -\gamma_0 |z_n|^2 -2\gamma_{\varphi} \sum_{m\neq n}|z_m|^2\right ] \] for n = 1,2, where \(\gamma_{\varphi}\) depends on the angle \(\varphi\ .\) In the case of a hexagonal lattice we can take \({\mathbf k}_n = q_c(\cos \varphi_n \sin\varphi_n)\) with \(\varphi_1= 0, \varphi_2= 2\pi/3, \varphi_1 = 4\pi/3 \) such that \[ \frac{dz_n}{dt} =z_n\left [\mu - \mu_c -\gamma_0 |z_n|^2 -\zeta z_{n-1}^*z_{n+1}^* \right ] -2\gamma_{\varphi_2} z_n \left (|z_{n-1}|^2+|z_{n+1}^2|\right ) \] where n =1,2,3 (mod 3). The fixed points of these ordinary differential equations can then be analyzed to determine which particular types of pattern are selected and to calculate their stability. The results can be summarized in a bifurcation diagram as illustrated in Figure 6c for the hexagonal lattice with \(\zeta > 0,\ 2\gamma_{\varphi_2} > \gamma_0\ .\) Note that in certain cases it may be necessary to go to higher orders to identify all of the basic patterns guaranteed by the underlying symmetry according to the Equivariant branching lemma (Golubitsky et al 1988, Hoyle 2006).
Analogous patterns have also been observed in fluids in the form of convection rolls and honeycombs as well as in animal coat markings in the form of stripes and spots. This indicates that although the physics may be very different, the interactions in all these phenomena are such that they can all be represented within the framework of the Turing mechanism (Turing 1952) and symmetric bifurcation theory.
Feature maps and the functional geometry of V1
The Ermentrout-Cowan theory of visual hallucinations is over-simplified in the sense that V1 is represented as if it were effectively a cortical version of the retina, in which cells only signal the level of local light illumination. However, most V1 cells are selective to a variety of features of a local stimulus such as the local orientation of a contrast edge or bar (Hubel and Wiesel 1974). Stimulus feature preferences such as orientation are distributed in an approximately periodic fashion across cortex forming a set of overlapping feature maps. Moreover, these feature maps are correlated with the distribution of long-range synaptic connections within cortex thus providing a substrate for the functional geometry of V1. The absence of orientation representation in the Ermentrout-Cowan model means that a number of the form constants cannot be generated by the model including form constants III and IV shown in Figure 1. These hallucinations are more accurately characterized as lattices of locally oriented contours or edges rather than in terms of contrasting regions of light and dark. This has motivated the development of a more general theory of visual hallucinations that takes into account the functional geometry of cortex (Bressloff et al 2001, 2002a,b).
Considerable information regarding the functional geometry of V1 has been revealed using the method of optical imaging (Blasdel and Salama 1986, Bonhoeffer and Grinvald 1991). The basic experimental procedure involves shining light directly on to the surface of the cortex. The degree of light absorption within each patch of cortex depends on the local level of activity. Thus, when an oriented image is presented across a large part of the visual field, the regions of cortex that are particularly sensitive to that stimulus will be differentiated. The topography revealed by these methods suggests that there is an underlying periodicity in the microstructure of V1 with a period of approximately 1mm (in cats and primates). The fundamental domain of this periodic tiling of the cortical plane is the hypercolumn (Hubel and Wiesel 1974), which contains two sets of orientation preferences per eye, organized around a pair of orientation singularities or pinwheels (Obermayer and Blasdel 1993). Combining optical imaging with a variety of stains and labels has made it possible to see how neurons within and between hypercolumns are interconnected, as is illustrated in Figure 7 for the tree shrew. These data suggest that there are at least two length-scales for cortico-cortical interactions:
local - cells less than a millimeter apart tend to make connections with most of their neighbors in a roughly isotropic fashion
long-range - cells make patchy connections every millimeter or so along their axons with cells in similar iso-orientation domains.
Figure 7 also shows that the long-range patchy connections, known as intrinsic lateral or horizontal connections and found mainly in superficial layers of V1 (Rockland and Lund 1983, Gilbert and Wiesel 1983), have a pronounced anisotropy. That is, the major axis of the distribution of patchy connections tends to run parallel to the visuotopic axis of their cell's orientation preference (Bosking et al 1997). Hence, the long-range connections between iso-orientation domains tend to form continuous contours following the topography of the retinocortical map. This anisotropy is the anatomical substrate of the Euclidean shift-twist symmetry underlying the hallucination model of Bressloff et al (2001).
Cortical model with shift-twist symmetry
Each location in the visual field is represented in V1, roughly speaking, by a hypercolumn-sized region containing all orientation preferences \(\theta\) (running from \(0^o\) to \(180^o\)). This suggests introducing an extended coordinate system in cortex given by \(({\mathbf r},\theta) \) such that all possible stimulus orientations are represented at each point in visual space under the inverse retinotopic map. Cortex is effectively treated as the continuum limit of a lattice of interconnected hypercolumns. The natural extension of the Ermentrout-Cowan model is then (Bressloff et al 2001) \[ \frac{\partial u({\mathbf r},\theta,t)}{\partial t} = - u({\mathbf r},\theta,t)+\int_{\R^2}\int_{0}^{\pi}w({\mathbf r},\theta|{\mathbf r}',\theta')f(u({\mathbf r}',\theta't))d\theta'd{\mathbf r}' \] The structure of the weight distribution \(w({\mathbf r},\theta|{\mathbf r}',\theta')\) connecting neurons with cortical coordinates \(({\mathbf r}',\theta')\) to neurons with coordinates \(({\mathbf r},\theta)\) is chosen so as to model the functional architecture of V1 encapsulated by Figure 7 and is illustrated schematically in Figure 8: \[ w({\mathbf r},\theta|{\mathbf r}',\theta')=w(\theta-\theta')\delta({\mathbf r}-{\mathbf r}')+J(|{\mathbf r}-{\mathbf r}'|)W(\theta-\theta')\Delta(\arg({\mathbf r}-{\mathbf r}' - \theta)) \] Here w represents the local isotropic connections within a hypercolumn and is taken to be a Mexican hat function as in the ring model (Ben-Yishai et al 1995, Somers et al 1995). The function J determines the distant-dependent variation of patchy connections between hypercolumns and the narrowly tuned function W ensures only neurons with similar orientations are linked via patchy connections. Finally the function \(\Delta\) ensures that only iso-orientation domains that approximately lie along the visuotopic axis connecting cells in different hypercolumns are linked, reflecting the anisotropy shown in Figure 7.
The total weight distribution is invariant under the so-called shift-twist action of the Euclidean group acting on the coordinates \(({\mathbf r},\theta,t)\ .\) That is, it is invariant with respect to translations \(({\mathbf r},\theta)\rightarrow ({\mathbf r}+{\mathbf s},\theta)\ ,\) reflections \((x,y,\theta)\rightarrow (x,-y,-\theta)\) and rotations \(({\mathbf r},\theta)\rightarrow ({\mathbf R}_{\phi}[{\mathbf r}],\theta+\phi)\ ,\) where \({\mathbf R}_{\phi}[{\mathbf r}]\) is the vector \({\mathbf r}\) rotated by the angle \(\phi\ .\) This form of the rotation operation follows from the anisotropy of the lateral weighting function and comprises a translation or shift of the orientation preference label \(\theta \) to \(\theta +\phi\ ,\) together with a rotation or twist of the position vector \({\mathbf r}\) by the angle \(\phi\ .\) Such a shift–twist operation (Zweck and Williams 2000) provides a novel way to generate the Euclidean group E(2) of rigid motions in the plane. The fact that the weighting function is invariant with respect to this form of E(2) has important consequences for the selection and stability of patterns that emerge via the Turing instability. In particular, close to the bifurcation point, the emerging patterns now belong to one of two classes, even (+) or odd (-), and can be represented as linear combinations of \(\theta\) modulated plane waves according to \[ u({\mathbf r}) =\sum_{n=1}^N\left [z_nu^{\pm}(\theta-\theta_n)\mbox{e}^{i {\mathbf k}_n \cdot {\mathbf r}}+ z_n^*u^{\pm}(\theta-\theta_n)\mbox{e}^{-i {\mathbf k}_n \cdot {\mathbf r}}\right ] \] Here \(\theta_n\) is the direction of the nth wavevector on the critical circle, \({\mathbf k}_n=q_c(\cos \theta_n,\sin\theta_n)\ ,\) and \(u^{\pm}(\theta)\) are even and odd functions of \(\theta\ ,\) respectively.
One can now distinguish between three basic classes of stripe pattern: non-contoured (\(u=1\)), even contoured (\(u=\cos \theta\)) and odd contoured (\(u=\sin \theta\)), see Figure 8. Linear combinations of the non-contoured stripes recover the patterns of the Ermentrout-Cowan model, whereas linear combinations of the contoured stripes generate new patterns that reproduce the edge-like hallucinatory form constants when mapped back to the visual field, as illustrated in Figure 9.
<review>A minor edit to associate this article with the correct category</review>
