3rd IFAC/IEEE CSS Workshop on Control of Systems Governed by 3rd IFAC/IEEE CSS Workshop on Control of Systems Governed by 3rd IFAC/IEEE CSS Workshop on Control of Systems Governed by Partial Partial Available online at www.sciencedirect.com Partial Differential Equation, and XI Workshop Control of Distributed 3rd IFAC/IEEE CSS Workshop on Control of Systems Governed by Differential Equation, Differential Equation, and and XI XI Workshop Workshop Control Control of of Distributed Distributed Parameter Systems Partial Parameter Systems Parameter Systems Oaxaca, Mexico, May 20-24, Differential Equation, and XI 2019 Workshop Control of Distributed Oaxaca, May Oaxaca, Mexico, Mexico, May 20-24, 20-24, 2019 2019 Parameter Systems Oaxaca, Mexico, May 20-24, 2019 IFAC PapersOnLine 52-2 (2019) 150–155
ScienceDirect
Control Mechanisms on the Initiation of Control Control Mechanisms Mechanisms on on the the Initiation Initiation of of Cortical Wave Propagation Preceding Control Mechanisms on the Initiation of Cortical Wave Propagation Preceding Cortical Seizure Wave Propagation Preceding Termination Cortical Seizure Wave Propagation Termination Seizure TerminationPreceding Seizure Termination an-Guerra ∗∗∗ Laura R. Gonz´ alez-Ram´ırez ∗∗∗ Rosalba Galv´
Laura a ırez Galv´ ∗ ∗ Rosalba ∗ a Laura R. R. Gonz´ Gonz´ alez-Ram´ lez-Ram´ ırez Rosalba Galv´ an-Guerra n-Guerra ∗∗ Juan E. Vel´ a zquez-Vel´ a zquez ∗ ∗ Juan E. Vel´ a zquez-Vel´ a zquez ∗ ∗ ∗ an-Guerra ∗ Laura R. Gonz´ alez-Ram´ Rosalba Galv´ Juan E. Vel´ aırez zquez-Vel´ azquez ∗ JuanNacional, E. Vel´ azquez-Vel´ azquez ∗ Interdisciplinaria ecnico Unidad Profesional ∗ Instituto Polit´ ∗ Nacional, Unidad Profesional Interdisciplinaria ∗ Instituto Polit´ e cnico ∗ Instituto Polit´ eCampus cnico Nacional, Unidad Profesional Interdisciplinaria de Ingenier´ ıa Hidalgo, Departmento de Formaci´ o n B´ a sica ıa Campus Hidalgo, Departmento de Formaci´ o a ∗ de Ingenier´ Instituto Polit´ cnico Nacional, Unidad Profesional Interdisciplinaria de Ingenier´ ıa eCampus Hidalgo, Departmento de Formaci´ on n B´ B´ asica sica Disciplinaria, M´ e xico (e-mail:
[email protected]) Disciplinaria, eexico
[email protected]) de Ingenier´ ıa CampusM´ Hidalgo, Departmento de Formaci´ on B´ asica Disciplinaria, M´ xico (e-mail: (e-mail:
[email protected]) Disciplinaria, M´exico (e-mail:
[email protected]) Abstract: Organized patterns of activity in the form of traveling waves are characteristic Abstract: Organized patterns of activity in the form traveling waves characteristic Abstract: Organized patterns of activity in thepreceding form of of seizure traveling waves are areUnderstanding characteristic features observed in clinical data recorded in vivo termination. features observed in clinical data recorded in vivo preceding seizure termination. Understanding Abstract: Organized patterns of activity in the form of traveling waves are characteristic features observed in clinical data recorded in vivo preceding seizure termination. Understanding the exact biological mechanisms leading to this type of cortical wave propagation is an open the exact biological mechanisms leading to this of cortical wave propagation is features in clinical datasimple recorded inphenomena vivo type preceding seizure termination. Understanding the exactobserved biological mechanisms leading to this type ofappear cortical wave propagation is an an open open challenge. Since these relatively brain rigth before seizure termination, challenge. Since these relatively simple brain phenomena appear rigth before seizure termination, the exact biological mechanisms leading to phenomena thisthat typeinitiate ofappear cortical wave propagation is Knowing an open challenge. Since these relatively the simple brain rigth before seizure termination, it is interesting to understand mechanisms such wave propagation. it to understand the mechanisms that initiate such wave propagation. Knowing challenge. Since these relatively simple brainhow phenomena appear rigth before seizure and termination, it is is interesting interesting to may understand the mechanisms that initiate such wave propagation. Knowing these mechanisms help to understand to shorten the length of aa seizure position these mechanisms may help to understand how to shorten the length of seizure and position it is interesting to understand the mechanisms that initiate such wave propagation. Knowing these mechanisms may help to understand how to shorten the length of a seizure and position the brain into a state of seizure termination sooner. the brain into a state of seizure termination sooner. these mechanisms mayofhelp to understand how tomechanism shorten the length of a seizure and position the brain into we a state seizure termination sooner. In this paper, present a minimum-time control for the initiation of traveling wave In this this paper, we present minimum-time control mechanism for for the the initiation initiation of of traveling wave wave the brain a state of seizure termination sooner. In paper, present aa minimum-time control mechanism activity ininto an we activity-based neural field model with features observed in clinical traveling recordings. activity in an activity-based neural field model with features observed in clinical recordings. In this paper, present a minimum-time control mechanism the initiation of traveling wave activity in an we activity-based neural field model with features for observed in clinical recordings. © 2019, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. activity in an activity-based neural field model with features observed in clinical recordings. Keywords: Keywords: Neural Neural Field, Field, Traveling Traveling Wave, Wave, Minimum-Time Minimum-Time Wave Wave Initiation, Initiation, Epilepsy Epilepsy Keywords: Neural Field, Traveling Wave, Minimum-Time Wave Initiation, Epilepsy Keywords: Neural Field, Traveling Wave, Minimum-Time Wave Initiation, Epilepsy 1. INTRODUCTION aa sufficiently sufficiently large large and and fixed fixed amount amount of of time. time. This This work work is is 1. isa 1. INTRODUCTION INTRODUCTION adevoted sufficiently large and fixed amount of time. This work devoted to the understanding of how much excitation in to the understanding of how much excitation in 1. INTRODUCTION a sufficiently large and fixed amount of time. This work devoted to the understanding of how much excitation in isa a minimum-time can initiate such activity. Epilepsy is a neurological condition characterized by reEpilepsy is is a neurological neurological condition condition characterized characterized by by rere- minimum-time can activity. devoted to the understanding of how much excitation in a minimum-time can initiate initiate such such activity. Epilepsy current and aunprovoked seizures. Epilepsy has been sucfield models describe mean activities of neuronal current seizures. Epilepsy has minimum-time can initiate such activity. Epilepsyand is aunprovoked neurological condition characterized bysucre- Neural Neural models mean of current and unprovoked seizures. Epilepsy has been been successfully described as a dynamic disease (Milton, 2003). Neural field field composed models describe describe mean activities activities of neuronal neuronal populations by excitatory and inhibitory cessfully described as a dynamic disease (Milton, 2003). current and unprovoked seizures. disease Epilepsy has been suc- populations composed by excitatory and inhibitory neuneucessfully described asmechanisms a dynamic (Milton, 2003). The exact biological that produce aa seizure Neural field models describe mean activities ofinhibitory neuronal populations composed by excitatory and inhibitory neuThe exact biological mechanisms that produce seizure rons. A neuron is considered to be excitatory or cessfully described asmechanisms a dynamic disease (Milton, 2003). A neuron is considered to be excitatory or inhibitory The exact biological that produce a seizure remain mostly unknown. However, there is evidence of rons. populations composed by excitatory and In inhibitory neurons. A neuron is considered toneighbors. be excitatory or inhibitory depending on its effect on its this work remain mostly unknown. However, there is evidence of The exact biological mechanisms produce a seizure remain mostly unknown. However, there of depending on its effect on its neighbors. In this work we abnormal neuronal activity in thethat wake ofisaevidence seizure inwe rons. A neuron is considered beofexcitatory or an inhibitory depending on its effect on itstoneighbors. In this work we focus on a neural field consistent exclusively excitaabnormal neuronal activity in the wake of a seizure inremain hyper-excitation mostly unknown. However, there of focus on a neural field consistent of exclusively an excitaabnormal neuronal activity in ofisaevidence seizure including (Mu˜ nozthe et wake al., 2007) and altered depending on its effect on its neighbors. In this work we focus on a neural field consistent of exclusively an excitatory population and an adaptation term. The adaptation cluding hyper-excitation (Mu˜ n oz et al., 2007) and altered abnormal neuronal activity in wake of a and seizure in- tory population and an adaptation term. The adaptation cluding hyper-excitation (Mu˜ nozthe et al.,inhibitory 2007) altered interactions between excitatory and neurons focus on a neural field consistent of exclusively an excitatory population and an adaptation term. The adaptation interactions between excitatory and inhibitory neurons prevents activity of population to cluding hyper-excitation (Mu˜ noz et 2007) and altered term term prevents the the activity of the the excitatory excitatory population to interactions between excitatory and inhibitory neurons (Lopantsev, 2009; Engel, 1996). Inal.,clinical recordings, tory population and an adaptation term. The adaptation term prevents Neural the activity of the excitatory population to stay excited. field models are usually considered (Lopantsev, 2009; Engel, 1996). In clinical recordings, interactions between excitatory and inhibitory neurons (Lopantsev, 2009; Engel, 1996). In clinical recordings, stay excited. Neural field models are usually considered epilepsy is manifested in brain voltage activity organized term prevents Neural the activity of the excitatory population to stay excited. field models are usually considered on unbounded domains where conditions for the existence epilepsy is manifested in brain voltage activity organized (Lopantsev, Engel, 1996). In clinical recordings, unbounded domains where conditions for the existence epilepsy is manifested inincluding brain voltage activity organized in the form of2009; patterns, traveling waves (Stacey, on stay excited.waves Neural field models are usually on unbounded domains where conditions for theconsidered existence of traveling are easily established. However, some in the form of patterns, including traveling waves (Stacey, epilepsy is manifested inincluding brain activity organized in the Viventi form of et patterns, traveling waves (Stacey, 2012; al., Gonz´ aavoltage lez-Ram´ ırez et 2015). traveling waves are established. However, some on unbounded domains where conditions for the existence of traveling waves are easily easily established. However, some 2012; Viventi et al., 2011; 2011; Gonz´ lez-Ram´ ırez et al., al., 2015). of work has been developed to study neural field models in the form of et patterns, traveling waves (Stacey, has been developed to study neural field models 2012; Viventi al., 2011; Gonz´ alez-Ram´ ırez et al., 2015). In (Gonz´ alez-Ram´ ırez etincluding al., 2015) important features con- work of traveling waves areincluding easily However, some work has been developed to established. study neural field models on bounded domains boundary conditions. In (Gonz´ a lez-Ram´ ırez et al., 2015) important features con2012; Viventi et al., 2011; Gonz´ alez-Ram´ ırez features ettermination al., 2015). In (Gonz´ alez-Ram´ ırez et al., 2015) important con- on bounded domains including boundary conditions. The The cerning wave propagation preceding seizure work has been developed to on study neural field models on bounded domains including boundary conditions. The effects of boundary conditions the existence of patterns cerning wave propagation preceding seizure termination In (Gonz´ alez-Ram´ ırez al.,preceding 2015)clinical important features conof boundary conditions on the existence of patterns cerning wave propagation seizure termination were determined fromet in vivo recordings. The effects on bounded domains including conditions. The effects of boundary conditions onboundary thewaves existence of patterns of activity in the form of traveling or labyrinthine were determined from in vivo clinical recordings. The cerning wave propagation seizure termination of activity in the form of traveling waves or labyrinthine were determined from features in preceding vivo was clinical recordings. The existence of such wave observed right at the effects of boundary conditions on the existence of patterns of activity in the form of traveling waves or labyrinthine existence of such wave features was observed right at the has been shown to be fundamental (Laing, 2005; wereofdetermined from features in vivothe clinical recordings. existence such wave was observed at The the patterns end the of seizure, determining beginning ofright the seizure patterns has shown be 2005; of activity the form waves or(Laing, labyrinthine patterns hasinbeen been shownofto totraveling be fundamental fundamental (Laing, 2005; Gokce, 2017). end of the seizure, determining the beginning of the seizure existence of such wave features was observed right at the end of the seizure, determining the beginning of the seizure Gokce, 2017). termination. Traveling wave patterns of different widths patterns has been shown to be fundamental (Laing, 2005; Gokce, 2017). termination. Traveling wave patterns of different widths end of the seizure, determining the beginning of the seizure termination. Traveling wave of different and observed for some of in In this 2017). work, we study the initiation of traveling wave and speeds speeds were were observed for patterns some amount amount of time timewidths in the the Gokce, this we the of traveling wave termination. Traveling wave patterns of different widths In this work, work, we study study thetheinitiation initiation oftime traveling wave and speeds were observed for some amount of time in the In clinical recordings. activity in aa neural field in minimum by exciting clinical recordings. activity in neural field in the minimum time by exciting and speeds were observed for some amount of time in the In clinical recordings. this work, we study the initiation oftime traveling wave activity in a neural field in the minimum by exciting the neuronal population at the boundary of a bounded the neuronal population the boundary of a Relatively little clinical recordings. in a neural field inat minimum time exciting the neuronal population atthe theissue boundary of by aasbounded bounded Relatively little attention attention has has been been paid paid to to the the mechanisms mechanisms activity domain. This is an important to address it domain. This population is an important issue to address it may may Relatively attention has beenwave paidpropagation. to the mechanisms behind initiation of Some the neuronal at the boundary of termination aas bounded behind the thelittle initiation of cortical cortical wave propagation. Some domain. This is an important issue to address as it may help determine artificial ways to start seizure Relatively little attention has been paid to the mechanisms behind the initiation of cortical wave propagation. Some determine artificial ways to start seizure termination of the work found in literature dealing with mathematical help domain. This isartificial an important address as it that may help determine ways toissue starttoseizure termination sooner. Studies described in (Milton, 2003) suggest of found in literature dealing with mathematical behind theofinitiation of cortical wave propagation. Some of the the work work found in wave literature dealing with mathematical Studies described in (Milton, 2003) suggest that modeling cortical propagation mainly focuses on sooner. help determine artificial ways to start seizure termination sooner. Studies described in inputs (Milton, 2003) suggest that applying controlled external in the cortex can help modeling of cortical wave propagation mainly focuses on of the work found in wave literature dealing with mathematical modeling of cortical propagation mainly focuses on applying biological conditions for the propagation of wave activity controlled external inputs in cortex can help sooner. Studies described in (Milton, 2003) suggest applying controlled external inputs in the the cortex can that help biological conditions for the propagation of wave activity achieve the termination of an epileptic seizure, being modeling of cortical wave propagation mainly focuses on achieve the termination of an epileptic seizure, being this this biological conditions for the propagation of wave activity and not deals with the initiation of such activity (Goulet, applying controlled external in the cortexbeing can help achieve the termination of to aninputs epileptic seizure, this and not deals with the initiation of such activity (Goulet, a less invasive procedure treat the disease. That biological conditions forinitiation the propagation of the wave activity a less invasive procedure to treat the disease. That is, and not deals with the of such In activity (Goulet, is, 2011; Laing, 2005; Ermentrout, 2001). numerical achieve the termination of to an the epileptic this a less invasive procedure treat the seizure, disease.being Thatthat is, it is important to determine exact mechanisms 2011; Laing, 2005; Ermentrout, 2001). In the numerical and not dealsdeveloped with the initiation of suchthe activity (Goulet, is important to determine the exact mechanisms that 2011; Laing, 2005; Ermentrout, 2001). In the numerical simulations, for these works, wave activity it a less invasive procedure to treat the disease. That is, it is important to determine exact mechanisms that set the brain dynamics into the termination state. The simulations, developed for these works, the wave activity 2011; Laing, by 2005; Ermentrout, 2001). In in the numerical set the brain dynamics into the termination state. The simulations, developed foraa these works, the wave activity is initiated exciting fixed interval space that it is important to determine exact mechanisms that set the brain dynamics into the termination state. The is initiated by exciting fixed interval in space that understanding of these mechanisms may help conceiving simulations, developed fora these works, theinwave activity is initiated by exciting fixed interval space that models spatially concentrated high neuronal activity for understanding of these mechanisms may help conceiving set the brain dynamics into the termination The understanding of these mechanisms maytermination help state. conceiving ways to set the brain dynamics into the state models spatially concentrated high neuronal activity for is initiated by exciting a fixed interval in space that models spatially concentrated high neuronal activity for ways to set the brain dynamics into the termination state The authors acknowledge support from SIP IPN 20195566 and understanding of these mechanisms may help conceiving ways to set the brain dynamics into the termination state once the seizures initiate. The authors acknowledge supporthigh from neuronal SIP IPN 20195566 models spatially concentrated activity and for once the seizures initiate. The authors acknowledge support from SIP IPN 20195566 and 20194994. ways to set the brain dynamics into the termination state once the seizures initiate. 20194994. 20194994. The authors acknowledge support from SIP IPN 20195566 and once the seizures initiate.
20194994. 2405-8963 © © 2019 2019, IFAC IFAC (International Federation of Automatic Control) Copyright 153 Hosting by Elsevier Ltd. All rights reserved. Copyright 2019 IFAC 153 Control. Peer review© responsibility of International Federation of Automatic Copyright © under 2019 IFAC 153 10.1016/j.ifacol.2019.08.027 Copyright © 2019 IFAC 153
2019 IFAC CPDE-CDPS Oaxaca, Mexico, May 20-24, 2019
Laura R. González-Ramírez et al. / IFAC PapersOnLine 52-2 (2019) 150–155
In (Bennighof and Boucher, 1992) minimum-time control mechanisms have been developed to study wave propagation in a distributed system finding bang off bang timedependent control mechanisms applied at the boundaries of a finite domain. These characteristics define an openloop controller. Using these settings the aim of this paper is to establish a minimum-time control mechanism for the initiation of cortical wave propagation. However, in comparison with (Bennighof and Boucher, 1992) we deal with a non-local system. The non-locality of the system plays a fundamental role in determining the existence and therefore the initiation of waves. In particular, the control input that initiates wave propagation is associated with algebraic equations determined by the non-locality of the system. The goal of this paper is to understand the mechanisms behind the initiation of wave propagation in a relatively simple neural field as this can help to further study minimum time-initiation of wave propagation on more realistic models established on a bounded two dimensional domain or initiation of other spatiotemporal wave phenomena like multi-bump solutions or spiral wave solutions. The paper is organized as follows: in section 2 the mathematical model is introduced and the conditions necessary for the problem are stated. In section 3 the formulation of the problem is stablished. In section 4 a minimum-time control for the initiation of cortical wave propagation is derived. Finally, section 5 concludes the paper. 2. MATHEMATICAL MODEL OF CORTICAL WAVE PROPAGATION In (Gonz´ alez-Ram´ırez et al., 2015) a descriptive analysis of wave propagation preceding seizure termination was developed. Important quantitative properties of cortical waves were obtained finding wave widths varying from 2000 µm to 5000 µm and wave speeds varying from 80 µm/ms to 500 µm/ms as observed in high density local field potential data (LFP). An additional feature was also found and labeled “reverberation time” referring to the time between an initial, large amplitude wave and a subsequent, smaller amplitude fluctuation or “reverberation”. Also, in (Gonz´ alez-Ram´ırez et al., 2015) an activity-based neural field model was developed to describe the mean activity of the excitatory neuronal population in the form of traveling wave solutions of this model. This mathematical model is given below and describes cortical wave propagation preceding seizure termination. In contrast to (Gonz´alezRam´ırez et al., 2015) a bounded domain is considered. ∂u (x, t) = − αu (x, t) − βq (x, t) ∂t |x−y | 1 )u(y, t)dy − k exp (− (1) + αH 2σ Ω σ ∂q (x, t) =δu (x, t) − δq (x, t) ∂t
In this model u(x, t) represents the mean neuronal excitatory activity. To simplify the mathematical modeling of cortical wave propagation we do not consider an inhibitory population. The term q(x, t) represents the mean adaptation at position x and time t. The objective of the adaptation is to depress the neuronal activity of the 154
151
excitatory population. This term is not modeling a specific biological contribution to this system, but rather diverse factors that do not permit that activity remains in the excited state. H(·) is the Heaviside function; when the spatial neuronal interaction is large enough (i.e., when the convolution of u(x, t) and an exponential kernel exceeds the fixed threshold k) the Heaviside function activates. The dynamics of q(x, t) depend linearly on both variables. The neuronal interactions are restricted to a bounded interval Ω ∈ R. Five parameters are considered: α is the decay rate of the neural activity; σ determines the spatial extent of the neuronal interaction; k is the activation threshold for the Heaviside function; δ is the adaptation decay rate, β is the adaptation strength that has been rescaled by the decay rate. 2.1 Traveling Wave Solutions A traveling wave solution is a solution moving with a fixed shape, fixed width w and fixed speed c across the domain. That is, u(x, t) is a traveling wave solution if is a solution of the form u(x, t) = u (x − ct). To establish a traveling wave solution for (1) the system is rewritten in a moving coordinate frame z = x − ct, with c > 0 and we look for stationary solutions of this system. That is ∂u ∂t = 0. In an unbounded domain we also consider conditions at infinity limz→±∞ u(z) = 0 and limz→±∞ ∂u ∂z (z) = 0, so the wave solution stays bounded. In a bounded domain we cannot consider this conditions however we consider traveling wave solutions as obtained in an unbounded domain that are moving in a bounded domain. We now write system (1) in a moving coordinate frame: α β ∂u (z) = u (z) + q (z) ∂z c c α |z−y | 1 )u(y)dy − k − H exp (− c 2σ Ω σ ∂q δ δ (z) = − u (z) + q (z) ∂z c c
(2)
Since we are interested in pulse-like solutions, it is assumed that the Heaviside function activates in exactly two points, being these: z = w and z = w0 . Since the wave solution is translationally invariant we simplify our parameters choices and set w0 = 0. The points z = 0 and z = w indicate when the activity starts to increase (beginning of the wave) and when activity starts to decrease (end of the wave), respectively. These points also determine “the matching conditions”. The matching conditions are the conditions that need to be satisfied to generate a traveling wave. We establish this conditions in (11) and (12). The constant coefficient part of system (1), is determined by ∂u α β (z) = u (z) + q (z) ∂z c c ∂q δ δ (z) = − u (z) + q (z) ∂z c c
(3)
and is directly solved. This gives rise to three different cases depending on the nature of the eigenvalues of the corresponding matrix of constant coefficients. These different cases are determined by the relationship between 2 β and (α−δ) 4δ . This work focuses on the case of complex eigenvalues, since this case better describes the feature
2019 IFAC CPDE-CDPS 152 Oaxaca, Mexico, May 20-24, 2019
Laura R. González-Ramírez et al. / IFAC PapersOnLine 52-2 (2019) 150–155
of the “reverberation” of activity observed in the clinical 2 data. This case is obtained when β > (α−δ) 4δ . Therefore, it
Activity 0.8
The traveling wave solution is assumed that β > for the neuronal population is:
0.6
u
(α−δ)2 4δ .
u (z) =
where
0 ve (z − w) (c3 cos(zφ) + c4 sin(zφ)) + ve (z) (c1 cos(zφ) + c2 sin(zφ))
ve (z) = φ= c1 =
α α+β
√
√
3000
wR
if z ≤ 0
4δβ−(α−δ)2
4δβ−(α−δ)2 2c
α+δ )× 4δβ − (α − δ)2 + exp(−w 2c (2β + α − δ) sin(wφ) − 4δβ − (α − δ)2 cos(wφ)
α+δ )× c2 = (2β + α − δ) + exp(−w 2c − (2β + α − δ) cos(wφ) − 4δβ − (α − δ)2 sin(wφ) c3 = (2β + α − δ) sin(wφ) −
0.2
if 0 < z < w
α exp ( α+δ 2c )z (α+β)
0.4
if z ≥ w
4δβ − (α − δ)2 cos(wφ)
c4 = − (2β + α − δ) cos(wφ) − 4δβ − (α − δ)2 sin(wφ).
(4) (5) (6) (7)
(8) (9) (10)
Throughout the manuscript we denote the traveling wave solution of (1) by u . Solution (4) is established for traveling wave solutions moving in the positive z direction and it is considering w > 0; but a similar expression can be determined for wave solutions moving in the negative z direction. That is, in the moving coordinate frame obtained under the transformation z = x − ct traveling waves are moving to the right in the positive z direction. The existence of the traveling wave solutions is determined by the matching conditions. That is, the wave solution u is evaluated in the input of the Heaviside function at the points where the Heaviside function activates/deactivates, z = 0 and z = w. The matching conditions are: |y| 1 u (y)dy = k exp − (11) 2σ Ω σ 1 |w−y | u (y)dy = k exp − (12) 2σ Ω σ Considering a fixed threshold k satisfying conditions (11) and (12) the relationship between wave speed c and wave width w is defined by: 2k c2 + cσ (α + δ) + δσ 2 (α + β) w = σln 1 − σα (c + δσ) (13)
2000
1000
w0 0.2
1000
zu
2000
Distance Μm
w
R
0.4
Fig. 1. A traveling wave solution of model (1) traveling to the right. The marked points indicate important features of the wave solution. The point w indicates the beginning of the wave solution. The point zu indicates the point at which the highest amplitude of the wave is achieved. The point w0 determines the width of the wave (w w0 ). The point of highest amplitude is labeled as u ˆ. That is, since it is assumed w0 = 0 the width of the wave is determined by w. The point wR indicates the end of the “reverberation” considering that one period of the oscillatory component of the wave solution has passed by after w. Parameters used in this plot: α = 1, δ = 0.1, β = 5 and σ = 200.
placed at the ends of the boundaries in a finite domain. This approach has been studied in a distributed system in (Bennighof and Boucher, 1992). However in our case we deal with a nonlocal system so further considerations need to be addressed. Consider system (1) with initial and boundary conditions in a bounded domain Ω established as the interval (0, L) with L >> 5000 µm together with Neumann boundary conditions: u(x, 0) = 0, and q(x, 0) = 0, (14) ∂u ∂u (0, t) = FL (t), and (L, t) = FR (t), (15) ∂n ∂n ∂q ∂q (0, t) = 0, and (L, t) = 0, (16) ∂n ∂n where n is the unit normal vector pointing outward from Ω. Model (1) is an excitable system with zero as its equilibrium point. To see this consider the space-clamped version of model (1), such that, the kernel of the space convolution is constant. In this case the system turns into an ODE: ∂u = −αu − βq + αH ( cu − k ) ∂t ∂q = δu − δq ∂t
(17)
3. FORMULATION OF THE PROBLEM
In this system small perturbations from the stable state decay towards the rest state at zero and sufficiently large perturbations from the rest state produce a long trajectory of the system before achieving the rest state. This implies that system (17) is an excitable system. Therefore, if the system (1) starts with conditions (14) it is expected that the system remains in the rest state unless a sufficiently large input (15) is applied at the boundaries.
The goal of this paper is to establish a control mechanism for carrying out initiation of wave propagation in a minimum time of system (1) with features observed in clinical recordings by means of two bounded control inputs
We now return to the main problem (1) with initial and boundary conditions determined by (14), (15) and (16). We consider a problem of minimum time initiation of cortical wave propagation. We consider a final condition at time T determined by the propagation of a single wave
This relationship is restricted to pairs (c, w) that satisfy conditions (11) and (12).
155
2019 IFAC CPDE-CDPS Oaxaca, Mexico, May 20-24, 2019
Laura R. González-Ramírez et al. / IFAC PapersOnLine 52-2 (2019) 150–155
of activity of speed c and width w through the interval L. We also consider that after this propagation of activity the whole interval is tending asymptotically to the equilibrium state and the system is regulated by control inputs at the ends of the interval. Since the domain Ω is in R equation (15) turns into: ∂u ∂u − (x, t) (x, t) = FL (t), = FR (t) (18) ∂x ∂x x=0 x=L
where FL and FR are Lebesgue integrable in [0, T ]. The controls consist in a generating and an absorbing component, such that, given the initial conditions and the mentioned assumptions, the initiation of wave propagation to the right is guaranteed without considering wave interactions. Therefore, FL only consists of a component that generates right-moving waves and FR only consists of a component that absorbs activity at the right end of the interval to avoid reflections. Due to the initial conditions, we set the absorbing component at the right end of the interval to zero for at least the time it takes for a single wave to propagate through the interval. That is: L (19) FR (t) = 0, for 0 < t < c
153
4. MINIMUM-TIME CONTROL FOR THE INITIATION OF CORTICAL WAVE PROPAGATION Since we consider the problem of wave propagation on a biological system we establish the assumption that the parameters determining brain activity described in model (1) remain fixed. Thus, two possible scenarios are considered. In the first scenario the brain is under a parameter configuration that supports wave propagation with a desired wave speed c and wave width w, that is, the matching conditions (11) and (12) hold true for that choice of parameters. For a sufficiently large external input the system naturally evolves into the wave solutions of interest. In this case we consider a minimum-time control input that initiates wave propagation. In the second scenario the brain is under a parameter configuration in which equation (11) and (12) do not hold true. The propagation of a wave with a given speed c and width w is not guaranteed. In this case we propose a control mechanism to initiate the propagation of the wave solution of interest including the relevant feature of the “reverberation”. However, in comparison to the previous case we cannot assure the propagation of the wave throughout the whole domain. 4.1 Cortical dynamics supporting wave propagation
3.1 Control Inputs The control inputs can be set to the boundary conditions as: ∂ FL (t) = − (u (x − ct)) (20) ∂x x=0 ∂ FR (t) = (u (x − ct)) (21) ∂x x=L
where u are the traveling wave solutions for the excitatory population. Thus, the wave activity generated at the left end of the interval travels to the right (under the moving coordinate frame z = x − ct) not changing its shape. Traveling waves travel across the domain until they reach the right end of the boundary in which they are absorbed. We consider the following relationship to establish the wave solution in terms of the control inputs: 1 ∂ ∂ (u (x − ct)) = − (u (x − ct)) (22) ∂x c ∂t The control input at the left end of the interval is then defined as: 1 ∂ FL (t) = (u (x − ct)) (23) c ∂t x=0
The last equation determines a way to establish a control mechanism for the initiation of wave propagation. In particular, the profile of the wave at time t and position x = 0 in terms of the control input can be defined as: 1 t u (x − ct) |x=0 = FL (y)dy + u (−ct) (24) c 0
t=0
The last term of the previous equation is a constant that describes the displacement of the left end at t = 0. Due to our initial conditions this term is set to zero. Depending on the parameter configuration of the system there are two different scenarios for wave initiation and propagation that are discussed in the next section.
156
Considering the linear stability analysis of traveling wave solutions developed in (Gonz´alez-Ram´ırez et al., 2015) and the fact that the system is in the right regime for wave propagation it is only necessary to apply a sufficiently large excitation at the boundaries of the domain in order to start wave propagation. That is, the system naturally evolves into a wave solution for a sufficiently strong initial input. In order to achieve this initial input we analyze the profile of traveling wave solutions determined by model (1) (see Figure 1) and observe the point at which the wave achieves its highest amplitude. The objective of the control is to produce the initiation of the traveling wave by taking the system to the value of highest amplitude u ˆ (see Figure 1). Maximizing the value of the solution u (z) in the interval w < z < 0 we obtain that the maximum value of the wave α + is at u ˆ = −ce ve ( nπ φ − w ) + α+β for a fixed n ∈ Z . If we set an initial positive input FL = B for B < ∞ and use (24) in 0 < t < T0 we obtain 1 T0 u ˆ= Bdy. (25) c 0 Thus, the time T0 necessary for initiation of activity is determined by T0 = cˆ u/B. Lemma 1. Consider the fixed parameters α, β, σ, δ. Assume that δ = α/10 and that the following equation is satisfied: |w−y| u (y)dy = 0 (26) exp − |y| σ u (y)dy − Ω exp − σ Ω
where c and w are also fixed and u is as defined in (4). In order to have a minimum-time initiation of propagation of activity of a traveling wave with wave speed of c and width width of w the following boundary conditions determined by (1),(14),(15) and (16) are necessary: FL (t) =
B 0
if 0 ≤ t < T0 w − T0 L + tr + if T0 ≤ t ≤ c c
(27)
Laura R. González-Ramírez et al. / IFAC PapersOnLine 52-2 (2019) 150–155
(28) where B is a fixed positive external input and tr = √ 4π is the time it takes to propagate approxi2 4δβ−(δ−α)
mately one period of the oscillatory component of the wave solution (see Figure(1)). Under the previous assumptions we are guarantedd the propagation of a traveling wave solution with a fixed wave speed of c and fixed wave width of w throughout the interval (0, L) without changing its shape. The minimum time to start such propagation is α T0 = cˆ u/B where u ˆ = −ce ve ( nπ φ − w ) + α+β for a fixed n ∈ Z + that maximizes the neural activity. The value for the control FL is obtained by using (25). Here, the end of the “reverberation” is determined at w +ctr . In (25) we have a generating component to initiate propagation of activity followed by a zero component that no longer generates activity until the reverberating wave propagates through the interval L. In a similar way and using (21) we can determine the control at the right end of the interval (FR ). According to the main result in this section the propagation of activity is started as observed in the clinical data applying a fixed positive input of B at the left end of the interval in a minimum time of T0 . Under this setting, we aim to initiate seizure termination in that time by artificially starting traveling wave propagation.
0.8 0.6 0.4 Activity
0 if 0 ≤ t < L/c −B if L/c ≤ t < T0 + L/c FR (t) = −v (z) (c cos(zφ) + c sin(zφ)) if T ≤ t ≤ w − T0 + t + L e 1 2 0 r c c
0.2 0.0 �0.2 �0.4 �8000 �6000 �4000 �2000 0 Distance Μm
2000
4000
t=25 0.8 0.6 0.4 Activity
2019 IFAC CPDE-CDPS 154 Oaxaca, Mexico, May 20-24, 2019
0.2 0 -0.2 -0.4 0
2000 4000 6000 8000 10000 12000 Distance
Fig. 2. (Top) Traveling wave profile using the analytic solution (4)
4.2 Cortical dynamics not supporting wave propagation In this scenario the parameter configuration is such that it does not support wave propagation with speed c and width w. That is, for that choice of parameter conditions (11) and (12) are not satisfied simultaneously so that the existence and stability of traveling wave solutions with speed c and with w cannot be established. In this setting the propagation of activity with speed c and width w cannot be naturally initiated like in the previous section. This system could be able to propagate traveling waves with different features of speed and width or could not propagate activity in the form of traveling waves. Therefore, mechanisms for propagation of activity imply more artificial inputs being applied at the boundary of the finite domain. This also implies that the duration of a control applied at the boundary cannot be minimized as a control signal is necessary during the generation of the traveling wave profile. Also, the propagation of the traveling wave cannot be discussed during the length of the interval as the system does not necessarily support propagation of such wave. For completion, we describe the initiation of a traveling wave in this setting. However, a further study of the propagation of this wave needs to be developed. Using (23) and (4) the generating component for the desired solution can be explicitly determined. This is stated in the next result. Lemma 2. Consider the fixed parameters α, β, σ, δ. Assume that δ = α/10. In order to start the initiation of a traveling wave with wave width of w and wave speed of 157
satisfying conditions (11) and (12) for a wave speed of c = 338 µm/ms and wave width of w = 2350 µm. (Bottom) Snapshot at t = 25 ms of two distinct numerical simulations of model (1) using control (27) and (28) for generation and absorption of traveling wave activity. Both traveling wave profiles are initiated at time t = 0 ms. a) (red curve) Initial input B = 200 applied during T0 = 1.25 ms. b) (blue dotted curve) Initial input B = 50 applied during T0 = 4.87 ms. Notice that activity is tending asymptotically towards the zero state behind the wave profile. Parameters used in these plots: α = 1, δ = 0.1, β = 3 and σ = 95 µm.
w the following boundary conditions (15) are needed for system (1): α d v (z − w) (c cos(zφ) + c sin(zφ)) + c e 3 4 dt α+β if 0 ≤ t < tw FL (t) = c d (ve (z) (c1 cos(zφ) + c2 sin(zφ))) dt if tw < t < tw + tR L 0 if tw + tR < t < tw + tR + c
(29)
0 if 0 ≤ t < L/c α d ve (z − w) (c3 cos(zφ) + c4 sin(zφ)) + −c dt α+β L L ≤ t < t + if FR (t) = w c c d (v (z) (c cos(zφ) + c sin(zφ))) −c e 1 2 dt L L if tw + < t < tw + tR + c c
(30)
2019 IFAC CPDE-CDPS Oaxaca, Mexico, May 20-24, 2019
Laura R. González-Ramírez et al. / IFAC PapersOnLine 52-2 (2019) 150–155
where tw = wc is the necessary time to propagate the width of the wave. Under this setting it is not possible to establish that the produced wave will propagate throughout the whole interval. However, wave propagation is initiated. The previous lemma is obtained by using (21). These equations also help to determine the absorbing component at the right end of the interval. However, under this scenario the wave propagation throughout the whole interval cannot be assured and a more thorough investigation needs to be pursued. 5. CONCLUSION In this work a minimum-time control mechanism for the initiation of cortical wave propagation with features observed preceding seizure termination is established. In this way, the cortical dynamics in the minimum-time are set into the termination state. These mechanisms are developed assuming that cortical dynamics remain fixed and naturally permit the propagation of such cortical wave activity with the desired features. These control mechanisms consist of a generating component that permits the initiation of traveling wave activity and an absorbing component that returns activity to its rest state. Mechanisms to generate traveling wave activity with the desired features are also established, assuming that cortical dynamics do not permit the existence of wave activity. However, these mechanisms are more difficult to implement when treated as external inputs in cortical activity. Also, initiation of wave activity under this setting does not necessarily establish the propagation of such activity. Establishing the exact mechanisms that permit the propagation of activity with the desired features needs to be further addressed. In this work we consider the problem of minimum-time initiation of cortical wave propagation in a relatively simple neural field. This simple model motivates to further study minimum-time initiation of more complex patterns such as spiral waves or multi-bump solutions that appear in more complex neural fields. ACKNOWLEDGEMENTS The authors thank the reviewers for their careful reading of the manuscript and their constructive remarks. REFERENCES J.K. Bennighof and R.L Boucher Exact minimum-time control of a distributed system using a traveling wave formulation J Optimiz Theory App, 1992. J.J. Engel Excitation and inhibition in epilepsy Can J Neurol Sci, 1996. L.R. Gonz´ alez-Ram´ırez, O.J. Ahmed, S.S. Cash, C.E. Wayne and M.A. Kramer. A biologically constrained, mathematical model of cortical wave propagation preceding seizure termination. PLoS Comput Biol., volume 11, 2015. V. Lopantsev, M. Both, A. Draguhn. Rapid plasticity at inhibitory and excitatory synapses in the hippocampus induced by ictal epileptiform discharges Eur J Neurosci, 2009. J. Milton, P. Jung. Epilepsy as a Dynamic Disease Springer Verlag, 2003. 158
155
J. Goulet, GB. Ermentrout. The mechanisms for compression and reflection of cortical waves. Biol Cybern, 2011. C. Laing Spiral Waves in Nonlocal Equations SIAM J. Appl. Dyn. Syst., 2005. GB. Ermentrout, D. Kleinfeld Traveling Electrical Waves in Cortex: Insights from phase dynamics and speculation on a computational role Neuron, 2001. A. Gokce, D. Avitabile, S. Coombes. The dynamics of neural fields on bounded domains: an interface approach for dirichlet boundary conditions. J Math Neuro, 2017. A. Mu˜ noz, P. M´endez, J. DeFelipe, F. Alvarez-Leefmans. Cation-chloride cotransporters and gabaergic innervation in the human epileptic hippocampus Epilepsia, 2007. W. Stacey Better resolution and fewer wires discover epileptic spiral waves Epilepsy Curr, volume 12, 2012. J. Viventi, D. Kim, L. Vigeland, E. Frechette, J. Blanco et al. Flexible, foldable, actively multiplexed, high-density electrode array for mapping brain activity in vivo Nature Neuroscience, volume 14, 2011.