Wilson–Cowan model

The Wilson–Cowan model considers a homogeneous population of interconnected neurons of excitatory and inhibitory subtypes. All cells receive the same number of excitatory and inhibitory afferents, that is, all cells receive the same average excitation, x(t). The target is to analyze the evolution in time of number of excitatory and inhibitory cells firing at time t, $E(t)$ and $I(t)$ respectively.

The equations that describes this evolution are the Wilson-Cowan model:

$E(t+\backslash tau)=\backslash left[1-\backslash int\_\{t-r\}^\{t\}E(t\text{'})dt\text{'}\backslash right]\; \backslash ;\; S\_e\backslash left(\; \backslash int\_\{-\backslash infty\}^\{t\}\backslash alpha(t-t\text{'})[c\_1E(t\text{'})-c\_2I(t\text{'})+P(t\text{'})]dt\text{'}\backslash right\; )$

$I(t+\backslash tau)=\backslash left[\; 1-\backslash int\_\{t-r\}^\{t\}I(t\text{'})dt\text{'}\backslash right]\; \backslash ;\; S\_i\; \backslash left(\; \backslash int\_\{-\backslash infty\}^\{t\}\backslash alpha(t-t\text{'})[c\_3E(t\text{'})-c\_4I(t\text{'})+Q(t\text{'})]dt\text{'}\backslash right)$

where:

• $S\_e\backslash \{\backslash \}$ and $S\_i\backslash \{\backslash \}$ are functions of sigmoid form that depends on the distribution of the trigger thresholds (see below)
• $\backslash alpha(t)$ is the stimulus decay function
• $c\_1$ and $c\_2$ are respectively the connectivity coefficient giving the average number of excitatory and inhibitory synapses per excitatory cell; $c\_3$ and $c\_4$ its counterparts for inhibitory cells
• $P(t)$ and $Q(t)$ are the external input to the excitatory/inhibitory populations.

If $\backslash theta$ denotes a cell's threshold potential and $D(\backslash theta)$ is the distribution of thresholds in all cells, then the expected proportion of neurons receiving an excitation at or above threshold level per unit time is:

$S(x)=\backslash int\_\{0\}^\{x\}D(\backslash theta)d\backslash theta$,

that is a function of sigmoid form if $D()$ is unimodal.

If, instead of all cells receiving same excitatory inputs and different threshold, we consider that all cells have same threshold but different number of afferent synapses per cell, being $C(w)$ the distribution of the number of afferent synapses, a variant of function $S()$ must be used:

$S(x)=\backslash int\_\{\backslash frac\{\backslash theta\}\{x\}\}^\{\backslash infty\}C(w)dw$

If we denote by $\backslash tau$ the refractory period after a trigger, the proportion of cells in refractory period is$\backslash int\_\{t-r\}^\{t\}E(t\text{'})dt\text{'}$ and the proportion of sensitive (able to trigger) cells is $1-\backslash int\_\{t-r\}^\{t\}E(t\text{'})dt\text{'}$.

The average excitation level of an excitatory cell at time $t$ is:

$x(t)\; =\; \backslash int\_\{-\backslash infty\}^\{t\}\backslash alpha(t-t\text{'})[c\_1\; E(t\text{'})-c\_2\; I(t\text{'})+P(t\text{'})]dt\text{'}$

Thus, the number of cells that triggers at some time $E(t+\backslash tau)$ is the number of cells not in refractory interval, $1-\backslash int\_\{t-r\}^\{t\}E(t\text{'})dt\text{'}$ AND that have reached the excitatory level, $S\_e(x(t))$, obtaining in this way the product at right side of the first equation of the model (with the assumption of uncorrelated terms). Same rationale can be done for inhibitory cells, obtaining second equation.

When time coarse-grained modeling is assumed the model simplifies, being the new equations of the model:

$\backslash tau\backslash frac\{d\backslash bar\{E\}\}\{dt\}=-\backslash bar\{E\}+(1-r\backslash bar\{E\})S\_e[kc\_1\backslash bar\{E\}(t)-kc\_2\backslash bar\{I\}(t)+kP(t)]$

$\backslash tau\text{'}\backslash frac\{d\backslash bar\{I\}\}\{dt\}=-\backslash bar\{I\}+(1-r\text{'}\backslash bar\{I\})S\_i[k\text{'}c\_3\backslash bar\{E\}(t)-k\text{'}c\_4\backslash bar\{I\}(t)+k\text{'}Q(t)]$

where bar terms are the time coarse-grained versions of original ones.

The determination of three concepts is fundamental to an understanding of hypersynchronization of neurophysiological activity at the global (system) level:{{cite journal | pmid = 18643304 | doi=10.1103/PhysRevE.77.061911 | volume=77 | title=From baseline to epileptiform activity: a path to synchronized rhythmicity in large-scale neural networks | year=2008 | journal=Physical Review E | page=061911 | last1 = Shusterman | first1 = V | last2 = Troy | first2 = WC| issue=6 | bibcode=2008PhRvE..77f1911S }}

1. The mechanism by which normal (baseline) neurophysiological activity evolves into hypersynchronization of large regions of the brain during epileptic seizures

1. The key factors that govern the rate of expansion of hypersynchronized regions

1. The electrophysiological activity pattern dynamics on a large-scale

A canonical analysis of these issues, developed in 2008 by Shusterman and Troy using the Wilson–Cowan model, predicts qualitative and quantitative features of epileptiform activity. In particular, it accurately predicts the propagation speed of epileptic seizures (which is approximately 4–7 times slower than normal brain wave activity) in a human subject with chronically implanted electroencephalographic electrodes.V. L. Towle, F. Ahmad, M. Kohrman, K. Hecox, and S. Chkhenkeli, in Epilepsy as a Dynamic Disease, pp. 69–81{{cite journal | last1 = Milton | first1 = J. G. | year = 2010 | title = Mathematical Review: From baseline to epileptiform activity: A path to synchronized rhythmicity in large-scale neural networks by V. Shusterman and W. C. Troy (Phys. Rev. E 77: 061911) | journal = Math. Rev. | volume = 2010 | page = 92025 }}

The transition from normal state of brain activity to epileptic seizures was not formulated theoretically until 2008, when a theoretical path from a baseline state to large-scale self-sustained oscillations, which spread out uniformly from the point of stimulus, has been mapped for the first time.

A realistic state of baseline physiological activity has been defined, using the following two-component definition:

(1) A time-independent component represented by subthreshold excitatory activity E and superthreshold inhibitory activity I.

(2) A time-varying component which may include singlepulse waves, multipulse waves, or periodic waves caused by spontaneous neuronal activity.

This baseline state represents activity of the brain in the state of relaxation, in which neurons receive some level of spontaneous, weak stimulation by small, naturally present concentrations of neurohormonal substances. In waking adults this state is commonly associated with alpha rhythm, whereas slower (theta and delta) rhythms are usually observed during deeper relaxation and sleep. To describe this general setting, a 3-variable $(u,I,v)$ spatially dependent extension of the classical Wilson–Cowan model can be utilized.{{cite journal | last1 = Wilson | first1 = H. R. | last2 = Cowan | first2 = J. D. | year = 1973 | title = A mathematical theory of the functional dynamics of cortical and thalamic nervous tissue | journal = Kybernetik | volume = 13 | issue = 2| pages = 55–80 | pmid = 4767470 | doi = 10.1007/BF00288786 | s2cid = 292546 }} Under appropriate initial conditions, the excitatory component, u, dominates over the inhibitory component, I, and the three-variable system reduces to the two-variable Pinto-Ermentrout type model{{cite journal | last1 = Pinto | first1 = D. | last2 = Ermentrout | first2 = G. B. | year = 2001 | title = Spatially Structured Activity in Synaptically Coupled Neuronal Networks: I. Traveling Fronts and Pulses | journal = SIAM J. Appl. Math. | volume = 62 | page = 206 | doi = 10.1137/s0036139900346453 | citeseerx = 10.1.1.16.6344 }}

$\{\backslash partial\; u\; \backslash over\; \backslash partial\; t\}=u-v+\; \backslash int\_\{R^2\}\backslash omega(x-x\text{'},y-y\text{'})f(u-\backslash theta)\backslash ,dxdy\; +\; \backslash zeta(x,y,t),$

$\{\backslash partial\; v\; \backslash over\; \backslash partial\; t\}=\backslash epsilon\; (\backslash beta\; u-v).$

The variable v governs the recovery of excitation u; $\backslash epsilon>0$ and $\backslash beta>0$ determine the rate of change of recovery. The connection function $\backslash omega(x,y)$ is positive, continuous, symmetric, and has the typical form $\backslash omega=Ae^\{-\backslash lambda\backslash sqrt\; \{-(x^2+y^2)\}\}$. In Ref. $(A,\backslash lambda)=(2.1,1).$ The firing rate function, which is generally accepted to have a sharply increasing sigmoidal shape, is approximated by $f(u-\backslash theta)=H(u-\backslash theta)$, where H denotes the Heaviside function; $\backslash zeta(x,y,t)$ is a short-time stimulus. This $(u,v)$ system has been successfully used in a wide variety of neuroscience research studies.{{cite journal | last1 = Folias | first1 = S. E. | last2 = Bressloff | first2 = P. C. | year = 2004 | title = Breathing pulses in an excitatory neural network | journal = SIAM J. Appl. Dyn. Syst. | volume = 3 | issue = 3| pages = 378–407 | doi = 10.1137/030602629 | bibcode = 2004SJADS...3..378F | hdl = 10044/1/107913 | hdl-access = free }}{{cite journal | last1 = Folias | first1 = S. E. | last2 = Bressloff | first2 = P. C. | date = Nov 2005 | title = Breathers in two-dimensional neural media | journal = Phys. Rev. Lett. | volume = 95 | issue = 20| page = 208107 | pmid = 16384107 | doi = 10.1103/PhysRevLett.95.208107 | bibcode = 2005PhRvL..95t8107F | hdl = 10044/1/107604 | hdl-access = free }}{{cite journal | last1 = Kilpatrick | first1 = Z. P. | last2 = Bressloff | first2 = P. C. | year = 2010 | title = Effects of synaptic depression and adaptation on spatiotemporal dynamics of an excitatory neuronal network | journal = Physica D | volume = 239 | issue = 9| pages = 547–560 | doi = 10.1016/j.physd.2009.06.003 | bibcode = 2010PhyD..239..547K | hdl = 10044/1/107212 | hdl-access = free }}{{cite journal | last1 = Laing | first1 = C. R. | last2 = Troy | first2 = W. C. | year = 2003| title = PDE methods for non-local models | journal = SIAM J. Appl. Dyn. Syst. | volume = 2 | issue = 3| pages = 487–516 | doi = 10.1137/030600040 | bibcode = 2003SJADS...2..487L }} In particular, it predicted the existence of spiral waves, which can occur during seizures; this theoretical prediction was subsequently confirmed experimentally using optical imaging of slices from the rat cortex.{{cite journal | last1 = Huang | first1 = X. | last2 = Troy | first2 = W. C. | last3 = Yang | first3 = Q. | last4 = Ma | first4 = H. | last5 = Laing | first5 = C. R. | last6 = Schiff | first6 = S. J. | last7 = Wu | first7 = J. Y. | date = Nov 2004 | title = Spiral waves in disinhibited mammalian neocortex | journal = J. Neurosci. | volume = 24 | issue = 44| pages = 9897–902 | doi = 10.1523/jneurosci.2705-04.2004 | doi-access = free | pmid = 15525774 | pmc = 4413915 }}

The expansion of hypersynchronized regions exhibiting large-amplitude stable bulk oscillations occurs when the oscillations coexist with the stable rest state $(u,v)=(0,0)$. To understand the mechanism responsible for the expansion, it is necessary to linearize the $(u,\; v)$ system around $(0,0)$ when $\backslash epsilon>0$ is held fixed. The linearized system exhibits subthreshold decaying oscillations whose frequency increases as $\backslash beta$ increases. At a critical value $\backslash beta^\{*\}$ where the oscillation frequency is high enough, bistability occurs in the $(u,v)$ system: a stable, spatially independent, periodic solution (bulk oscillation) and a stable rest state coexist over a continuous range of parameters. When $\backslash beta\backslash ge\backslash beta^\{*\}$ where bulk oscillations occur, "The rate of expansion of the hypersynchronization region is determined by an interplay between two key features: (i) the speed c of waves that form and propagate outward from the edge of the region, and (ii) the concave shape of the graph of the activation variable u as it rises, during each bulk oscillation cycle, from the rest state u=0 to the activation threshold. Numerical experiments show that during the rise of

u towards threshold, as the rate of vertical increase slows down, over time interval $\backslash Delta\; t,$ due to the concave component, the stable solitary wave emanating from the region causes the region to expand spatially at a Rate proportional to the wave speed. From this initial observation it is natural to expect that the proportionality constant should be the fraction of the time that the solution is concave during one cycle." Therefore, when $\backslash beta\backslash ge\backslash beta^\{*\}$, the rate of expansion of the region is estimated by

$Rate\; =(\backslash Delta\; t/T)*c~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~(1)$

where $\backslash Delta\; t$ is the length of subthreshold time interval, T is period of the periodic solution; c is the speed of waves emanating from the hypersynchronization region. A realistic value of c, derived by Wilson et al.,{{cite journal | last1 = Wilson | first1 = H. R. | last2 = Blake | first2 = R. | last3 = Lee | first3 = S. H. | year = 2001 | title = Dynamics of travelling waves in visual perception | journal = Nature | volume = 412 | issue = 6850| pages = 907–910 | doi = 10.1038/35091066 | pmid = 11528478 | bibcode = 2001Natur.412..907W | s2cid = 4431136 }} is c=22.4 mm/s.

How to evaluate the ratio $\backslash Delta\; t/T?$ To determine values for $\backslash Delta\; t/T$ it is necessary to analyze the underlying bulk oscillation which satisfies the spatially independent system

$\{\{du\}\; \backslash over\; \{dt\}\}=u-v+H(u-\backslash theta),$

$\{\{dv\}\; \backslash over\; \{dt\}\}=\backslash epsilon\; (\backslash beta\; u-v).$

This system is derived using standard functions and parameter values $\backslash omega=2.1e^\{-\backslash lambda\backslash sqrt\; \{-(x^2+y^2)\}\}$, $\backslash epsilon=0.1$ and $\backslash theta=0.1$ Bulk oscillations occur when $\backslash beta\; \backslash ge\; \backslash beta^\{*\}=12.61$. When $12.61\; \backslash le\; \backslash beta\; \backslash le\; 17$, Shusterman and Troy analyzed the bulk oscillations and found $0.136\; \backslash le\; \backslash Delta\; t/T\; \backslash le\; 0.238$. This gives the range

$3.046\; mm/s\; \backslash le\; Rate\; \backslash le\; 5.331\; mm/s~~~~~~~~~~~~(2)$

Since $0.136\; \backslash le\; \backslash Delta\; t/T\; \backslash le\; 0.238$, Eq. (1) shows that the migration Rate is a fraction of the traveling wave speed, which is consistent with experimental and clinical observations regarding the slow spread of epileptic activity.{{page needed|date=January 2014}}Epilepsy as a Dynamic Disease, edited by J. Milton and P. Jung, Biological and Medical Physics series Springer, Berlin, 2003. This migration mechanism also provides a plausible explanation for spread and sustenance of epileptiform activity without a driving source that, despite a number of experimental studies, has never been observed.

The rate of migration of hypersynchronous activity that was experimentally recorded during seizures in a human subject, using chronically implanted subdural electrodes on the surface of the left temporal lobe, has been estimated as

$Rate\; \backslash approx\; 4\; mm/s$,

which is consistent with the theoretically predicted range given above in (2). The ratio $Rate/c$ in formula (1) shows that the leading edge of the region of synchronous seizure activity migrates approximately 4–7 times more slowly than normal brain wave activity, which is in agreement with the experimental data described above.

To summarize, mathematical modeling and theoretical analysis of large-scale electrophysiological activity provide tools for predicting the spread and migration of hypersynchronous brain activity, which can be useful for diagnostic evaluation and management of patients with epilepsy. It might be also useful for predicting migration and spread of electrical activity over large regions of the brain that occur during deep sleep (Delta wave), cognitive activity and in other functional settings.

{{reflist|30em}}

{{DEFAULTSORT:Wilson-Cowan model}}

Category:Computational neuroscience

Category:Epilepsy

Category:Neural circuitry

Category:Neurophysiology

Category:Electrophysiology