BioSystems 58 (2000) 249 – 257 www.elsevier.com/locate/biosystems
Dynamics of a hybrid system of a brain neural network and an artificial nonlinear oscillator Norihiro Katayama *, Mitsuyuki Nakao, Hiroaki Saitoh, Mitsuaki Yamamoto Laboratory of Neurophysiology and Bioinformatics, Graduate School of Information Sciences, Tohoku Uni6ersity, Aobayama 05, Sendai 980 -8579, Japan
Abstract In the brain, many functional modules interact with each other to execute complex information processing. Understanding the nature of these interactions is necessary for understanding how the brain functions. In this study, to mimic interacting modules in the brain, we constructed a hybrid system mutually coupling a hippocampal CA3 network as an actual brain module and a radial isochron clock (RIC) simulated by a personal computer as an artificial module. Return map analysis of the CA3-RIC system’s dynamics showed the mutual entrainment and complex dynamics dependent on the coupling modes. The phase response curve of CA3 was modeled regarding the CA3 as a nonlinear oscillator. Using the phase response curves of CA3 and RIC, we reconstructed return maps of the hybrid system’s dynamics. Although the reconstructed return maps almost agreed with the experimental data, there were deviations dependent on the coupling mode. In particular, we noted that the deviation was smaller under the bidirectional coupling conditions than during the one-way coupling from RIC to CA3. These results suggest that brain modules may flexibly change their dynamical properties through interaction with other modules. © 2000 Elsevier Science Ireland Ltd. All rights reserved. Keywords: Hippocampal neurons; Field potential; Coupled oscillators; Phase resetting; Mutual entrainment; Synchronization
1. Introduction The brain is a complex system that consists of many functional modules. Higher-order brain functions are thought to be realized though inter-
* Corresponding author. Tel.: +81-22-2177179; fax: + 8122-2639438. E-mail address:
[email protected] (N. Katayama).
actions among these modules. For example, functional modules in the visual cortex have been revealed to work in parallel and hierarchically to process complex information (Felleman and Van Essen, 1991), and synchronous neuronal activities in the different modules have been anticipated to play an important role in solving the ‘binding problems’ in the brain (Milner, 1974; von der Malsburg and Schneider, 1986). This suggests practical aspects of the inter-module interactions. Therefore, in order to understand the mechanisms underlying the brain’s information processing, it
0303-2647/00/$ - see front matter © 2000 Elsevier Science Ireland Ltd. All rights reserved. PII: S 0 3 0 3 - 2 6 4 7 ( 0 0 ) 0 0 1 2 9 - 5
250
N. Katayama et al. / BioSystems 58 (2000) 249–257
is essential to know the dynamical characteristics of each module, as well as the nature of interactions among the modules. In this study, as a framework for investigating the dynamics of a multi-module system, we constructed a hybrid system consisting of a neural network in the hippocampal CA3 region as an actual brain module and a nonlinear oscillator as an artificial module working together as a minimal unit representing interacting modules in the brain. In this framework, we analyzed the dynamics of the hybrid system under varied coupling conditions. Considering the close relationship between the dynamical properties and the function of the modules, variation of the dynamics of the hybrid system implies diversity in the modules’ functions. For considering these variations, we investigated how interactions determine the function of each module in a multi-module system.
2. Materials and methods
2.1. Brain module: hippocampal CA3 region We adopted the neural network in the CA3 region of a hippocampal slice preparation (abbreviated CA3) as the brain module. Transverse slices (0.4 mm thickness) of hippocampi were obtained from the brains of guinea pigs (300920 g) by the conventional method (Ito et al., 1989) and kept in normal artificial cerebrospinal fluid (ACSF) containing (in mM) NaCl 124, KCl 3.5, NaH2PO4 1.25, MgSO4 2.0, CaCl2 2.0, NaHCO3 22.0, and glucose 10.0 (pH 7.4), which was equilibrated with a gas mixture of 5% CO2 in O2. All the experiments were conducted in a perfusing fluid chamber kept at 349 1°C. In order to detect neural activities, we measured the extracellular field potential at the CA3 pyramidal cell layer with a single glass pipette filled with the ACSF. The potential was amplified through the bandpass filter (0.3–25 Hz) and sampled with 200 Hz by a personal computer (IBM-PC compatible, CPU AMD K6-233; abbreviated PC) with an A/D and D/A converter card (National Instruments, AT-MIO-16E-2, 12-bit precision). The PC simultaneously simulated the artificial module (see
below), which provided electric stimuli to the mossy fiber (MF) bundle with a bipolar tungsten electrode. Before starting experiments, pulse-current stimuli with an amplitude of 0.5–1 mA and duration of 200 ms were applied every 30 s to MF. Only the slices showing stable responses to those stimuli for 30 min were used for the following experiments. The CA3 neurons are known to have recurrent excitatory synaptic connections, which tend to cause a rhythmic neural activity. We adopted a CA3 neural network as a model of a rhythmicacting brain module. Actually, when the ACSF was replaced with an experimental ACSF containing 8.5 mM KCl and 2 mM penicillin chloride, the CA3 neurons exhibited spontaneous bursts (Hayashi and Ishizuka, 1995). The CA3 always responded with a burst to a stimulus stronger than a threshold level (0.03–0.07 mA). When threshold-level stimuli were applied, the neural responses were strongly dependent on the temporal pattern of the stimuli. In order to study the organizing process of the response, we used threshold-level stimuli.
2.2. Artificial module: radial isochron clock As the artificial module, we used a radial isochron clock (RIC) model to mimic a rhythmicacting brain module, because it is the simplest model of a biological oscillator (Winfree, 1980; Glass and Mackey, 1988; Nomura et al., 1994). Fig. 1b shows the dynamics of the RIC, which is described by the following equations: r; = Kr(1− r) u: =
1 TRIC
(1a) (1b)
where r= x 2 + y 2, u is the normalized ‘phase’ of the RIC’s state point ranging between 0 and 1, TRIC is the period of the limit cycle, and K is a constant. In this study, we set TRIC = 5 s and K . When the RIC is stimulated with an impulse, the RIC’s state point (indicated by P in Fig. 1b) is moved toward the direction of the x-axis responsible for the stimulus intensity A (P%). Then, the state point instantaneously ap-
N. Katayama et al. / BioSystems 58 (2000) 249–257
proaches the point on the limit cycle (P%) along with the isochron because of K . The relationship between the old phase u at P and the new phase u% at P% is described below: u% =
1 sin 2pu arctan 2p A+ cos 2pu
u +DRIC(u)
251
inhibitory one (Nomura et al., 1994). In the case of A= 0, corresponding to one-way coupling from RIC to CA3, periodic stimulation with an interval of TRIC is realized.
3. Results (2)
where DRIC is called the phase response curve (PRC) of the RIC. In practice, RIC was simulated by the PC, which executed real-time processing of Eq. (1)b Eq. (2). When the RIC’s phase passed across zero, the RIC was considered to have fired. At this moment, the PC generated a pulse by the D/A converter to stimulate the CA3. On the other hand, when a CA3 bursting was detected, the RIC was stimulated by an impulse with an amplitude of A. The occurrence of a CA3 bursting was detected when the band-pass-filtered field potential crossed a threshold level (about 0.1 mV). Time series of the CA3 burstings and RIC firings were recorded under varied coupling modes, each for 10 min, where intensity of the stimulus to the RIC was changed from A = 0.95, −0.95 and then A = 0. If one views the RIC as a single pacemaker neuron, a positive A corresponds to an excitatory coupling, whereas a negative A corresponds to an
3.1. Dynamics of inter-bursting inter6als Fig. 2a–c show the bursting activities of the CA3 and the firings of the RIC under the varied coupling conditions. The CA3-RIC system exhibited diverse dynamics that were dependent on the coupling modes. First, we investigated the dynamics of the inter-bursting intervals of the CA3 ICA3 and the RIC IRIC (Fig. 3, see also Fig. 2d for definitions). In the case of one-way coupling from RIC to CA3 (A= 0, Fig. 3b), IRIC remained unchanged while ICA3 fluctuated significantly (interval= 3.499 0.67 s, mean9 S.D). When the stimulus intensity for the CA3 was increased well above the threshold, CA3 burstings were synchronized to the RIC firings with a ratio of 1:1 (1:1 forced-entrainment, not shown here). This seemed to indicate that the stimulus intensity used was not sufficiently strong to entrain the CA3 burstings.
Fig. 1. (a) Schematic representation of the hybrid system consisting of the brain module and the artificial module (the CA3-RIC system). Field potentials in the CA3 region were amplified and sent to a personal computer (PC), which simultaneously simulated a nonlinear oscillator (RIC). When a CA3 bursting potential was detected by the PC, the RIC was stimulated. On the other hand, the RIC generated an electric stimulus for the CA3 at the moment when the RIC fired (u = 0). (b) Dynamical properties of the RIC on the phase plane.
252
N. Katayama et al. / BioSystems 58 (2000) 249–257
Fig. 2. Observed dynamics of the CA3 and the RIC. (a–c) Field potentials of the CA3 (upper) and the output sequences of the RIC (lower) under the varied coupling modes. (d) Definitions of the inter-bursting intervals of the CA3 ICA3 and the inter-firing interval of the RIC IRIC. (e) Definition of the time interval c.
When CA3 and RIC were mutually coupled (i.e. A" 0), fluctuations of ICA3 seemed to be reduced. As shown in Fig. 3a and d, ICA3 for A = 0.95 was more stable than the other cases (interval= 2.9190.36 s for A =0.95 and 3.19 9 0.38 s for A= −0.95). Data points on the return map of ICA3 were shown to mostly concentrate around period= 3 s, but they occasionally deviated from the period. On the other hand, IRIC changed frequently between two distinct intervals, about 4 and 5.5 s. Analysis of the return map did not allow us to determine the switching cycle between the intervals. In the case of A = −0.95, RIC fired with period 2 (Fig. 3c and i), which was known from the return map IRIC( j+ 2) −IRIC( j ) consisting of two clusters on the diagonal line (not shown here). In contrast, periodicity of CA3 burstings was not so clear. Nevertheless, by comparing ICA3 and IRIC, we obtained the following integer ratio: mean firing interval of CA3:3.19 s 2 : . mean firing interval of RIC:4.83 s 3 In addition, data points on the return map of ICA3(i +3) versus ICA3(i ) were concentrated on a diagonal line (not shown). These results suggest that CA3 and RIC were mutually entrained with a period ratio of 3:2. As we mentioned in the case of one-way coupling, the stimulus for CA3 was not sufficiently strong for inducing forced-entrainment. However, even with such a weak stimulus, it is possible that CA3 and RIC cooperate to
synchronize their own activities through mutual interactions.
3.2. Dynamics of the latency of CA3 burstings from RIC firings In order to determine the relationship between the CA3’s dynamics and the RIC’s ones, we analyzed the time interval between the RIC firing and the latest CA3 bursting (see Fig. 2e for a diagram). Fig. 4a (basic plot of c) shows that c’s dynamics are seemingly complex in the case of A= 0.95. On the corresponding return map (Fig. 4d), one can clearly recognize five clusters. Since the clusters approximately form a single-valued function, the complex dynamics would be regarded as deterministic. In the case of A= − 0.95, c’s dynamics are almost period 3 because the return map consists mainly of three clusters (Fig. 4f). This result also supports 3:2 mutual entrainment between the CA3 and the RIC. In contrast, c’s dynamics for A=0 seem to be quasi-periodic (Fig. 4b). Nevertheless, as shown in Fig. 4e, distribution of the data points on the return map are qualitatively similar to the case of A = −0.95.
3.3. Modeling the phase response cur6e of CA3 To model the phase response curve of the CA3, we explored the stimulus-to-response relationship
N. Katayama et al. / BioSystems 58 (2000) 249–257
of the CA3 to understand the dynamics of c. We assumed the CA3 as a nonlinear oscillator, because it generated bursts repetitively with an approximately fixed period TCA3 if no stimulus was applied. In addition, even if the bursting interval was changed by some stimulus, the interval was restored by itself. Therefore, we investigated the phase response curve (PRC) of the CA3, which is an analyzing tool of internal structure of a nonlinear oscillator. Fig. 5 shows the relationship between the bursting interval T and the stimulus timing d, which
253
were obtained from the data presented in Fig. 3. As shown in the figure, overall response characteristics were consistent through the coupling modes. Regardless of the coupling modes, the bursting interval T(d) changed dramatically, depending on d. That is, if a stimulus was applied within TR after the last CA3 bursting (dBTR), the stimulus did not affect the bursting interval. In contrast, if a stimulus was applied after the refractory period TR, CA3 generated a burst immediately TD after the stimulus. We estimated (TCA3, TR, TD) as (3.10, 1.50, 0.10) for A=0.95,
Fig. 3. Dynamics of the CA3-RIC system under the varied coupling modes. (a – c) Basic plots of the inter-bursting intervals of CA3 (ICA3, ) and the inter-firing interval of the RIC (IRIC, ) versus event time. (d – f) Return maps of ICA3 corresponding to the data (a–c). (g – i) Return maps of IRIC corresponding to the data (a – c).
N. Katayama et al. / BioSystems 58 (2000) 249–257
254
Fig. 4. Return maps of the time-interval c under the varied coupling modes. (a – c) Basic plots of c. (d – f) Return maps of c. Dots indicate the experimental data. Thick lines are the reconstructed return maps; thin lines represent the stable orbits allowed in the reconstructed return maps. The entrainment ratio for the stable orbit is indicated beside the orbit.
(3.70, 1.80, 0.10) for A =0, and (3.42, 1.80, 0.10) for A= − 0.95 using the least-square method. Based on these results, we modeled the CA3’s PRC DCA3 as follows (Fig. 5e): DCA3(f)=
=
(TCA3 −T(TCA3f)) TCA3
!
0, (TCA3 − TD)/TCA3 − f
05f BTR/TCA3 otherwise,
where f(=d/TCA3) is the normalized ‘phase’ of the CA3 ranging between 0 and 1.
3.4. Reconstructed return map of c based on the phase response cur6es If the CA3’s PRC fully describes the dynamical properties of CA3 under the varied coupling conditions with RIC, the dynamics of c could be reconstructed based on the PRCs of both oscillators. Conversely, if there were disparities between the experimentally obtained and the reconstructed return maps; they indicate a nature of the CA3RIC system’s dynamics that could not be de-
scribed within the framework of a coupled two-oscillator system. The reconstructed return maps of c are presented in Fig. 4d–f together with the experimental data. The mathematical expression of the maps is described in Appendix A. As shown in the graphs, the reconstructed maps almost agree with the experimental data. The stable orbits for the reconstructed maps indicate that the CA3 and RIC would mutually entrain with period ratios of 11:7 for A= + 0.95, and of 3:2 for A= 0 and for A= − 0.95. Regardless of coupling conditions, the reconstructed orbit almost covers the experimental data. Therefore, the actual orbit is considered to wander in the neighborhood of the stable orbit predicted by the reconstructed map. Nevertheless, there are deviations between the reconstructed map and the experimental data, especially in the case of one-way coupling, which is depicted by an oval in Fig. 4e. We believe the deviation indicates fluctuations in the intrinsic bursting period TCA3 just after the CA3 bursting immediately evoked by the stimulus from the RIC, because we know the concrete cluster is
N. Katayama et al. / BioSystems 58 (2000) 249–257
255
Fig. 5. Stimulus-to-response characteristics of the CA3 under the varied coupling modes. (a – c) Plots of inter-bursting interval T as a function of the stimulus timing d. (d) Definition of the inter-bursting interval T and the stimulus timing d. (e) Modeled phase response curve of the CA3 (DCA3).
constructed by the point (c(n), c(n + 1)) = (TD,TCA3 +TD) according to the theory (see Appendix A). Actually, the fluctuation was significantly large in the case of one-way coupling. These results suggest that the dynamical properties of fluctuation in TCA3 changed dependent on the coupling mode.
4. Discussion We studied the dynamics of a hybrid system composed of an actual brain module (CA3 neural network) and an artificial module (RIC) under varied coupling conditions. Within our experimental conditions, the CA3 and RIC were mutually entrained when they were connected bidirectionally, whereas the CA3 was not firmly entrained to the RIC under the one-way coupling. The dynamics of the hybrid system were reconstructed by combining the phase response curves of the CA3 and the RIC. The reconstructed dynamics almost agreed with the experimental data except for the case of the one-way coupling. We
believe the exceptional disparity was due to fluctuations in the intrinsic bursting period of the CA3. Bursting activity of the CA3 in response to periodic stimulation corresponding to the oneway coupling has been investigated by Hayashi and Ishizuka (1995), who observed diverse behaviors such as phase-locking and chaotic activities dependent on the stimulus intensity and frequency. Our experiments showed that CA3 burstings were not entrained to the periodic stimuli if the stimulus intensity was not sufficiently strong for inducing forced-entrainment. On the other hand, even with the weak stimulation, mutual entrainment between the CA3 and RIC was realized under the bidirectional coupling condition. In this case, the CA3 burstings would modulate the phase of the RIC, thereby controlling the phase of stimuli coming from the RIC. This might be useful for synchronizing the modules. In addition, complex dynamical properties inherent in the CA3 network could also underlie these phenomena, although we treated the CA3 as a simple oscillator in this study. In any case, these results suggest that cooperative interactions naturally
256
N. Katayama et al. / BioSystems 58 (2000) 249–257
emerge in a framework such as the one we used, and such cooperation may be substantial in inter-module interactions in the brain. (Felleman and Van Essen, 1991). The hippocampus is known to couple with many modules in the brain, including the septum and entorhinal cortex. Therefore, the artificial module coupled with the CA3 should have dynamical properties similar to the properties of those modules. In spite of this, we cannot compare the RIC that we adopted in this study with those modules, because properties such as phase response curves have not been clarified yet. However, it is likely the RIC shares some dynamical properties with those brain modules, because many neural systems (e.g. thalamocortical neurons (Soltesz and Crunelli, 1992) and pacemaker networks (Ayers and Selverston, 1977) were reported to have PRCs similar to that of the RIC under certain conditions. We await future studies that will reveal the dynamical properties of the modules. From the viewpoint of biosystem control, our framework is useful. Recently, the chaos control technique has shown to be successful for decreasing as well as increasing periodicity of the CA3 bursting behavior (Schiff et al., 1994). This technique, like ours, is based on the mutual interaction between biological tissue and an artificial module. In our study, the rhythmic activity of a brain module was modulated by interacting with an RIC. Since the RIC is simpler than the chaos control algorithm, our framework may be advantageous for biosystem control.
Acknowledgements Experiments were authorized by the committee for experiment on animal of Graduate School of Information Sciences, Tohoku University. This work was partially supported by Grants-in Aids for Scientific Research from the Ministry of Education, Science, Sports, and Culture of Japan, No.10750329, No.11145207, and No.10480080.
Appendix A. Reconstruction of the CA3-RIC system dynamics In order to describe the dynamics of the CA3RIC system simply, we introduce several definitions as follows. The dynamics of the CA3-RIC system are represented by C for the CA3 bursting and R for the RIC firing. For example, when a firing of RIC is followed by the double bursts of CA3, this sequence is represented by RCC. To describe the return map concisely, we define the following functions: fRIC(t) t+ TRICDRIC(t/TRIC), fCA3(d) d+TCA3DCA3(d/TCA3) =
!
d, TCA3 − TD,
05 dB TR otherwise,
where TRIC and TCA3 are the intrinsic firing and the bursting periods, respectively, t and d represent the time phase, i.e. t= TRICu for RIC and d= TCA3f for CA3, where u and f are the normalized ‘phases’ of the respective module ranging between 0 and 1. DRIC and DCA3 represent the PRCs of the respective modules. Evidently, 05 fRIC B TRIC and 05 fCA3 B TCA3 hold. fCA3 is a non-decreasing function because of TCA3 − TD \ TR according to the experimental results. fRIC is monotonically increasing, because we set A B 1 in this study. fRIC and fCA3 give the new time phase of the respective system just after a stimulus coming from the other system. For example, when the RIC is stimulated at t=TRICu, the new time phase t%= TRICu% is instantaneous after the stimulus is given by t%= fRIC(t). Hence the time interval between the stimulus and the following RIC firing is described by (TRIC − fRIC(t)), if additional stimuli were not given within the interval. Using these definitions, we can describe the return maps of c as follows: c(n+ 1)=c(n)+ TCA3, 1 05 c(n)Bf − RIC(TRIC − TCA3)
(a.1)
N. Katayama et al. / BioSystems 58 (2000) 249–257
257
Fig. 6. Patterns observed in the output sequence of the CA3-RIC system. Vertical arrows represent the events (firings or burstings) of the module. (a) Sequence of RCC, (b) RCRC, and (c) RCCRC.
c(n +1)=TCA3 − fCA3(m1), 1 f− RIC(TRIC −TCA3) 5c(n) BTCA3
(a.2)
c(n + 1)=TCA3 − fCA3(m2), 1 TCA3 5 c(n)B TCA3 +f − RIC(TRIC −TCA3)
(a.3) 1 where f − RIC represents the inverse function of fRIC, m1 = TRIC − fRIC(c(n)) and m2 =TRIC −fRIC ( fRIC(c(n)− TCA3)+ TCA3). The dynamics of the CA3-RIC system for the respective cases Eqs. (a.1), (a.2) and (a.3) are depicted in Fig. 6. According to the reconstructed return maps, within our experimental conditions, each output sequence contains only RC and RCC patterns; in other words, combinations such as RR and CCC do not appear in the sequence if no behavior deviating from the map occurs. Due to paucity of space, we cannot explain the reason here.
References Ayers, J.L., Selverston, A.I., 1977. Synaptic control of an endogenous pacemaker network. J. Physiol.(Paris) 73,
.
453 – 461. Felleman, D.J., Van Essen, D.C., 1991. Distributed hierarchical processing in the primate cerebral cortex. Cerebral Cortex 1, 1 – 47. Glass, L., Mackey, M.C., 1988. From Clocks to Chaos: The Rhythms of Life. Princeton University Press. Hayashi, H., Ishizuka, S., 1995. Chaotic response of the hippocampal CA3 region to a mossy fiber stimulation in vitro. Brain Res. 686, 194 – 206. Ito, K.H., Miyakawa, H., Izumi, Y., Kato, H., 1989. Longterm potentiation induced by repeated paired stimulation of the stratum radiatum and alveus in CA1 pyramidal neurons of hippocampal slices. Biomed. Res. 10, 111 – 124. Milner, P.M., 1974. A model for visual shape recognition. Psychol. Rev. 81, 521 – 535. Nomura, T., Sato, S., Doi, S., Segundo, J.P., Stiber, M., 1994. A modified radial isochron clock with slow and fast dynamics as a model of pacemaker neurons. Biol. Cybern. 72, 93 – 101. Schiff, S.J., Jerger, K., Duong, D.H., Chang, T., Spano, M.L., Ditto, W.L., 1994. Controlling chaos in the brain. Nature 370, 615 – 620. Soltesz, I., Crunelli, V., 1992. A role for low-frequency, rhythmic synaptic potentials in the synchronization of cat thalamocortical cells. J. Physiol. (Lond.) 457, 257 – 276. von der Malsburg, C., Schneider, W., 1986. The neural cocktail-party processor. Biol. Cybern. 54, 29 – 40. Winfree, A.T., 1980. The Geometry of Biological Time. Springer, New York.