Next Article in Journal
Rigid Shape Registration Based on Extended Hamiltonian Learning
Next Article in Special Issue
Chaos Control and Synchronization of a Complex Rikitake Dynamo Model
Previous Article in Journal
Entropy, Free Energy, and Work of Restricted Boltzmann Machines
Previous Article in Special Issue
Neural Computing Enhanced Parameter Estimation for Multi-Input and Multi-Output Total Non-Linear Dynamic Models
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Birhythmic Analog Circuit Maze: A Nonlinear Neurostimulation Testbed

by
Ian D. Jordan
1,2 and
Il Memming Park
1,2,3,*
1
Department of Applied Mathematics and Statistics, Stony Brook University, Stony Brook, NY 11794, USA
2
Institute for Advanced Computing Science, Stony Brook University, Stony Brook, NY 11794, USA
3
Department of Neurobiology and Behavior, Stony Brook University, Stony Brook, NY 11794, USA
*
Author to whom correspondence should be addressed.
Entropy 2020, 22(5), 537; https://doi.org/10.3390/e22050537
Submission received: 29 March 2020 / Revised: 4 May 2020 / Accepted: 9 May 2020 / Published: 11 May 2020

Abstract

:
Brain dynamics can exhibit narrow-band nonlinear oscillations and multistability. For a subset of disorders of consciousness and motor control, we hypothesized that some symptoms originate from the inability to spontaneously transition from one attractor to another. Using external perturbations, such as electrical pulses delivered by deep brain stimulation devices, it may be possible to induce such transition out of the pathological attractors. However, the induction of transition may be non-trivial, rendering the current open-loop stimulation strategies insufficient. In order to develop next-generation neural stimulators that can intelligently learn to induce attractor transitions, we require a platform to test the efficacy of such systems. To this end, we designed an analog circuit as a model for the multistable brain dynamics. The circuit spontaneously oscillates stably on two periods as an instantiation of a 3-dimensional continuous-time gated recurrent neural network. To discourage simple perturbation strategies, such as constant or random stimulation patterns from easily inducing transition between the stable limit cycles, we designed a state-dependent nonlinear circuit interface for external perturbation. We demonstrate the existence of nontrivial solutions to the transition problem in our circuit implementation.

1. Introduction

Multistability is a widespread phenomenon in the field of dynamical systems where a system exhibits multiple stable states or more generally attractors [1]. Appearing in nearly all disciplines of natural science and engineering, including biology [2,3], chemistry [4,5], electronics [6,7], fluid mechanics [8,9], genetics [10,11], and physics [12,13], researchers have shown substantial interest in such behavior [1]. In neuroscience, computations are thought to be implemented by multistable dynamical systems, and recent experimental and methodological advances have generated renewed interests [14,15,16,17,18,19]. The multiple attractors within these dynamical systems seem to underlie a wide array of functions, including sensory perception [20], motor function [21], and cognition [22,23], as well as dysfunctions, such as movement disorders [24], epilepsy [25], and disorders of consciousness [26,27]. We hypothesized that multistability underlie some dynamical neurological diseases such that manifested symptoms are fundamentally due to the inability to naturally transition from one basin of attraction to another. Under this hypothesis, neurostimulation techniques provide a means to perturb neural systems to assist transition between attractors as a treatment option.
Open-loop electrical stimulation therapies have shown remarkable successes, most notably with Parkinson’s disease [28]. However, open-loop strategies are likely to be insufficient for the general induction of attractor transitions that manifest complex nonlinear dynamics and non-trivial stimulus induced perturbations. For example, high-amplitude low-frequency signals, such as those dominant in disorders of consciousness [27], suggest the existence of strong attractor dynamics which may require a sophisticated feedback control system to transition out [29]. This sets the stage for the next-generation closed-loop neural stimulators that can intelligently learn to induce attractor transitions. The added complexity of the closed-loop stimulation systems calls for a platform to develop and test their efficacy.
In this paper, we aimed to develop a hardware platform by which these intelligent stimulation algorithms could be tested and validated on. Clinically, this platform can serve as an initial stepping stone, prior to the use of animal models and human patients, in demonstrating the successful use of a stimulation regime under our hypothesis for a given application. As is a common method for realizing a dynamical system physically, an analog electronic circuit exhibiting the desired dynamics is constructed [30,31,32]. Due to the nature of widely used electrical stimulation for neurological implants, this medium will serve well as a testbed. While such a system can be developed as a software package, neurostimulators have an analog interface, where typically current is injected into the neural tissue. By keeping our platform entirely analog, we maintain this interface. Furthermore, the analog nature of the brain ensures the existence of inherent noise and imperfections in the system dynamics. We would like to implement these imperfections here as a result of the medium by which the platform is designed on, rather than needing to over-complicate our model. This added variability will require a stimulation algorithm to be more robust.
To reduce unnecessary complexity, we constructed our system to demonstrate the simplest form of oscillatory multistability, birhythmicy. More specifically, the system will simultaneously exhibit two self-exciting limit cycles of notably different frequencies. Given that the state of the system is sufficiently close to one of the two attractors, an intelligent stimulation algorithm can be tested by trying to perturb the system into the other basin of attraction. While this may seem like a substantial simplification of global brain dynamics for typical neurological function, when viewing brain activity at different spatiotemporal scales under specific conditions far fewer attracting states may be present. Two examples are patients under deep anesthesia or within a comatose state. Under such conditions, there appears to be coherent synchronization across large brain regions [26,27]. Such homogeneity would vastly simplify the underlying dynamic behavior, resulting in far fewer attractors. Under the right lens, our developed model becomes more immediately applicable.
In the following section, we derive the system from the general continuous-time dynamical system underlying the gated recurrent unit (GRU), a commonly used recurrent neural network architecture [33,34]. In Section 3, we discuss the details of the circuit design and present the results of the physical realization. In Section 4, we discuss the addition and design of a nonlinear circuit, state dependent on the system described in Section 2 and Section 3, by which external stimulation is interfaced. The addition of this nonlinear stimulator circuit will ensure random or periodic stimulation patterns will be ineffective in inducing transitions between the two attractor states.

2. Birhythmic Dynamics in 3-Dimensions

Our goal was to find a simple bistable dynamical system where each attractor corresponds to a periodic orbit. We draw from the recurrent neural network literature on simple forms of stable limit cycles. Specifically, we utilize the autonomous continuous-time gated recurrent unit (ct-GRU) [33,34] formulation, which can be represented as follows:
h ˙ = ( 1 s ( t ) ) ( T ( t ) h ( t ) ) ( hidden state ) ,
s ( t ) = σ ( U s h ( t ) + b s ) ( update gate ) ,
r ( t ) = σ ( U r h ( t ) + b r ) ( reset gate ) ,
T ( t ) = tanh ( U h ( r ( t ) h ( t ) ) + b h ) ,
where h ( t ) R d is the state of the system, U s , U r , U h R d × d are the parameter matrices, b s , b r , b h R d are the bias vectors, ⊙ represents the Hadamard product, and σ ( z ) = 1 / ( 1 + e z ) is the element-wise logistic sigmoid function. For a given set of parameters, fixed points of the system exist where h ˙ = 0 . Since 1 s ( t ) > 0 , s , this term does not influence the roots of the right-hand side of Equation (1). Therefore, s ( t ) can only affect the speed of the flow and, in turn, can be neglected when choosing a set of parameters for the system to enact a desired structure of attractors [34]. Note that, if the parameters of r ( t ) have been set to zero, the ct-GRU architecture simplifies to the classic ct-tanh-RNN if the parameters of s ( t ) are also set to zero.
In previous work [34], we have shown that, for d = 2 , the ct-GRU is capable of expressing a single limit cycle (attracting closed orbit) in phase space under the following set of parameters:
U r , b r , b h = 0 , U h = 3 cos α sin α sin α cos α ,
where α S 0 and S 0 ( π 21 , π 3.8 ) . The phase portrait depicting this behavior for α = π 5 can be seen in Figure 1, where h x y T .
Extending the system to 3-dimensions allows for the simultaneous existence of two limit cycles in phase space under a single set of parameters. More specifically, the addition of a third dimension enables us to mirror any attractor structure representable for d = 2 across an unstable manifold on the plane defined by the original two dimensions in R 3 . This behavior is depicted in Figure 2A, where now h x y z T , and the parameters are set as follows:
U r , b r , b h = 0 , U h = 3 cos π 5 sin π 5 0 sin π 5 cos π 5 0 0 0 1 .
As stated before, s ( t ) only acts to adjust the speed of phase flow. If U s , b s = 0 , the periods of both limit cycles are equal. As a means to easily decouple the two frequencies of oscillation, the velocity of flow may be made dependent on its vertical position with respect to the z-axis. While the range of the logistic sigmoid function has the benefit of always being defined on ( 0 , 1 ) , it may produce inaccurate results when physically realized along its tails. Furthermore, any function that is strictly positive and sufficiently well-behaved on the region of phase space we are interested in will work under this context.
For simplicity, we redefine s ( t ) linearly as s ( t ) = U s h ( t ) + b s . We then note that Equation () is asymptotically bound to [ 1 , 1 ] d . To account for potential error in the electronic realization, we set s ( t ) such that its output remains strictly positive on ( 1.5 , 1.5 ) × ( 1.5 , 1.5 ) × R . The results of this linear update-gate are demonstrated in Figure 2B,C and achieved under the following set of parameters in conjunction with those of Equation (6):
b s = 0.5 0.5 0.5 , U s = 0 0 1 0 0 1 0 0 0 .
To better grasp the dynamical system depicted in Figure 2 to be later realized, we can rewrite Equations (1)–(4) explicitly in terms of x , y and z with our chosen parameters from Equations (6) and (7) as follows:
x ˙ = z 3 2 x tanh x · 3 2 cos π 5 y · 3 2 sin π 5 ,
y ˙ = z 3 2 y tanh x · 3 2 sin π 5 + y · 3 2 cos π 5 ,
γ z ˙ = 1 2 z tanh 3 2 z ,
where γ R is an added time constant that will be implemented in the circuit realization to adjust the difficulty of transitioning between attracting states. For our implementation of Equations (8)–(10), we let γ = 10 6 .

3. Electronic Physical Realization

Within most applications, smooth continuous-time systems can be realized as electronic circuits comprised of inexpensive components and integrated circuits [30]. In this section, we introduce a comprehensive circuit design to realize Equations (8)–(10) and construct the system on a breadboard. Experimental recordings of trajectories of interest are then compared with the theoretical system derived in Section 2 as a means to validate the realization. All basic operational amplifiers used are TL082CP, and all individual transistors are MPS2222. In addition, two analog multiplier chips are used, which are the standard AD633 four quadrant multipliers. Note that all schematics shown assume unity gain associated with each multiplier. A complete list of all component values in the following schematics can be found in Appendix A.

3.1. Nonlinear Activation Function Circuit

To properly realize Equations (8)–(10) as an analog circuit, we first must account for the nonlinearity in the system; the hyperbolic tangent function. Previous work has allowed us to easily realize this nonlinearity by means of a simple op-amp and transistor circuit as depicted in Figure 3 [30]. Allowing V i n and V o to represent the input and output voltages of the circuit, respectively, Duan and Liao [35] showed that the input-output relation takes the following form:
V 0 = tanh R 2 2 R V T V i n ,
where V T 26 mV is the thermal voltage of the transistors at room temperature. Allowing R 2 = 520   Ω , R 3 = R 4 = 1 k Ω , R 11 = 11 k Ω , all other resistors set to R = 10 k Ω , V C C = 15 V, and V E E = 15 V, reduces the coefficient in front of V i n in Equation (11) to 1, thereby successfully implementing the hyperbolic tangent function as an analog circuit. Further information regarding the error associated with our constructed hyperbolic tangent units can be found in Appendix B.

3.2. Schematics of Electronic Birhythmic RNN

Analog circuits can successfully perform addition and subtraction with just operational amplifiers and appropriately connected resistors. The additional use of capacitors allows for an analog implementation of integration [36]. In regards to our realization of Equations (8)–(10), the analog implementation of the hyperbolic tangent function can be achieved with the schematic demonstrated in Figure 3. In the following schematics, such nonlinear operations are represented by boxes labeled “-tanh”, where terminals IOP1 and IOP2 are the input and output voltage to each hyperbolic tangent unit, respectively. Furthermore, two analog multiplication chips can be used to introduce the phase flow speed dependence on z ( t ) in Equations (8) and (9). Thus, the entire system can be constructed entirely from simple analog components.
We begin with the circuit realization of Equation (10), where the schematic is shown in Figure 4. Note that Equation (10) only comprises of the state variable z ( t ) , and can therefore be build independently of Equations (8) and (9). The variable z ( t ) is represented as the voltage across capacitor C1 in the provided circuit diagram. The input reading stim_out represents the output of the nonlinear stimulator circuit and will be discussed in Section 4. The schematic of the electronic realization of Equations (8) and (9) is shown in Figure 5, where the state variables, x ( t ) and y ( t ) , are represented by the voltages across capacitors C1 and C2, respectively. The “z” input to analog multipliers, M1 and M2, is taken from the integrator output labeled “z” in Figure 4.

3.3. Circuit Construction and Experimental Results

A circuit following the schematics depicted in Figure 4 and Figure 5 was constructed on a breadboard. This system was photographed and is depicted in Figure 6. The blue boxes indicate the three separate hyperbolic tangent units. The magenta box indicates the realization of Equation (10), as described by Figure 4, and the green box highlights the two analog multiplier chips and their configuration.
The circuit was constructed so that each limit cycle maintains an approximately constant z-component value when sufficiently close to one of the two attractors. Using an four channel oscilloscope, we observed the behavior of the system for different initial conditions as a means to demonstrate stability across runs. These included ten trials initialized randomly within a subset of each of the two basins of attractions. These subsets are ( 1.5 , 1.5 ) × ( 1.5 , 1.5 ) × ( 0 , 1.5 ) and ( 1.5 , 1.5 ) × ( 1.5 , 1.5 ) × ( 1.5 , 0 ) for the slow and fast limit cycles, respectively. Each trajectory was given five seconds to relax, and then the state of the system was recorded for five seconds. For the slow limit cycle, the sample mean of the average z-component values across all recorded time steps was calculated to be 0.858 V across trials, with a standard deviation of 5.800 × 10 3 V across those trial averages. For any single trial, the expected standard deviation of the z-component values across all time steps was calculated to be 1.110 × 10 16 V, suggesting the upper limit cycle is extremely stable. Similarly, for the fast limit cycle, the sample mean of the average z-component values across all recorded time steps was calculated to be 0.907 V across trials, with a standard deviation of 9.248 × 10 3 V across those trial averages. For any single trial, the expected standard deviation of the z-component values across all time steps was calculated to be 2.133 × 10 3 V, suggesting the lower limit cycle is also highly stable.
An example of these recordings are depicted in Figure 7 for two different initializations; one within the basin of attraction for each of the two stable limit cycles. Figure 7A,C show the asymptotic behavior of the three state variable voltages in time [ x ( t ) , y ( t ) , z ( t ) ] . As intended, the period of the each limit cycles are visibly different, and qualitatively match the asymptotic behavior demonstrated in Figure 2B,C. Furthermore, Figure 7B,D depict the trajectories projected onto the x-y plane, shown in Figure 7A,C, respectively. These recordings indicate the same asymptotic behavior as demonstrated in Figure 1. As such, we can conclude that the analog implementation of Equations (8)–(10) was successfully realized from the model derived in Section 2.

4. Nonlinear Stimulator Circuit Design and Discussion

When looking at the system derived in Section 2, as expressed by Equations (8)–(10), we notice the geometric simplicity of the global dynamics. As shown in Figure 2A, the two periodic attractors are mirror images of one another across a planar unstable manifold on the x-y plane. This symmetry enabled us to easily decouple the frequencies of each limit cycle by introducing a strictly positive variable time constant dependent only on our z-coordinate. However, this introduces a clear problem from the point of view of developing a neurostimulation testbed. In order to transition between basins of attraction, we only need to worry about the component z ( t ) . In other words, transition between attractors can be achieved with constant stimulation on z. This solution to attractor transition is trivial, and does not require the use of an intelligent algorithm to solve, rendering it inadequate as a testbed. In order to negate this issue, we develop a state-dependent nonlinear stimulator interface circuit by which all stimulus must pass through. In addition, we allow only one location of stimulation within the circuit previously designed, as marked by the node labeled stim_out in Figure 4, where the output to the stimulator circuit will be fed in and summed with the current value of z ˙ in the system. By extending the system properly in this way, we can prevent straightforward stimulation patterns (i.e., constant, random, periodic, etc.) from inducing attractor transitions.
Ideally, we want to develop a stimulation circuit such that a pulse stimulated at a random time may cause the output, stim_out, to either increase or decrease. Furthermore, the amount by which the signal can change should exist on a continuum, rather than outputting a voltage from a finite set of values. The final requirement we will enforce for such a circuit will be that a stimulation pulse delivered at a random time should have equal probability of increasing or decreasing the output signal. This last requirement ensures that if one stimulates randomly or continuously, the expected value of stim_out averaged across all time will be zero as time approaches infinity.
Due to the sinusoidal nature of the x and y components in Equations (8) and (9), Equation (12) will serve as the input/output relation of the stimulator circuit at a given time t.
S o u t ( t ) = S i n ( t ) x ( t ) y ( t ) ,
where S o u t is the output voltage of the stimulator circuit, and S i n is the input stimulus (note that S i n = 0 when no stimulus is applied). We note that this system satisfies all of our requirements, as z is independent of x and y. To transition from the slow limit cycle to the fast limit cycle, stimulation must be applied primarily when x ( t ) and y ( t ) are of opposite signs. By doing this, the rate of change in z ( t ) can be made negative on a subset of these intervals (dependent on the current value of z ( t ) ). Assuming that this downward forcing of z ( t ) can both overcome the pullback from the attractor during these intervals and travel further downward than the accumulated upward travel during periods of no stimulation, then, after a sufficient amount of time, the system will jump over the unstable manifold at z = 0 and into the basin of attraction of the fast limit cycle. Since the speed of oscillation is dependent on z ( t ) , any added stimulation will change the time window by which proper stimulation should be applied, ensuring that a periodic stimulation regime will fail to transition between the attractors. By the same argument, to transition from the fast limit cycle to the slow limit cycle, stimulation must be applied primarily when x ( t ) and y ( t ) are of the same sign. Figure 8 depicts the schematic for the described stimulator circuit, where S i n is labeled as Stimulus.
A physical realization of the stimulator circuit was created in conjunction with the birhythmic system and tested with a stimulus that could be turned on (4 V) or off (0 V). While this voltage range is certainly larger than those seen during neurostimulation, we note that this testbed acts to validate the underlying logic of a stimulation algorithm, and that this range is made proportional, in a practical sense, to the voltage range of the constructed dynamics. For interfacing with a traditional neural stimulation device, an additional signal amplifying circuit will need to be constructed to take in stimulus (typically low current biphasic pulse injections), and control a properly scaled voltage output which will be fed into the input of our designed nonlinear stimulator circuit.
In order to test the validity of our constructed system we recorded twenty trials for each of several stimulation methods, where half of the trials had the state of the system initialized on the slow limit cycle, and the other half of the trials were initialized on the fast limit cycle. Each trial was recorded for five seconds. We show that neither 4 V constant stimulation nor 4 V manually random stimulation induce state transition, and that these results appear invariant to the initial phase angle of the oscillations in x ( t ) and y ( t ) when stimulation is applied. Furthermore, we implement a stimulation pattern that can successfully transition between the stable limit cycles.
For constant stimulation, we statistically quantify the expectation of the maximum euclidean distance z ( t ) moves away from the limit cycle that the system is initialized on for a given trial. The z coordinates of both limit cycles empirically derived in Section 3.3 will be used. In the case of the slow limit cycle, the sample mean of the maximum distances achieved on each trial is 0.230 V, with a standard deviation of 0.026 V across trials, implying that each trial was nearly identical in its effectiveness to transition between states. For the fast limit cycle, this sample mean is 0.147 V, and no variance was detectable in these recordings in 24-bit resolution. A 1.5 s segment of z ( t ) for all trials is depicted in Figure 9, where constant stimulation is denoted in turquoise. Figure 10A,B depict an example of the resultant behavior of the system under constant 4 V stimulation initialized on the slow and fast limit cycles, respectively. The yellow, blue, and pink curves represent the x, y, and z components of the system, and the green curve depicts the voltage over time of stim_out.
For random stimulation, we do the same thing. In the case of the slow limit cycle, the sample mean of the maximum distance achieved on each trial is 0.242 V, with a standard deviation of 0.022 V across trials. For the fast limit cycle, this sample mean is 0.177 V, with standard deviation of 0.013 V across trials. A segment of z ( t ) for all of these trials is depicted in orange in Figure 9. Figure 11A,B depict an example of the resultant behavior for trajectories initialized on the slow and fast limit cycles, respectively.
Despite executing stimulation at a random initial phase angle in terms of the oscillation in x and y, the sample statistics suggest a qualitative homogeneity across trials for the two stimulation patterns, which can be seen in Figure 9. In addition, as expected, both of these stimulation regimes do not escape the basin of attraction by which the system is initialized in, as to do so would require the stimulation when the system state is in the proper two quadrants of the x-y plane to overtake the stimulation in the other two for a prolonged period of time. By restricting our stimulation to the desired intervals, we can achieve state transition. In this case, we simply place a 1N4148 diode in the appropriate orientation just prior to the stim_out node in Figure 8 and apply constant stimulation. All of the trials successfully traveled over the separating unstable manifold at z = 0 as depicted by the red and blue curves in Figure 9, representing the trials initialized on the slow and fast limit cycles, respectively. An example of the behavior for all system coordinates and the stimulator circuit output is shown in Figure 12 for both trajectories initialized on the slow and fast limit cycles.
Such a demonstration indicates the existence of a stimulation pattern capable of transitioning between basins of attraction in the system. An input stimulus pattern of this form will have to be mimicked by an intelligent algorithm to a sufficient degree of accuracy without the added rectifier, thus aiding in the validation of the use of that algorithm. However, note that a stimulation system would not have direct access to the state variables, the dynamical system, nor how the stimulus modulates the states. Any such algorithm will need to discover the latent dynamics and learn to control the states from observations [32] at the same time it learns to transition out of the current basin of attraction. The continuously changing intervals of when stimulation should be applied, while geometrically simple, are a highly nontrivial relation to uncover within a dynamical system in this manner and should prove challenging for a general algorithm required to learn it.

5. Conclusions

For this paper, an electronic testbed for intelligent neurostimulation methods was developed from a physical realization of the dynamical system underlying the architecture of an artificial gated recurrent neural network. Using simple analog components, the system is fabricated such that it exhibits birhythmic behavior, but the stimulation pattern required to transition between attracting states is made nontrivial. As such, standard open-loop stimulation regimes will be unable to induce attractor transitions. We hypothesized that an inability to perform analogous state transitions may underlie some neurological diseases, and prior evidence suggests that, under this hypothesis, global cortical brain function may give rise to very few attractors under the proper spatiotemporal scale of viewing for comatose patients and patients under deep anesthesia. If correct, this enables our system to aid in validating the efficacy of next-generation neurostimulation algorithms upon successfully jumping between basins of attraction. We suspect multiple realizations of this system can be appropriately coupled together to form a more complicated testbed with a wider array of dynamics. Such an extension may be more immediately applicable to more complex neurological functions and dysfunctions. We leave this to future work.

Author Contributions

I.D.J. and I.M.P. conceptualized the project. I.D.J. conducted experiments. I.D.J. and I.M.P. wrote the manuscript. All authors have read and agreed to the published version of the manuscript.

Funding

This work is supported by Stony Brook University’s Discovery Award, NIH EB-026946, and NSF IIS-1845836. I.D.J. was supported partially by the Institute for Advanced Compuational Science Jr. Researcher Fellowship, at Stony Brook University.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A. Circuit Component Values

Appendix A.1. Component Values for Figure 4

V C C = 15 V, V E E = 15 V
R 1 = R 3 = R 4 = R 5 = R 6 = R 12 = R 13 = 100 k Ω , R 2 = 1 M Ω , R 7 = 50 k Ω , R 8 = R 9 = R 11 = 10 k Ω , R 10 = 5 k Ω
C 1 = 1 μ F

Appendix A.2. Component Values for Figure 5

V C C = 15 V, V E E = 15 V, V D D = 1 V
R 1 = R 3 = R 5 = R 6 = R 7 = R 9 = R 10 = R 11 = R 12 = R 14 = R 16 = R 17 = R 19 = R 21 = R 22 = R 23 = R 24 = R 25 = R 26 = R 27 = R 28 = R 29 = R 31 = R 32 = R 33 = R 34 = R 35 = R 36 = R 37 = R 38 = R 39 = R 40 = R 42 = R 43 = R 44 = R 45 = R 46 = R 47 = R 48 = R 49 = 10 k Ω , R 2 = R 4 = R 13 = R 15 = 100 k Ω , R 8 = R 20 = 8820 Ω , R 18 = R 50 = 2315 Ω , R 30 = R 41 = 5 k Ω .
C 1 = C 2 = 10 μ F , C 3 = C 4 = C 5 = C 6 = 1 μ F

Appendix A.3. Component Values for Figure 8

V C C = 15 V, V E E = 15 V
R 1 = R 2 = R 3 = 100 k Ω
C 1 = C 2 = C 3 = C 4 = C 5 = 1 μ F

Appendix B. Hyperbolic Tangent Implementation

Our system depicted in Equations (8)–(10) requires the construction of three separate hyperbolic tangent units. Following construction, 21 input voltages, equally spaced on [−2, 2] V, were applied to each hyperbolic tangent unit as input, and each respective output voltage was recorded and compared with the hyperbolic tangent function for that given input. These recordings and comparisons of accuracy of the three hyperbolic tangent units are displayed in Figure A1. The average mean squared error of the built hyperbolic tangent units with respect to the hyperbolic tangent function was approximately 86.5 mV.
Figure A1. Input/output relation of three physically realized hyperbolic tangent circuits, interpolated through 21 points, compared with the analytic hyperbolic tangent function.
Figure A1. Input/output relation of three physically realized hyperbolic tangent circuits, interpolated through 21 points, compared with the analytic hyperbolic tangent function.
Entropy 22 00537 g0a1

References

  1. Pisarchik, A.N.; Feudel, U. Control of multistability. Phys. Rep. 2014, 540, 167–218. [Google Scholar] [CrossRef]
  2. Abou-Jaoudé, W.; Chaves, M.; Gouzé, J.L. Mechanisms for coexistence of two limit cycles in a biochemical model1. IFAC Proc. Vol. 2011, 44, 8378–8383. [Google Scholar] [CrossRef]
  3. Decroly, O.; Goldbeter, A. Birhythmicity, chaos, and other patterns of temporal self-organization in a multiply regulated biochemical system. Proc. Natl. Acad. Sci. USA 1982, 79, 6917–6921. [Google Scholar] [CrossRef] [Green Version]
  4. Hudson, J.L.; Mankin, J.C. Chaos in the Belousov–Zhabotinskii reaction. J. Chem. Phys. 1981, 74, 6171–6177. [Google Scholar] [CrossRef]
  5. Ngonghala, C.N.; Feudel, U.; Showalter, K. Extreme multistability in a chemical model system. Phys. Rev. E 2011, 83, 056206. [Google Scholar] [CrossRef] [Green Version]
  6. Biswas, D.; Banerjee, T.; Kurths, J. Control of birhythmicity through conjugate self-feedback: Theory and experiment. Phys. Rev. E 2016, 94, 042226. [Google Scholar] [CrossRef] [Green Version]
  7. In, V.; Kho, A.; Neff, J.D.; Palacios, A.; Longhini, P.; Meadows, B.K. Experimental observation of multifrequency patterns in arrays of coupled nonlinear oscillators. Phys. Rev. Lett. 2003, 91, 244101. [Google Scholar] [CrossRef] [Green Version]
  8. Ravelet, F.; Marié, L.; Chiffaudel, A.; Daviaud, F. Multistability and Memory Effect in a Highly Turbulent Flow: Experimental Evidence for a Global Bifurcation. Phys. Rev. Lett. 2004, 93, 164501. [Google Scholar] [CrossRef] [Green Version]
  9. Shiau, Y.H.; Peng, Y.F.; Hwang, R.R.; Hu, C.K. Multistability and symmetry breaking in the two-dimensional flow around a square cylinder. Phys. Rev. E 1999, 60, 6188–6191. [Google Scholar] [CrossRef] [Green Version]
  10. Koseska, A.; Volkov, E.; Zaikin, A.; Kurths, J. Inherent multistability in arrays of autoinducer coupled genetic oscillators. Phys. Rev. E 2007, 75, 031916. [Google Scholar] [CrossRef] [Green Version]
  11. Labavić, D.; Meyer-Ortmanns, H. Networks of coupled circuits: From a versatile toggle switch to collective coherent behavior. Chaos Interdiscip. J. Nonlinear Sci. 2014, 24, 043118. [Google Scholar] [CrossRef] [Green Version]
  12. Ding, W.C.; Xie, J.H.; Sun, Q.G. Interaction of Hopf and period doubling bifurcations of a vibro-impact system. J. Sound Vib. 2004, 275, 27–45. [Google Scholar] [CrossRef] [Green Version]
  13. Kiss, I.Z.; Hudson, J.L. Chaotic cluster itinerancy and hierarchical cluster trees in electrochemical experiments. Chaos Interdiscip. J. Nonlinear Sci. 2003, 13, 999–1009. [Google Scholar] [CrossRef]
  14. Little, W.A. The existence of persistent states in the brain. Math. Biosci. 1974, 19, 101–120. [Google Scholar] [CrossRef]
  15. Breakspear, M. Dynamic models of large-scale brain activity. Nat. Neurosci. 2017, 20, 340–352. [Google Scholar] [CrossRef]
  16. Fairhall, A.; Machens, C. Editorial overview: Computational neuroscience. Curr. Opin. Neurobiol. 2017, 46, A1–A5. [Google Scholar] [CrossRef]
  17. Sussillo, D. Neural circuits as computational dynamical systems. Curr. Opin. Neurobiol. 2014, 25, 156–163. [Google Scholar] [CrossRef]
  18. La Camera, G.; Fontanini, A.; Mazzucato, L. Cortical computations via metastable activity. Curr. Opin. Neurobiol. 2019, 58, 37–45. [Google Scholar] [CrossRef]
  19. Zhao, Y.; Park, I.M. Interpretable Nonlinear Dynamic Modeling of Neural Trajectories. arXiv 2016, arXiv:q-bio.QM/1608.06546. [Google Scholar]
  20. Li, H.H.; Rankin, J.; Rinzel, J.; Carrasco, M.; Heeger, D.J. Attention model of binocular rivalry. Proc. Natl. Acad. Sci. USA 2017, 114, E6192–E6201. [Google Scholar] [CrossRef] [Green Version]
  21. Churchland, M.M.; Cunningham, J.P.; Kaufman, M.T.; Foster, J.D.; Nuyujukian, P.; Ryu, S.I.; Shenoy, K.V. Neural population dynamics during reaching. Nature 2012, 487, 51–56. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  22. Mante, V.; Sussillo, D.; Shenoy, K.V.; Newsome, W.T. Context-dependent computation by recurrent dynamics in prefrontal cortex. Nature 2013, 503, 78–84. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  23. Jazayeri, M.; Afraz, A. Navigating the Neural Space in Search of the Neural Code. Neuron 2017, 93, 1003–1014. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  24. Wagenaar, R.C.; van Emmerik, R.E.A. Dynamics of movement disorders. Hum. Mov. Sci. 1996, 15, 161–175. [Google Scholar] [CrossRef]
  25. Iasemidis, L.D. Epileptic seizure prediction and control. IEEE Trans. Biomed. Eng. 2003, 50, 549–558. [Google Scholar] [CrossRef]
  26. Lewis, L.D.; Piantoni, G.; Peterfreund, R.A.; Eskandar, E.N.; Harrell, P.G.; Akeju, O.; Aglio, L.S.; Cash, S.S.; Brown, E.N.; Mukamel, E.A.; et al. A transient cortical state with sleep-like sensory responses precedes emergence from general anesthesia in humans. Elife 2018, 7, e33250. [Google Scholar] [CrossRef]
  27. Schiff, N.D.; Nauvel, T.; Victor, J.D. Large-scale brain dynamics in disorders of consciousness. Curr. Opin. Neurobiol. 2014, 25, 7–14. [Google Scholar] [CrossRef] [Green Version]
  28. Little, S.; Brown, P. Focusing Brain Therapeutic Interventions in Space and Time for Parkinson’s Disease. Curr. Biol. 2014, 24, R898–R909. [Google Scholar] [CrossRef] [Green Version]
  29. Hocker, D.; Park, I.M. Myopic control of neural dynamics. PLoS Comput. Biol. 2019, 15, e1006854. [Google Scholar] [CrossRef] [Green Version]
  30. Bao, B.; Qian, H.; Wang, J.; Xu, Q.; Chen, M.; Wu, H.; Yu, Y. Numerical analyses and experimental validations of coexisting multiple attractors in Hopfield neural network. Nonlinear Dyn. 2017, 90, 2359–2369. [Google Scholar] [CrossRef]
  31. Rahman, A.; Jordan, I.; Blackmore, D. Qualitative models and experimental investigation of chaotic NOR gates and set/reset flip-flops. Proc. R. Soc. Math. Phys. Eng. Sci. 2018, 474, 20170111. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  32. Zhao, Y.; Nassar, J.; Jordan, I.; Bugallo, M.; Park, I.M. Streaming Variational Monte Carlo. arXiv 2019, arXiv:1906.01549. [Google Scholar]
  33. Cho, K.; van Merrienboer, B.; Gulcehre, C.; Bahdanau, D.; Bougares, F.; Schwenk, H.; Bengio, Y. Learning Phrase Representations using RNN Encoder-Decoder for Statistical Machine Translation. arXiv 2014, arXiv:1406.1078. [Google Scholar]
  34. Jordan, I.D.; Sokol, P.A.; Park, I.M. Gated recurrent units viewed through the lens of continuous time dynamical systems. arXiv 2019, arXiv:1906.01005. [Google Scholar]
  35. Duan, S.; Liao, X. An electronic implementation for Liao’s chaotic delayed neuron model with non-monotonous activation function. Phys. Lett. A 2007, 369, 37–43. [Google Scholar] [CrossRef]
  36. Wilmshurst, T.H. Analog Circuit Techniques: With Digital Interfacing; Elsevier: Amsterdam, The Netherlands, 2001. [Google Scholar]
Figure 1. Planar Limit Cycle with 2D continuous-time gated recurrent unit (ct-GRU) depicted in phase space: The red dot indicates an unstable fixed point at the origin unstable, while orange and pink lines represent the x and y nullclines, respectively. Purple lines indicate various trajectories of the hidden state. Direction of the flow is determined by the black arrows, where the colormap underlying the figure depicts the magnitude of the velocity of the flow in log scale.
Figure 1. Planar Limit Cycle with 2D continuous-time gated recurrent unit (ct-GRU) depicted in phase space: The red dot indicates an unstable fixed point at the origin unstable, while orange and pink lines represent the x and y nullclines, respectively. Purple lines indicate various trajectories of the hidden state. Direction of the flow is determined by the black arrows, where the colormap underlying the figure depicts the magnitude of the velocity of the flow in log scale.
Entropy 22 00537 g001
Figure 2. Birhythmicy in 3-dimensions: (A): light blue manifold on the x y plane separates the basins of attraction of the upper and lower limit cycles. Trajectories are colored either dark blue or purple, depending on which basin of attraction they are initialized in. Red dots indicate fixed points, and black arrows depict the direction of flow. (B,C): x, y, and z components of trajectories initialized in the basins of attractions for the top and bottom limit cycles, respectively. Solid colored lines indicate x ( t ) , dashed lines indicate y ( t ) , and black lines indicate z ( t ) .
Figure 2. Birhythmicy in 3-dimensions: (A): light blue manifold on the x y plane separates the basins of attraction of the upper and lower limit cycles. Trajectories are colored either dark blue or purple, depending on which basin of attraction they are initialized in. Red dots indicate fixed points, and black arrows depict the direction of flow. (B,C): x, y, and z components of trajectories initialized in the basins of attractions for the top and bottom limit cycles, respectively. Solid colored lines indicate x ( t ) , dashed lines indicate y ( t ) , and black lines indicate z ( t ) .
Entropy 22 00537 g002
Figure 3. Electronic circuit realization of the hyperbolic tangent function, as implemented in Reference [30]. V i n and V o represent the input and output signals, respectively.
Figure 3. Electronic circuit realization of the hyperbolic tangent function, as implemented in Reference [30]. V i n and V o represent the input and output signals, respectively.
Entropy 22 00537 g003
Figure 4. Circuit schematic of z ˙ for the birhythmic system. The system block labeled -tanh represents the circuit depicted in Figure 3, where I01 and I02 correspond to V i n and V o , respectively. The terminal labeled stim_out represents the output to the stimulator circuit, as discussed in Section 4.
Figure 4. Circuit schematic of z ˙ for the birhythmic system. The system block labeled -tanh represents the circuit depicted in Figure 3, where I01 and I02 correspond to V i n and V o , respectively. The terminal labeled stim_out represents the output to the stimulator circuit, as discussed in Section 4.
Entropy 22 00537 g004
Figure 5. Circuit schematic of x ˙ and y ˙ for the birhythmic system. The system blocks labeled -tanh1 and -tanh2 represent the circuit depicted in Figure 3, where I01 and I02 correspond to V i n and V o , respectively, for both blocks. The two multiplier chips, M1 and M2, are assumed to operate with unity gain.
Figure 5. Circuit schematic of x ˙ and y ˙ for the birhythmic system. The system blocks labeled -tanh1 and -tanh2 represent the circuit depicted in Figure 3, where I01 and I02 correspond to V i n and V o , respectively, for both blocks. The two multiplier chips, M1 and M2, are assumed to operate with unity gain.
Entropy 22 00537 g005
Figure 6. Physical birhythmic circuit constructed on a breadboard. Blue boxes represent hyperbolic tangent units. The magenta box indicates the subsection of the circuit generating z ˙ , and the green box indicates the analog multipliers.
Figure 6. Physical birhythmic circuit constructed on a breadboard. Blue boxes represent hyperbolic tangent units. The magenta box indicates the subsection of the circuit generating z ˙ , and the green box indicates the analog multipliers.
Entropy 22 00537 g006
Figure 7. Experimental recordings of the birhythmic circuit: (A,C): x (yellow), y (blue), and z (pink) with respect to time of trajectories within the basin of attraction of the fast and slow limit cycles, respectively. (B,D): Projection of the corresponding trajectories in (A,C) onto the x-y plane, respectively.
Figure 7. Experimental recordings of the birhythmic circuit: (A,C): x (yellow), y (blue), and z (pink) with respect to time of trajectories within the basin of attraction of the fast and slow limit cycles, respectively. (B,D): Projection of the corresponding trajectories in (A,C) onto the x-y plane, respectively.
Entropy 22 00537 g007
Figure 8. Schematic for nonlinear stimulator circuit, with input labeled as Stimulus. The output, labeled stim_out, is fed into the terminal with the same name presented in the circuit diagram shown in Figure 4. The two multiplier chips, M1 and M2, are assumed to operate with unity gain, and the x and y terminals are fed into the equivalently named terminals depicted in Figure 5.
Figure 8. Schematic for nonlinear stimulator circuit, with input labeled as Stimulus. The output, labeled stim_out, is fed into the terminal with the same name presented in the circuit diagram shown in Figure 4. The two multiplier chips, M1 and M2, are assumed to operate with unity gain, and the x and y terminals are fed into the equivalently named terminals depicted in Figure 5.
Entropy 22 00537 g008
Figure 9. A time window of z ( t ) for all experimental trials: Red and blue trajectories demonstrate resultant behavior from stimulation patterns designed to transition states from the slow and fast limit cycles, respectively. Turquoise trajectories depict trials of constant stimulation, and orange trajectories show trials of random stimulation.
Figure 9. A time window of z ( t ) for all experimental trials: Red and blue trajectories demonstrate resultant behavior from stimulation patterns designed to transition states from the slow and fast limit cycles, respectively. Turquoise trajectories depict trials of constant stimulation, and orange trajectories show trials of random stimulation.
Entropy 22 00537 g009
Figure 10. x (yellow), y (blue), z (pink), and stim_out (green) with respect to time of trajectories within the basin of attraction of the slow (A) and fast (B) limit cycles, under constant stimulation. Note that this stimulation regime does not successfully transition between the two attracting states in either direction.
Figure 10. x (yellow), y (blue), z (pink), and stim_out (green) with respect to time of trajectories within the basin of attraction of the slow (A) and fast (B) limit cycles, under constant stimulation. Note that this stimulation regime does not successfully transition between the two attracting states in either direction.
Entropy 22 00537 g010
Figure 11. x (yellow), y (blue), z (pink), and stim_out (green) with respect to time of trajectories within the basin of attraction of the slow (A) and fast (B) limit cycles, with random stimulation. Note that this stimulation regime does not successfully transition between the two attracting states in either direction.
Figure 11. x (yellow), y (blue), z (pink), and stim_out (green) with respect to time of trajectories within the basin of attraction of the slow (A) and fast (B) limit cycles, with random stimulation. Note that this stimulation regime does not successfully transition between the two attracting states in either direction.
Entropy 22 00537 g011
Figure 12. Examples of stimulation patterns capable of inducing transition between states: x (yellow), y (blue), z (pink), and stim_out (green) with respect to time of trajectories initialized within the basin of attraction of the slow (A) and fast (B) limit cycles. As z ( t ) changes with stimulation so does the frequency of oscillation. As such, the time window to stimulate shifts continuously.
Figure 12. Examples of stimulation patterns capable of inducing transition between states: x (yellow), y (blue), z (pink), and stim_out (green) with respect to time of trajectories initialized within the basin of attraction of the slow (A) and fast (B) limit cycles. As z ( t ) changes with stimulation so does the frequency of oscillation. As such, the time window to stimulate shifts continuously.
Entropy 22 00537 g012

Share and Cite

MDPI and ACS Style

Jordan, I.D.; Park, I.M. Birhythmic Analog Circuit Maze: A Nonlinear Neurostimulation Testbed. Entropy 2020, 22, 537. https://doi.org/10.3390/e22050537

AMA Style

Jordan ID, Park IM. Birhythmic Analog Circuit Maze: A Nonlinear Neurostimulation Testbed. Entropy. 2020; 22(5):537. https://doi.org/10.3390/e22050537

Chicago/Turabian Style

Jordan, Ian D., and Il Memming Park. 2020. "Birhythmic Analog Circuit Maze: A Nonlinear Neurostimulation Testbed" Entropy 22, no. 5: 537. https://doi.org/10.3390/e22050537

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop