Walter J. Freeman Journal Article e–Reprint

Neural mechanisms underlying destabilization of cortex by sensory input

Reprinted From: Physica D 75 (1994) 151-164

Walter J. Freeman

Abstract

A mathematical model of the mechanism of olfaction has been developed in the form of coupled nonlinear integrodifferential equations. The model was first constructed to simulate the performance of the separated bulb or in response to electrical impulse stimulation of its input and output nerves. With sufficiently high internal feedback gains, the model entered a stable limit cycle state at frequencies characteristic of the EEG bursts. The instability was diminished when the connection densities and input amplitude were spatially nonuniform, though these factors could he compensated by increasing the internal gains (connection strengths). This model failed to simulate three important characteristics of olfactory activity: (1) the dependence of spatial patterns of bulbar output on internal connectivity rather than on input patterns: (2) the exquisite selective sensitivity of the bulb to learned inputs; and (3) its broad spectrum and aperiodic wave forms. The model was changed by including an asymmetric sigmoid nonlinearity and by introducing long feedback connections from the cortex to the bulb that had dispersive delays. These changes sufficed to provide chaotic solutions of the equations, and to selectively sensitize the system to destabilization by input, thereby improving the simulations of the biological activity. The conclusion is put forth that for pattern recognition to be done, the chaotic activity of the sensory system must arise from a macroscopic attractor. Local fluctuations are everywhere damped by the distributions of return times and strengths owing to biological variability in neuron size and connection density. Thus high-dimensional local stability is combined with low-dimensional global instability so as to emphasize the relationships among local features in a holistic sensory input pattern.

1. Introduction

The purpose of this report is to describe a way in which sensory information may be processed in a distributed mode by one of the simpler forms of cerebral cortex. The report is based on a mathematical model of the neural interactions of the mammalian olfactory bulb and cortex, which has been derived from anatomical and physiological studies of these two structures. These studies have been focused on the electroencephalographic (EEG) waves generated by the bulb and cortex in waking cats, rabbits, and other species during the performance of goal-directed behavior. The EEG in this condition shows a prominent sinusoidal burst at 40 to 80 Hz that lasts 100 to 200 ms and tends to recur with each act of inspiration. Because the bursts are absent or of low amplitude in sleeping or waking but non-motivated animals, and because the correlation between the EEGs of the bulb and cortex often rises to a high level during the EEG burst (Bressler, 1988), it is reasonable to believe that the processing of olfactory information by neural activity in the bulb, and its transmission to and reception by the neurons by the cortex, are closely related to the spatiotemporal patterns of EEG bursts.

It has been shown (Freeman, 1975) that the probability of pulse firing by single neurons in the bulb transmitting to the cortex oscillates at the frequency of the EEG burst activity when bursts are present. Each burst reflects a pattern of spatial coherence and organization of the activities of many neurons in the bulb during stimulation with odorants. Two questions are considered here. What are the bases for the neural dynamics leading to the burst activity? What is the possible role of the dynamics in the process of odor perception?

The approach was to construct a net of nonlinear integrodifferential equations in order to simulate the dynamics of the olfactory neural masses. The principal features of this model were taken from a lumped piece-wise linear model (Freeman, 1975) that was developed to describe the main neural interactions in the bulb and cortex. The nonlinear model was developed in stages, first as a lumped system (without the spatial variable), then as a uniform distributed system, and thereafter as a spatially nonuniform distributed system.

Fig. 1. The main types of neurons and connections in the olfactory bulb and cortex are shown schematically. See text for symbols. From Freeman (1975).

2. The lumped nonlinear model

The principal cell types and connections in the bulb and cortex are shown in Fig. 1. The receptors (R, excitatory) transmit by axons in the primary olfactory nerve (PON) to the mitral (M) and the tufted (T) cells in the bulb in synaptic clusters called glomeruli. M and T are excitatory to the main interneurons in the bulb, the granule cells (G), and they are inhibitory to the mitral and tufted cells (Shepherd, 1972). The interaction between M, T, and G forms a negative feedback loop. Additionally there are interactions among mitral and tufted cells (mutual excitation) and among granule cells (mutual inhibition). These main connections are shown schematically in Fig. 2. The output to the bulb is from the mitral and tufted cells by axons in the lateral olfactory tract (LOT). The output of the granule cells is reflected in the EEG generated by them.

The neurons of each type are considered to form a set, and in a local neighborhood a spatially homogeneous subset. Each subset is assumed to perform four operations in series. (a) It receives pulses and converts them to synaptic currents by an amplitude-dependent nonlinear process G_d(p). (b) The currents are summed and weighted over time by a linear time-invariant function f_d(t). (c) The resultant synaptic current is converted back to pulses by a nonlinear function G_a(v). (d) The pulses are transmitted by dispersive delay to the other subsets by a linear function F(t). For the present purposes G_d(p) is sufficiently linear to be replaced by a fixed coefficient k, and F_a(t) has the form F_a(t - T), where T is sufficiently small to be omitted, T= 0. Then F(t) and G(v) are to be defined.

Fig. 2. The connectivity of excitatory and inhibitory neurons in the bulb is shown schematically In the minimal form that includes mutual excitation, mutual inhibition, and negative feedback in a mixed population. From Freeman (1975).

Fig. 3. This nonlinear function for wave to pulse conversion at trigger zones places the maximal gain at the rest point of v = 0. This is incompatible with destabilization by input, because the gain is decreased with any input leading to output. The K's are fixed gains for lumped piece-wise linearapproximations. From Freeman (1975).

The time-dependent properties of the two main types of neurons in the open loop state have been determined by measurement of the responses of these neurons to single shock electrical stimulation. Precise simulation of these impulse responses requires at least a third-order linear differential equation. For present purposes, a second-order equation suffices as an approximation. For the ith set of neurons (Freeman, 1975),

where a = 220/sec + 720/sec and b = 220 x 720/sec². These rate constants correspond to the time constant for the passive membrane of the bulbar neurons (4.5 msec) and dendritic cable delay (1.4 msec). The nonlinear function G(v) has been found by determining the dependence of pulse probability of mitral and tufted cells on the wave amplitude of granule cells. This empirical relation is approximated by exponential function (Fig. 3) in which the coefficient b is included,

where P_o is a fixed coefficient representing the mean pulse density in the subset. Then from Fig. 2 and Eqs. (1) and (2) the dynamics of the bulbar set for an impulse input I∂(t) to the mutual and tufted cells is given by

It is important to note that in these equations the nonlinear operation precedes the summation at each node. The output of the model for low values of P_o(<3.8) consists of an oscillation that can be fitted with a damped sine wave,

For P_o = 0 the output gives the open loop response predicted by the solution to Eq. (1). For increasing values of P_o the frequency w_i increases and the decay rate a_i decreases, owing to the increase in the negative feedback gain (denoted "K" in Fig. 4). Examples of the root loci (-180^°) in the upper left complex plane near the origin are shown there with the gain contours given in decibels for a lumped piece-wise linear model (Freeman, 1975). For sufficiently high values of P_o (>1.8) the output is a limit cycle, which is stable owing to the saturation bilaterally of the sigmoid curve (Fig. 3). The amplitude varies with P_o and the frequency is almost invariant at 54Hz.

The adequacy of the model was tested by measuring the impulse response solutions in the left half of the complex plane (the zero equilibrium domain) with varying intensity of input I∂(t). The lumped nonlinear system behaved in the same way with respect to varying input intensity as did the lumped piece-wise linear model and the bulb and cortex to which that model conformed. This validated the use of the functions F(v) and G(v) for modeling the Hopf bifurcation.

Fig. 4. A root locus diagram is shown in the left upper quadrant of the complex plane for the Laplacian operator s = a ± jw in a linearized model of a negative feedback loop with two double poles at -220/sec and -720/sec. The root loci for negative feedback are at –180^° selected gain contour and the root locus. The threshold for instability lies at the intersection of the root locus and the jw axis at K = 0.52.

3. The spatially uniform distributed nonlinear model

The main problem in making the distributed model was to present local interactions within each subset in an array of subsets without introducing self-excitation and self-inhibition. These have been excluded because neurons are refractory during impulse transmission and cannot respond to their own input. Moreover, the likelihood of a neuron forming a synaptic connection onto itself in the cortex is vanishingly small, owing to the large number of neurons and their sparseness of connectivity. A viable solution is to construct an array of subsets, each of which has the connections shown in Fig. 2. Each excitatory element labeled "1e" represents a subgroup of neurons in a local neighborhood that receives input from the neighboring excitatory and inhibitory neurons (Fig. 5). A similar set connections is formulated for each inhibitory element labeled "1i". The strength of connection of two neighborhoods located at x and x' is considered to decrease with the distance at x - x' between those elements such that

In the computations d(x, x') has been approximated by a triangle with a half-width of 2s. Three coefficients are needed: s_n for negative feedback, s_e for mutual excitation, and (T_i for mutual inhibition. Experimentally the value for s_n for transmission between mitral-tufted and granule cells was approximately 0.4 mm. In the model this half-width is 5 elements, so the width of each element is 0.16 mm, and an array of 50 elements is 8 mm in width. (For computation the size of the element is ∆x' = 1.) These dimensions approximate the average interglomerular distance and the gross dimensions of the bulb. The values for s_e and s_i are estimated to be 5 from anatomical measurements.

Fig. 5. Input connectivity is shown schematically for an excitatory element in an array: a comparable connectivity holds for inhibitory elements as well. See Eqs. (7).

In order to avoid boundary effects the elements are connected in a circular array. Then from Eqs. (1)-(3) and (5) and the pattern of connections in Fig. 5,

where k_ij is a fixed coupling coefficient, and R = 50 is defined by the connectivity, including the circular boundary condition. Then,

For spatially uniform input I(x, t) = I_o∂(t) at every element and for P_O < 2 the output (impulse response) is again a damped sine wave. The root locus for increasing is essentially the same as for a lumped KII set, except that zero equilibrium solutions are found at higher values of P_o, and the difference increases with large values for u.

Above a certain value for P_o (depending on s) the solution is a stable limit cycle for which the frequency is almost invariant near 54 Hz and the amplitude depends on P_o. For excessive values of P_o the limit cycle becomes non-sinusoidal.

The spectral properties of the uniform system were determined in the equilibrium domain by giving an impulse input for which the amplitude varied with the spatial frequency f_x in cycles/min. The spatial cosine function was superimposed on a unit step function to simulate pure excitation,

The amplitude W_i of the damped sine wave output

of each element varied with x at the input spatial frequency. The attenuation of W_i with increasing spatial frequency was rapid. At f_x = 0.3/mm (about 3 cycles in an array of 50 elements) the spatial variation in output v_i was 1% of the output amplitude W_i at f_x = 0.

In the limit cycle domain the output amplitude was constant for all elements, so that the spatial spectrum was an impulse, whereas the temporal spectrum showed a sharp peak at 40 Hz. This finding showed that in the low-gain receiving mode, the pass bands of the spatial and temporal spectra are broad, so that high spatial and temporal frequencies of input corresponding to patterns of action potentials are accepted. However, after a state transition from equilibrium to limit cycle, a Hopf bifurcation based on increased internal gain, the pass bands are closed to all but the lowest spatial and temporal frequencies. In effect, the model receives and accepts patterned input only in the low-gain state, and it transmits its patterned output only in the high-gain oscillatory state. This is because the intracortical associational input in the diastolic low-gain state is less than the external input, allowing the latter to dominate, and, further, the disorganized and incoherent activity that typically results from input is smoothed out under spatial integration by the output pathway. In the high-gain state the spatial coherence is elevated, resulting in more effective transmission through the spatially divergent output path, and the controlling input to each transmitting neuron is more heavily intracortical (associational) than extracortical (sensory) in the systolic high-gain state.

4. The nonuniform distributed nonlinear model

Spatial variations in the gains and connectivity functions were introduced by varying P_o and k_ij as functions of x. Eqs. (2) and (6) were re-defined:

Changing P_o(x) caused a change at x in the slope of the nonlinear function for wave to pulse conversion at trigger zones as shown in Fig. 3.

The values for P_o(x) were varied sinusoidally across the array. The principal effect was to increase the amplitude w_i(x) and phase F_i(x) of the limit cycle activity for elements with high P_o(x) and to decrease amplitude and phase for elements with low P_o(x). The temporal frequency was constant across the array at w_o, and phase f was measured with respect to the temporal phase of the ensemble average of the outputs of all the elements. The mean amplitude of output V_o was close to the amplitude when P_o(x) was uniform and equal to the mean value for P_o(x)=^∆=P_o. The level for P_o at the transition from the zero equilibrium to the limit cycle domain was approximately the same as for P_O in the uniform case. The magnitude of variation in w_i(x) decreased with increasing spatial frequency of P_o(x) and approached 5% of V_o at 0.6/mm. That is, the limit cycle activity tended to become uniform when the spatial period of the variation in P_o(x) was less than twice the dispersion coefficient s_e.

P_o(x) was also varied by assigning random numbers and then normalizing the values to an appropriate mean P_o and variance. The limit cycle amplitude V_o and frequency w_o were not changed. Values for w_i(x) and j_i(x) depended on the variance of the setting of P_o(x) and had the shape of a smoothed version of P_o(x). The correlation between input P_o(x) and output w_i(x) ranged from 0.6 to 0.9, demonstrating the strong dependence of output on input for this model.

5. Modulation of the limit cycle activity by input

Characteristically with each inspiration the receptors deliver a volley of impulses to the bulb, and the background impulse activity P of the bulb increases and decreases in a phasic manner with the respiratory cycle. This fluctuation in P and its average value P_o are enhanced when the animal is motivated as by hunger or fear. When that occurs, the sinusoidal burst appears in the EEG. Because the interaction strengths (feedback gains) among elements in the model depended on the slope of the nonlinear sigmoid curve at v_i = 0 given by Eq. (2), and because that slope depended on P_o, it has been postulated (Freeman, 1975) that the bulbar system might be shifted with each inspiration from an equilibrium domain to a limit cycle domain, and that the sinusoidal EEG burst manifested that limit cycle. Presumably it was quenched by the fall in background impulse activity during each expiration. This postulate was expressed as the dependence of P_o(x) on the integrated receptor pulse input i(x, t) from the start of inhalation at t = 0 to the onset of the limit cycle activity at t = t_o, where t_o marked the time at which the Hopf bifurcation occurred:

where k_o was a fixed coefficient representing the strengths of the synaptic input connections to the excitatory neurons.

On the basis of single unit recording (Adrian, 1950; Lettvin and Gesteland, 1965; Eeckman and Freeman, 1990), there is good reason to suppose that the density of receptor input to the bulb is spatially nonuniform, and that it varies spatially in relation to particular odors being presented. If the input to the model i_i(x) for a given inspiration was not uniform, and if P_o(x) depended on i_i(x), then the threshold for the Hopf bifurcation to a limit cycle was not affected, and the output in the form of w_i(x) and j_i(x) were simple functions of the input.

When the values for k_ee(X) or k_i(x) were varied randomly in the presence of uniform i_i(x) and P_o(x), the domain of equilibrium stability was enlarged. That is, the system required a higher value of P_o to enter the limit cycle state than it did when the connectivities were uniform. The implication is that variation in synaptic coupling strengths, which are expected in a large population on the basis of biological variability, cause a distribution of characteristic frequencies to emerge among excitatory and inhibitory neurons coupled by negative feedback into local oscillators. The distribution has the effect of local damping of fluctuations. This is an important finding, because it means that an instability can arise only by global interactions. It is this feature that underlies the large-scale spatial coherence of activity as it is observed in multichannel records from the olfactory system (Freeman, 1991).

The opposite effect, destabilization, was achieved in the model when the connectivities k_ee(X) were spatially varied in a correlated way with the input function P_o(x). In this case the limit cycle state occurred for lower values of P_o. This was done in the model in two ways.

First, the mutually excitatory connections k_ee(x) were changed in relation to P_o(x) according to an empirical function:

For appropriate values of the coefficients a₁-a₃ the characteristic frequencies of the local subsets were tuned to a common resonance frequency w_o given by

Families of root loci from Eqs. (7) and (12) were calculated while both k_ee(x) and P_o(x) were changed uniformly across the array. For increasing k_ee(x) the characteristic temporal frequency w did not change and the decay rate a decreased, implying a decrease in stability.

Second, the inhibitory connection strengths k_i(x) were changed uniformly in relation to P_o(x) according to another empirical function,

so that the characteristic frequencies again conformed to a resonance locus given by Eq. (13). In this case the magnitude of k_i(x) decreased with larger P_o(x). The root loci for k_i showed that the characteristic frequency t_o, decreased with decreasing k_i. There was no change in stability. This result implies that the mutually excitatory connections are much more important in relation to destabilization than are the mutually inhibitory connections. It can further explain how small changes in k_ee(x) with learning (Freeman, 1975, 1987) can destabilize a neural mass and encourage the onset of oscillations.

When a certain pattern was fixed in the connectivities according to a set of random numbers for i_i(x) and Eqs. (11) and (12) or (14) and the same input was given, the system readily entered into a limit cycle at a relatively low frequency, for example, 40 Hz, for appropriate values of the a_ij's. For constrained values of k_ee(X) there was a strong positive correlation between Ii(x) and output V_i(x), whereas for constrained values of k_i(x) there was a strong negative correlation between I_i(x) and V_i(x). If an input I_j(x) unrelated to the connectivity pattern in the array was given, the frequency w_i either reverted toward the characteristic frequency of 54 Hz for the uniform system, or the tendency to limit cycle activity was reduced or abolished, and the correlations between input and output decreased. These properties showed that the system responded maximally to a pattern of input that correlated with a pre-existing connectivity pattern of strengthened reciprocal connections in the array, but the output pattern was dominated by the input pattern in all cases, and the output could be construed as completing the pattern of the interconnected elements in the array in response to partial inputs. This finding seemed to explain the prominence of "40 Hz" in olfactory cortex. An alternative conception has been advanced to explain neocortical 40 Hz oscillations observed at the macroscopic level in the human magnetoencephalogram and at the microscopic level in single neurons in the primate thalamus, in which single-cell oscillators are coupled to "resonate" over large areas of the brain (Llinas and Ribary, 1993). However, neither approach has sufficed to explain the olfactory electrophysiological data.

6. Improvement by modification in the sigmoid curve

Experiments with recording EEGs from spatial arrays of electrodes on the olfactory bulb and cortex of animals have shown that the spatial patterns of amplitude modulation of oscillatory bursts occur in relation to specific odorants that the animals have been trained to discriminate. The results have shown, further, that the spatial patterns of bulbar and cortical output are not determined by the spatial patterns of input, but by the internal connectivities of the bulb and cortex. An example is shown in Fig. 6. On the left is a set of 64 traces of a single burst on one inhalation showing the common wave form of the oscillation and the differing amplitudes across the array. The contour plots on the right show the spatial pattern differences between a control odorant state with no reinforcement and an odorant amyl acetate followed by a weak aversive electric shock to the skin, disagreeable but not painful. Although the odorants were unchanged from the 1st to the 3rd weeks of training, both patterns were modified, as the animal learned to adapt to the context of an aversive stimulus in the recording box, which it had not previously had to deal with.

Fig. 6. The 64 traces of a single burst unaveraged during an inhalation were recorded simultaneously from the surface of the olfactory bulb with an 8 X 8 array of electrodes spaced at 0.5 min as determined by the Nyquist frequency of the spatial spectrum (Freeman and Baird, 1987). The frames at the right show contour plots of the average root mean square (RMS) amplitudes of sets of ten bursts, which were recorded during training of the rabbit to identify the odor of amyl acetate (banana oil) under mild aversive reinforcement. The changes are summarized over the first three weeks of training. From Freeman and Schneider (1982).

This finding and many others like it showed that internal factors and not the particular stimuli determined the spatial patterns of bulbar output in the waking, motivated animals. A search for the mechanism revealed that Eqs. (2) previously derived for the sigmoid curve (Freeman, 1975) were incorrect, because the maximal slope was required by those equations to be at the rest state. Re-evaluation of the experimental data showed that this was not the case. Instead, the maximal slope in the data was always displaced to the excitatory side of the rest point of the sigmoid curve. A new sigmoid curve was derived from first principles (Freeman, 1979a), which fit the data more accurately, and which confirmed the finding of the displacement of the maximal gain (Fig. 7). According to the new asymmetric sigmoid curve the input not only causes an increase in pulse density P_o but it also increases the forward gain of each element and therefore the feedback gain between elements. The gain is therefore input-dependent. The average pulse density P_o does not change with inhalation, but it does change with the degree of motivation, as, for example, hunger increases with the duration of food deprivation. This property holds for all parts and cell types of the olfactory system (Eeckman and Freeman, 1991).

Fig. 7. The derivation of the asymmetric sigmoid curve is shown in terms of the intervening variable, m, which denotes the voltage-dependent sodium permeability of axons that gives the property of an exponential increase in sensitivity with depolarization of the membrane potential. From Freeman (1979a).

The finding was already established (Freeman, 1975) that when an animal learned to identify a new stimulus, the synaptic change was a selective increase in k_ee(i,j) according to a modified Hebb rule, where i and j represented two elements that were simultaneously excited by the input, and whose outputs were therefore highly correlated. An example is shown in Fig. 8 of simulated learning. Three of n elements received input in a training trial. The outputs of the oscillatory bursts were correlated, and the connections between those pairs having high amplitudes and high correlations were strengthened (upper frame). Thereafter a stereotypic output pattern of spatial amplitude modulation occurred irrespective of which element in the interconnected subset (a Nerve Cell Assembly in the Hebbian terminology) received the input. With this change in the model it conformed to the experimental findings, in which output was not a simple function of input, but it depended on the connectivities (Freeman, 1979c), more specifically on the selection of a Nerve Cell Assembly by input from among those stored previously by learning.

Fig. 8. Selective sensitization is achieved by increasing the connection strengths k_ee between elements having high correlation of their outputs during the formation of a burst on input. Thereafter the stereotypic output spatial pattern reflects the Nerve Cell Assembly and not the current input. From Freeman (1979b,c).

Further studies of the performance of the model during learning (Yao et al., 1991; Shimoide et al., 1993) have shown that associative Hebbian learning is insufficient, and that it must be accompanied by habituation, the attenuation of responsiveness of cortical elements to input that is undesired as being irrelevant, ambiguous or confusing. This non-Hebbian learning affects all prior stored Nerve Cell Assemblies in the process of removing overlap, so that the patterns of output are not invariant with respect to the classes of input. This property conforms well with the lack of invariance of the spatial patterns from the EEGs with respect to serial conditioning (Freeman, 1991).

The increase in the strength of the long-range, reciprocal synaptic connections between the excitatory neurons plays a crucial role in the selection of the basin of attraction to which the system transits upon induction of a Hopf bifurcation. This has been shown in a study by Edelman and Freeman (1990), in which the threshold for bifurcation was determined as a function of the strength of k_ee (Fig. 9). The selective sensitivity of the system to learned input is greatly enhanced by the operation of the Nerve Cell Assembly in conjunction with the asymmetric sigmoid curve, owing to the input-dependent increase in feedback gain. By this mechanism a learned stimulus embedded in background giving a low signal-to-noise ratio will trigger an explosive growth of activity just as the system is approaching a separatrix that will take its trajectory into a new basin of attraction. This study also revealed that the sensitization was specific for the excitatory connections to the exclusion of modifiable inhibitory connections, and, moreover, that the use was inadmissible of self-excitation to compensate for the omission of mutual excitation or for the use of a first-order ODE in place of a second- or higher-order ODEs to include the effects of dendritic cable delay as in Eq. (1).

Fig. 9. The output is plotted as a function of input for the circuit shown in Fig. 2 for two values of k _ee. The bimodal character of the response at high values of the excito-excitatory gain parameter is characteristic of this system. Its importance lies in the fact that during learning, this is the parameter that is increased (Freeman, 1968, 1975, 1979b,c). It is responsible for the selective sensitivity of the chaotic system to learned input. When the Nerve Cell Assembly receives weak input that is masked by background input to give a low signal-to-noise ratio, the early onset of explosively growing activity in the assembly carries the bulb across a separatrix at the onset of the state transition, and guides the entire bulb into the relevant wing of the global attractor. From Edelman and Freeman (1990).

7. Improvement by modifications in the connectivity

This model failed to simulate the aperiodic oscillations of the EEGs of the olfactory system (Fig. 10) and the forms of their spectra, which are broad rather than peaked and tend to the "1/f" pattern characteristic of intermittent chaos. The limit cycle forms of EEG activity are observed in the bulb or cortex only after they have been surgically isolated from each other and been given a sustained excitatory chemical stimulus. Three modifications by addition of new elements and pathways (Fig. 12) were required in order to correct this deficiency (Freeman, 1987).

Fig. 10. Examples are shown of the chaotic background EEG activity that is generated by the olfactory bulb (OB) and prepyriform cortex (PC) and the simulations from the model of these outputs, and those of the mitral (M) cells and the AON as well. From Freeman (1987).

First, a mutually excitatory population at the input of the bulb was required to simulate the actions of interneurons there that are known as the periglomerular cells (P) owing to their juxtaposition with the synaptic nests called glomeruli (Fig. 1), in which the sensory axons make their connections onto the mitral and tufted cells of the bulb. The set of P cells is self-stabilizing and provides a strong and continuing excitatory bias that maintains the olfactory system in an oscillatory dynamic range (Freeman, 1975). The P cells also receive input from sensory receptors, which excites them and enhances the excitatory bias, so they play a role in the destabilization of the bulb by enhancing the increase in gain through the asymmetric sigmoid function (Freeman, 1979b). Furthermore, one of the targets to which the bulb projects, the anterior olfactory nucleus (AON), sends axons back to the P cells forming a long excitatory feedback loop. This loop is critical for the induction of chaotic activity by the olfactory system (Freeman, 1987). particularly the broad "1/f" type spectrum that emphasizes the amplitude in the low end of the spectrum, which gives rise to low-frequency, high-amplitude fluctuations as exemplified in the time series, simulated below and recorded above, in Fig. 10, and the large, irregular loops in phase portraits (Fig. 11). The model suffices to replicate state transitions to burst activity on the presentation of input that simulates the normal sensory input with repeated inhalations (Freeman, 1987; Kay et al., 1993).

Second, the anterior olfactory nucleus (AON) and prepyriform cortex (PC) were included in an enlarged model (Fig. 12), in recognition that the olfactory system requires that these three parts be intact and in direct feedback relations, in order to produce the broad spectrum activity that is normally observed (Freeman, 1975). Each has its characteristic frequency as measured from its impulse response in experimental conditions where the interactions become negligible, namely, in forming a time ensemble average by repetitive stimulation at a low rate to yield an averaged evoked potential (AEP). These frequencies are amplitude-dependent, but on extrapolation into the low-amplitude, small-signal range of linearity, it becomes apparent that the three frequencies are incommensurate, meaning that they are not in integer ratios. This is an important property, because it means that the three parts cannot agree on a single frequency, even though they are locked into a multiple loop feedback system.

Fig. 11. Phase portraits are shown of the type of activity displayed in Fig. 10 to compare the real and simulated EEGs. Rotation is counterclockwise in each frame, because the mitral (M) activity precedes the others by one quarter cycle of the oscillations on the average. The traces are 500 msec in duration, digitized at 1 msec intervals. From Freeman (1987).

Fig. 12. Modifications of the connectivity include the introduction of the periglomerular (P) cells to provide a modulatory excitatory bias control; incorporation of three oscillators comprising the olfactory system, the bulb (OB), anterior olfactory nucleus (AON), and prepyriform cortex (PC) having noncommensurate characteristic frequencies; and the necessity for distributed delays in the long feedback paths that interconnect the three oscillators and the P cells. These elements and their connections sufficed to simulate the chaotic properties of the EEGs of the three structures. From Yao and Freeman (1990).

Third, delays are introduced into the four main feedback paths, which are due to the slow conduction velocities of the small neurons in the feedback pathways and amount to between a tenth and a half cycle duration of the characteristic frequencies. Owing to variations in axon length and diameter, the arrival times are dispersed, so that each path acts as a low pass filter, thereby enhancing the low-frequency activity that is fed back.

Among the feedback pathways there are three from excitatory to inhibitory neurons, which constitute stabilizing negative feedback, because the delays up to half a cycle are not sufficiently dense to lead to reinforcement of a negative going wave. The fourth path is from excitatory neurons in the AON to excitatory P cells and constitutes a positive feedback pathway, which, if strong enough, can be destabilizing. Analysis of the KIII model shows that this loop is the most likely site of a positive Lyapunov exponent calculated from the real and simulated EEGs. It therefore appears that the aperiodic activity of the olfactory system is a global property, and it is not due to local instabilities that are prevented by the distribution of feedback gains within each of the three main parts. This is important in respect to the state change that takes place with each inhalation, because it means that the entire olfactory system participates as a macroscopic entity by jumping from one wing of the global attractor to another. Only in this way could the entire system perform pattern recognition on the spatially and temporally distributed input patterns that constitute sensory stimuli. It means, further, that the state transition is not a Hopf bifurcation, but it is a step in a continuing evolution that Tsuda (1991) has described in terms of "chaotic itinerancy".

The solutions of the equations that serve to simulate the statistical, spectral, and dynamical properties of the EEGs have been derived by numerical integration. The introduction of the distributed delays is required by the physiology being modeled. Unfortunately, this type of differential delay equations cannot be solved analytically by presently known techniques, and in the chaotic domains the numerical solutions fail to be robust, owing to sensitivity to small changes in the initial conditions and in the parameter values. Seemingly trivial changes in parameters may lead to unexpected outcomes. Hence further development of neurodynamics for the cerebral cortex may depend heavily on prior developments in the requisite mathematics.

Acknowledgement

This work was supported by grants from the National Institute of Mental Health and the Office of Naval Research. A preliminary report has been given from which this work has been adapted (Freeman and Ahn, 1976).

References

E.D. Adrian, The electrical activity of the mammalian olfactory bulb, Electroenceph. Clin. Neurophysiol. 2 (1950). 377-388.

S.B. Bressler, Changes in electrical activity of rabbit olfactory bulb and cortex to conditioned odor stimulation., Behav. Neurosci. 102 (1988) 740-747.

J.A. Edelman and W.J. Freeman, Simulation and analysis of a model of mitral-granule cell population interactions in the mammalian olfactory bulb, in: Proc. IJCNN 1 (1990), pp. 62-65.

F.H. Eeckman and W.J. Freeman, Correlations between unit firing and EEG in the rat olfactory system, Brain Res. 528 (1990) 238-244.

F.H. Eeckman and W.J. Freeman, Asymmetric sigmoid nonlinearity in the rat olfactory system, Brain Res. 557 (1991), 13-21.

W.J. Freeman, Mass Action in the Nervous System (Academic Press, New York, 1975).

W.J. Freeman, Nonlinear gain mediating cortical stimulusresponse relations, Biol. Cybern. 33 (1979a) 237-247.

W.J. Freeman, Nonlinear dynamics of paleocortex manifested in the olfactory EEG, Biol. Cybern. 35 (1979b) 21-37.

W.J. Freeman, EEG analysis gives model of neuronal template-matching mechanism for sensory search with olfactory bulb, Biol. Cybern. 35 (1979c) 221-234.

W.J. Freeman, Simulation of chaotic EEG patterns with a dynamic model of the olfactory system, Biol. Cybern. 56 (1987), 139-150.

W.J. Freeman, The physiology of perception, Sci. Amer. 264 (1991) 78-85.

W.J. Freeman and S.M. Ahn, Spatial and temporal characteristic frequencies of interactive neural masses, in: Proc. IEEE Int. Conf. on Cybernetics and Society Washington, D.C., November 1976 (1976) pp. 279-284.

W.J. Freeman and B. Baird, Relation of olfactory EEG to behavior: Spatial analysis: Behav. Neurosci. 101 (1987) 393-408.

W.J. Freeman and W. Schneider, Changes in spatial patterns of rabbit olfactory EEG with conditioning to odors, Psychophysiol. 19 (1982) 44-56.

L. Kay, S. Shimoide and W.J. Freeman, Comparison of EEG time series from rat olfactory system with model composed of nonlinear coupled oscillators (1993), in preparation.

J.Y. Lettvin and R.C. Gesteland, Speculations on smell, Cold Spring Harbor Symp, Quant. Biol. 30 (1965) 217-225.

R. Llinas and U. Ribary, Coherent 40-Hz oscillation characterizes dream state in humans, Proc. Nat. Acad. Sci. USA 90 (1993) 2078-2081.

G.M. Shepherd Synaptic organization of the mammalian olfactory bulb, Physiol. Rev. 52 (1972) 864-917.

K. Shimoide, M.C. Greenspon and W.J. Freeman, Modeling of chaotic dynamics in the olfactory system and application to pattern recognition, in: Eeckman F.H., Neural Systems Analysis and Modeling (Kluwer, Boston, 1993) pp. 365-372.

I. Tsuda, Chaotic itinerancy as a dynamical basis of hermeneutics in brain and mind, World Futures 32 (1991) 167-184,

Y. Yao and W.J. Freeman, Model of biological pattern recognition with spatially chaotic dynamics, Neural Networks 3 (1990) 153-170.

Y. Yao, W.J. Freeman, B. Burke and Q. Yang, Pattern recognition by a distributed neural network: An industrial application, Neural Networks 4 (1991) 103-121.