Excitatory-Inhibitory Homeostasis and Bifurcation Control in the Wilson-Cowan Model of Cortical Dynamics

Excitatory-Inhibitory Homeostasis and Bifurcation Control in the Wilson-Cowan Model of Cortical Dynamics | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results Excitatory-Inhibitory Homeostasis and Bifurcation Control in the Wilson-Cowan Model of Cortical Dynamics View ORCID Profile Francisco Páscoa dos Santos , View ORCID Profile Paul FMJ Verschure doi: https://doi.org/10.1101/2024.06.05.597312 Francisco Páscoa dos Santos 1 Eodyne Systems SL , Barcelona, Spain 2 Department of Information and Communication Technologies, Universitat Pompeu Fabra (UPF) , Barcelona, Spain Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Francisco Páscoa dos Santos For correspondence: f.pascoadossantos{at}gmail.com Paul FMJ Verschure 3 Donders Institute for Brain, Cognition and Behavior, Radboud University , Nijmegen, The Netherlands Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Paul FMJ Verschure Abstract Full Text Info/History Metrics Data/Code Preview PDF Abstract Although the primary function of excitatory-inhibitory (E-I) homeostasis is the maintenance of mean firing rates, the conjugation of multiple homeostatic mechanisms is thought to be pivotal to ensuring edge-of-bifurcation dynamics in cortical circuits. However, computational studies on E-I homeostasis have focused solely on the plasticity of inhibition, neglecting the impact of different modes of E-I homeostasis on cortical dynamics. Therefore, we investigate how oscillations and edge-of-bifurcation dynamics are shaped by the diverse mechanisms of E-I homeostasis employed by cortical networks. Using the Wilson-Cowan model, we explore how distinct modes of E-I homeostasis maintain stable firing rates in models with varying levels of input and how it affects circuit dynamics. Our results confirm that E-I homeostasis can be leveraged to control edge-of-bifurcation dynamics and that some modes of homeostasis maintain mean firing rates under higher levels of input by modulating the distance to the bifurcation. Additionally, relying on multiple modes of homeostasis ensures stable activity while keeping oscillation frequencies within a physiological range. Our findings tie relevant features of cortical networks, such as E-I balance, the generation of gamma oscillations, and edge-of-bifurcation dynamics, under the framework of firing-rate homeostasis, providing a mechanistic explanation for the heterogeneity in the distance to the bifurcation found across cortical areas. In addition, we reveal the functional benefits of relying upon different homeostatic mechanisms, providing a robust method to regulate network dynamics with minimal perturbation to the generation of gamma rhythms and explaining the correlation between inhibition and gamma frequencies found in cortical networks. Author summary We study how excitatory-inhibitory (E-I) homeostasis controls edge-of-bifurcation dynamics in cortical networks and how it impacts the generation of gamma oscillations. Importantly, while previous studies have limited E-I homeostasis to the plasticity of inhibition, we explore the wide range of mechanisms employed by cortical networks and, more importantly, how they interact. Here, we derive the mathematical solution for the Wilson-Cowan model under distinct modes of homeostasis and study how they shape model dynamics and the generation of gamma oscillations. That said, we demonstrate that E-I homeostasis, particularly of excitation and intrinsic excitability, modulates model dynamics relative to the bifurcation between damped and sustained oscillations in a manner previously unaccounted for, providing a mechanism for the implementation of heterogeneous distances to the bifurcation across cortical areas. Furthermore, our results stress the functional benefits of relying on multiple modes of homeostasis, allowing for the control of firing rates and circuit dynamics while ensuring that gamma oscillations remain within a physiological range and explaining the relationship between inhibition and gamma frequencies found in empirical data. With these results, we unify E-I balance, edge-of-bifurcation dynamics, and gamma oscillations under the lens of firing-rate homeostasis. Introduction The human cortex comprises repeated motifs of excitatory-inhibitory (E-I) interactions [ 1 – 3 ]. For this reason, cortical excitatory pyramidal (PY) neurons are constantly bombarded by excitatory and inhibitory inputs, mostly from neurons in local circuits [ 4 – 8 ]. More importantly, there is extensive evidence indicating that the excitatory and inhibitory inputs to single neurons are, by and large, tightly balanced [ 9 – 12 ] and that ensuring balanced E-I interactions is essential for the operations that cortical networks perform [ 13 – 16 ]. In addition, there is strong evidence suggesting that cortical networks homeostatically regulate E-I balance, ensuring the stability of firing rates in PY neurons [ 17 , 18 ]. The precise mechanisms of E-I homeostasis can take many forms, from synaptic scaling of excitatory [ 19 , 20 ] and inhibitory synapses [ 21 , 22 ] to the plasticity of the intrinsic excitability of PY neurons [ 23 – 25 ]. In addition, the mechanisms that maintain the stability of PY firing rates have been found to contribute to the self-regulation of cortical networks toward criticality [ 18 , 26 ]. This is particularly relevant since E-I balance is linked to edge-of-bifurcation dynamics [ 26 – 29 ], which are thought to be a foundational mechanism underlying the emergence of collective dynamics, such as metastability, in the cortex [ 27 , 30 – 34 ]. Furthermore, recent studies suggest that distinct cortical areas might be poised in different dynamical regimes relative to the bifurcation [ 35 ], reflecting their position in the cortical hierarchy [ 7 ]. Importantly, this heterogeneity in the distance to the bifurcation of cortical areas has been suggested to contribute to the emergence of metastable large-scale dynamics [ 32 ], suggesting a functional relevance beyond the local circuits. However, it is not clear how local cortical networks self-regulate not only toward the edge of bifurcation, but also how this set-point is determined across cortical areas to ensure heterogeneity in their dynamical regime relative to the bifurcation. That said, understanding how cortical networks regulate edge-of-bifurcation dynamics through E-I homeostasis might be essential to understanding large-scale cortical dynamics. In this context, several studies have endeavored to study E-I balance in large-scale brain models, demonstrating its role in supporting criticality [ 36 ] and metastability [ 37 ], and in shaping functional resting-state networks [ 37 – 39 ]. Furthermore, our results suggest that E-I balance is not only relevant for the emergence of spatiotemporal dynamics across timescales [ 40 ], but that E-I homeostasis can contribute to the recovery of functional networks following perturbations to the structural connectome [ 41 ]. There are, however, some caveats common to all these approaches. The first is that they commonly rely only on the regulation of feedback inhibition [ 36 – 39 ], even though it is known that E-I homeostasis is not limited to inhibitory synapses alone [ 17 , 18 , 42 ]. More specifically, it can also operate through the modulation of excitatory synapses [ 19 ] or intrinsic excitability [ 23 ]. Furthermore, it has been speculated that the concurrent action of multiple homeostatic mechanisms is relevant to avoiding drastic changes in neuronal dynamics [ 43 ] and to maintaining network dynamics close to criticality [ 18 , 42 ]. Importantly, our results suggest that homeostasis of inhibition alone is not sufficient for recovering network dynamics after lesions to the connectome [ 41 ]. However, to our knowledge, other modes of E-I homeostasis have not yet been validated in large-scale models. Furthermore, another gap in the literature is the lack of detailed exploration of the impacts of homeostasis on model dynamics. For example, manually regulating the synaptic couplings of the Wilson-Cowan (WC) model, a common neural-mass model relying on E-I interactions [ 44 , 45 ], affects the frequency of oscillations and bifurcation dynamics of the model [ 46 ], suggesting that processes that modulate the same parameters, such as E-I homeostasis, might have similar effects. However, this possibility remains unexplored in the literature. Given the foundational relationship between the dynamics of individual cortical areas and the emergent behavior of large-scale cortical networks, we aim to elucidate the impact of E-I homeostasis on the dynamics and bifurcation portrait of cortical areas. Therefore, in this paper, we study the behavior of the WC model under different forms of homeostasis by deriving analytically the fixed points of the different modes of E-I homeostasis and studying their impact on node dynamics using bifurcation theory. In addition, we devise a method to compute the model’s fixed points when multiple modes of E-I homeostasis are implemented in conjunction, allowing for the study of the WC model under multiple modes of homeostasis and how their concerted action shapes oscillatory dynamics. Furthermore, we analyze the effects of each mode of homeostasis on node frequency, suggesting that the cooperation of multiple modes of homeostasis allows for minimal impact on node parameters, maintaining the oscillation frequency within the gamma range. More importantly, our results reveal how edge-of-bifurcation dynamics can be controlled in cortical networks through the homeostatic regulation of mean firing rate. In addition, some modes of homeostasis can modulate the distance to the bifurcation, depending on the magnitude of external inputs, while maintaining the same mean firing rate. Therefore, we advance E-I homeostasis as a mechanism not only for the maintenance of stable activity in cortical neurons, but also to ensure edge-of-bifurcation dynamics, as previously suggested [ 18 , 26 ], offering an explanation as to why different cortical areas might be poised in different regimes relative to criticality [ 32 , 35 ]. In summary, the work presented here not only provides a mathematical basis for implementing multiple modes of E-I homeostasis in the WC model, but also elucidates how cortical networks can maintain dynamics at the edge of bifurcation through E-I homeostasis and how the combination of multiple processes of E-I homeostasis can ensure the maintenance of stable firing rates with minimal impact on the intrinsic rhythms of cortical networks. Methods Wilson-Cowan Model The WC model represents the dynamics of reciprocal interactions between populations of excitatory and inhibitory neurons in the cortex [ 44 , 45 ]. Given the relevance of excitatory synaptic scaling as a mechanism of E-I homeostasis [ 17 , 19 , 21 , 22 , 47 , 48 ], we propose an updated form of the canonical WC model by adding a term G E to control the gain of all excitatory inputs (i.e. recurrent excitation and excitatory external inputs). Therefore, we consider the following model equations, describing the firing rates of coupled excitatory ( r E ( t )) and inhibitory ( r I ( t )) neural masses: where c EE represents the recurrent coupling in excitatory populations, c IE represents the coupling from excitatory to inhibitory populations and c EI represents the coupling from inhibitory to excitatory populations. I ext represents an external excitatory input to the excitatory population, used here to represent the input received from other cortical areas, which is by and large excitatory [ 49 – 51 ]. The sigmoid activation (or input-output) functions F E ( x ) and F I ( x ) can be defined as: where µ and σ represent the firing threshold and sensitivity of the neural masses, respectively. Here, we consider that the activation functions of excitatory and inhibitory populations can be dissimilar and, thus, we define F E ( x ) with parameters µ E and σ E for the excitatory neural mass and F I ( x ) with parameters µ I and σ I for its inhibitory counterpart. This modification allows for the implementation of plasticity of intrinsic excitability of the excitatory population only [ 23 , 47 ]. The default model parameters follow the implementation of [ 39 – 41 ] ( Table 1 ). Unless stated otherwise, all parameters have arbitrary units. The time constants τ E and τ I were tuned so that our default model generates intrinsic oscillations at ∼ 40Hz near the bifurcation. This was done to approximate the resonance of cortical networks at gamma rhythms through the interaction between PY cells and fast-spiking inhibitory interneurons [ 52 – 55 ], in line with previous modeling studies [ 40 , 41 , 56 , 57 ]. On another note, it is relevant to point out that the original definition by Wilson and Cowan [ 44 ] also accounts for self-inhibition of the inhibitory population. However, to ensure that the model displays a Hopf-Bifurcation behavior between damped and sustained oscillations, conforming with the edge-of-bifurcation hypothesis [ 26 – 29 ], it must be guaranteed that the inhibitory-to-inhibitory coupling is close to 0 or, at least, considerably smaller than the rest of the couplings in the model [ 45 ]. Here, similarly to previous approaches [ 39 – 41 ], we set that parameter to 0. All the code necessary to run the models can be consulted in https://gitlab.com/francpsantos/wc_node_homeostasis . View this table: View inline View popup Download powerpoint Table 1. Default Parameters of the Wilson-Cowan Model. Derivation of Nullclines and Fixed Points The fundamental function of E-I homeostasis is to ensure the stability of firing rates in PY neurons by adapting their synapses or excitability [ 17 ]. In dynamical models, this translates into adjusting the parameters until the system reaches a stable state corresponding to the target activity levels. Therefore, we start by deriving the expression for the fixed points of the WC model. The first step is to obtain the expression for the system nullclines, which describe all the values of ( r E , r I ) for which either or for which are equal to zero. The intersection of the system nullclines defines the system’s fixed points, where both derivatives equal 0. The expression for the r E nullcline, obtained from Eqs 1 and 2, is the following: Similarly, we can derive the nullcline of r I by equating to 0, yielding: Here, log denotes the natural logarithm and exp the exponential function. For more details on the derivation of the system’s nullclines, refer to S1 Appendix. In Figure 1a , we plot the nullclines for three examples of the WC model under different levels of external input I ext . Download figure Open in new tab Fig 1. Phase portraits of the Wilson-Cowan model and detection of fixed points. Considering the r E and r I nullclines, the fixed points of the system can be derived by computing the intersection of both lines, yielding the values of ( r E , r I ) for which both and . Therefore, the fixed point r E can be obtained by solving the following equality: Due to the highly non-linear character of the equation, it is not possible to find a closed-form solution. Therefore, we find the fixed points by evaluating Eq (5) for all values of r E between 0 and 1, with a step of 10 − 6 , and detecting the points where the function crosses 0 ( Fig 1b ). Then, the corresponding values of can be obtained by substituting r E by in Eq (4). Analytical Approach to E-I Homeostasis: Computation of Model Parameters at the Fixed-Point To approximate the effect of E-I homeostasis in cortical networks, the parameters of the WC model should be adapted to maintain the activity of excitatory populations at a target level, regardless of perturbations in other parameters [ 17 , 26 , 48 ]. Here, to evaluate the conditions under which E-I homeostasis can maintain activity stable, we focus on changes at the level of the external input I ext , as is commonly done in studies of E-I homeostasis through the use of sensory deprivation [ 26 , 48 ]. Here, we explore four different modes of homeostasis derived from empirical studies ( Fig 2 ) and, for each, we derive the expressions for obtaining the values of the parameters in question as a function of I ext and the target excitatory activity . Download figure Open in new tab Fig 2. Different Implementations of Excitatory-Inhibitory Homeostasis in the Wilson-Cowan Model For each mode of homeostasis, we present a diagram of the Wilson-Cowan model, on which we highlight the model components that are modulated by E-I homeostasis. For plasticity of intrinsic excitability, we explore two methods which either modulate the threshold ( µ E ) of the input-output function, or adapt its threshold and slope ( σ E ) in a coordinated manner. Starting with synaptic scaling of excitatory synapses, it is known that, following perturbations in the normal firing rates of PY cells, their excitatory synapses are uniformly scaled to return firing rates to pre-perturbation levels [ 19 , 21 ] through post-synaptic changes in excitatory receptors or synaptic spines [ 22 , 47 , 48 ]. In our implementation of the WC model, the strength of excitatory synapses on pyramidal populations can be scaled by adapting the parameter G E ( Fig 2 ). For a given combination of r E and I ext , we can then derive the homeostatic solution of the system - the value of G E necessary to ensure that r E corresponds to a fixed point under I ext - by isolating G E in Eq (5): Another common form of E-I homeostasis is the synaptic scaling of fast-spiking inhibitory synapses onto pyramidal neurons [ 12 , 17 , 21 ] ( Fig 2 ). Although it is not clear if synapses are adjusted in a pre-or post-synaptic manner [ 17 ] or if there is an overall reduction in the number of synapses [ 21 , 22 , 58 ], the function of this type of plasticity is also to maintain the firing rates of pyramidal neurons [ 17 ]. Similarly to G E , the equation for the homeostatic value of c EI , representing the inhibitory-to-excitatory coupling in the WC model, can be derived from Eq (5): Beyond synaptic scaling, there are mechanisms of E-I homeostasis that adapt the intrinsic excitability of PY neurons to counteract perturbations in firing rates, such as sensory deprivation or activity blockade [ 20 , 23 – 25 ]. While some studies suggest that homeostasis of intrinsic excitability modulates firing thresholds [ 20 ] ( Fig 2 ), others show, instead, concurrent modulation of the firing threshold and slope of input-output curves, so that the activity of a neuron with no external input remains the same [ 23 – 25 ] ( Fig 2 ). Regardless, it is not clear to which extent cortical networks rely on these mechanisms of homeostasis [ 21 ], with the possibility that they only come into play following substantial perturbations in neuronal activity [ 17 ] or during development [ 25 ]. That said, homeostasis of intrinsic excitability can be implemented in two manners in the WC model. The first is at the level of the firing threshold of excitatory populations µ E [ 20 ] ( Fig 2 ). In this case, reordering Eq (5) gives: In the second case, the implementation of E-I homeostasis should ensure that both the firing threshold µ E and slope σ E change in a coordinated manner so that F E (0) remains the same, in line with empirical results [ 23 – 25 ] ( Fig 2 and S1 Fig ). Starting with the default parameters for the uncoupled node ( Table 1 ) ( and ), we have: where K is a constant equal to .By substituting σ E for Kµ E in Eq (8) and isolating µ E , we obtain the expression for this mode of homeostasis as: In addition to these forms of homeostatic plasticity, we also explored the synaptic scaling of excitatory synapses in inhibitory interneurons [ 43 , 59 ], through modulation of the c IE parameter. However, due to the lack of robustness to higher levels of external input and the scarce empirical evidence, this mode of homeostasis is not explored in detail in this paper. Nonetheless, the value of c IE and model dynamics under different combinations of r E and I ext in the WC model with homeostasis of excitatory-to-inhibitory synapses can be consulted in S2 Appendix. Implementation of Multiple Modes of Homeostasis Evidence suggests that different modes of homeostasis contribute simultaneously to the maintenance of E-I balance in the human cortex [ 17 , 18 , 20 , 21 , 42 ]. Therefore, we explore the possibility of the mathematical derivation of the system’s fixed points under more than one type of E-I homeostasis. From Eq (7), it follows that the system can have multiple solutions of ( G E , c EI ) that satisfy the condition of maintaining or a given I ext . Here, we present a heuristic to select a combination of G E and c EI relying on the assumption that both forms of homeostasis have similar timescales τ homeo , which is in line with empirical results [ 20 , 21 , 26 ]. First, we define the putative equations describing the dynamics of G E and c EI as: where τ homeo is substantially larger than the time constants τ E and τ I , representing the slow timescale of E-I homeostasis [ 17 , 18 ] and ensuring separation of timescales between neural dynamics and homeostasis. Note that, when the firing rate r E is higher than the target ρ, G E will tend to decrease and c EI will tend to increase to bring firing rates toward lower levels, and vice-versa . Given that the magnitude of both derivatives is the same, G E and c EI have the same rate of variation at any point in time, albeit in opposite directions. Therefore, the absolute difference between steady-state G E and c EI and their initial values is always the same. For a more detailed mathematical proof and examples, refer to S3 Appendix. Note that, as opposed to the common implementation of c EI homeostasis [ 14 , 36 , 37 , 39 ], our mechanism does not depend on r I . This is relevant because, in this case, the assumption of equal variation from initial values would no longer hold. However, in S4 Appendix, we provide a detailed explanation as to why our implementation of c EI homeostasis is more sensible for this version of the WC model. That said, under the assumption of equal timescales for the plasticity of G E and c EI , we arrive at the following expression relating the steady state values of both parameters, accounting for the fact that they change in opposite directions: Therefore, by inserting this expression into Eq (7), substituting G E and reorganizing the terms, we arrive at the following equations for c EI and G E under synaptic scaling of both excitation and inhibition: where and represent the initial values of G E and c EI , corresponding to the default parameters of the WC model, presented in Table 1 . Furthermore, in S3 Appendix, we present an adaptation of our approach to accommodate different timescales for each homeostatic mechanism. Assuming that homeostasis of intrinsic excitability operates in a timescale similar to the other homeostatic mechanisms, it is possible to derive an analytical expression for all the parameters under simultaneous G E , c EI , and µ E homeostasis. In this case, we consider that: In simpler terms, if firing rates are above the target, the threshold of the activation function increases to adjust activity toward the target setpoint. By combining this expression with Eq (11), we arrive at: Then, we can combine Eq (15) with Eq (8), yielding: Similarly, for plasticity of µ E and σ E , we substitute (15) in (10), obtaining: Computation of Homeostatic Parameters as a Function of the Target Firing Rate In cortical networks, homeostatic plasticity mechanisms regulate the average firing rates of excitatory neurons [ 17 , 19 , 20 , 26 , 47 , 48 ]. In the WC model, when the solution of the system is a stable fixed point or a stable spiral, the average firing rate of the system is precisely equal to . However, when the target state corresponds to a limit cycle, the average firing rate differs from ( Fig 3 ). For this reason, it is necessary to derive a method to obtain the solution of the system as a function of the target firing rate ρ , instead of . There is, however, no analytical method able to provide the mean value of a limit-cycle. Therefore, to estimate the mean firing rate ρ corresponding to a given , the equations must be solved numerically and the mean r E computed from the data. Download figure Open in new tab Fig 3. Mean firing rate as a function of fixed point r E for the Wilson-Cowan model under homeostasis of G E and I ext = 2.. While the solid blue line represents the mean firing rate as a function of , the dashed black line corresponds to . Red dots show the iterations of the procedure described in section 4.4.1. On the right, we include a zoom-in on the plot around the chosen value of corresponding to ρ . That said, to obtain an estimate of the steady-state model parameters as a function of ρ and I ext , we apply the following iterative procedure, illustrated in Fig 3 , with ϵ = 10 − 10 : Define upper and lower limits for ( and ). For , compute the homeostatic parameters of the system under I ext using the equations derived in the previous sections. Solve the system numerically for 5 seconds and compute the average firing rate ⟨ r E ⟩ from the last 2 seconds of activity. If |⟨ r E ⟩ − ρ | ρ , restart the procedure from step 1, updating Conversely, if ⟨ r E ⟩ < ρ , restart from step 1 with . This method allows for the iterative approximation of the model parameters corresponding to a mean firing rate of ρ , within an arbitrarily small level of precision dictated by ϵ . Note that, for the lower values of , for which the system solution corresponds to a stable fixed point, the mean firing rate is precisely the same as the fixed point ( Fig 3 ). However, as the system enters the limit-cycle regime, the mean firing rate starts to grow faster than the fixed point and it is in this regime that the numerical estimation of parameters described above is necessary. Here, unless stated otherwise, the system was solved numerically using the Euler method with a time step of 0.1 ms. Derivation of Jacobian Matrix and Classification of System Dynamics After obtaining the fixed point solutions for each type of E-I homeostasis, we characterize network dynamics around the fixed points. The characterization of dynamics can be done with linear stability analysis [ 60 , 61 ], by linearizing the system around a given fixed point using the first-order approximation of its Taylor expansion. To simplify the expressions, let us consider and . With this notation, the first-order Taylor expansion of the WC system can be written as: The matrix containing the partial derivatives of f ( x ) and g ( x ) at the fixed point is the Jacobian ( J ) of the system and, through the analysis of the eigenvalues of J , one can evaluate the stability and dynamics of the system [ 60 , 61 ]. A detailed derivation of the Jacobian in the WC model can be consulted in S5 Appendix. The eigenvalues λ ± of J can be expressed as a function of the trace ( Tr ( J )) and determinant ( Det ( J )) of the matrix as follows: That said, the fixed point is stable if the real part of both eigenvalues of J is negative (i.e. Tr ( J ) 0) and the dynamics have periodicity if the eigenvalues have a non-zero imaginary component (i.e. Tr ( J ) 2 − 4 Det ( J ) < 0). In addition, it is relevant to point out that fixed points in the “unstable spiral” regime do not necessarily represent an unstable solution of the system. As pointed out in [ 44 ], the WC model is always constrained by a bounding surface due to the sigmoid non-linearity, which, in our case, limits values between 0 and 1. Therefore, considering the Poincaré-Bendixson theorem [ 60 ], when there is only one fixed point, corresponding to an unstable spiral, and a boundary over which the flux points inwards, the system will necessarily settle on a periodic orbit. Given the number of fixed points and the trace and determinant of J , we classify node dynamics in the following manner: Stable Fixed Point : only one fixed point, corresponding to , Tr ( J ) < 0 and Stable Spiral : only one fixed point, corresponding to , Tr ( J ) 0 and . Bistable : two stable equilibria, one of which is the target . Unstable : more than one fixed point with the target corresponding to an unstable solution of the system (unstable fixed point, unstable spiral, or saddle point). In this study, we focus particularly on how different modes of E-I homeostasis modulate the Hopf-bifurcation between the stable spiral and limit cycle regimes, given the relevance of edge-of-bifurcation dynamics [ 27 , 30 – 34 ], related to the bifurcation between stable activity and sustained oscillations [ 26 – 28 ]. Analysis of Intrinsic Oscillation Frequency To analyze the impact of homeostasis on oscillation frequency, we solve the system numerically for each combination of ρ and I ext . Here, we simulate node dynamics for 5 seconds and use the last 3 seconds to compute the power spectrum of node activity, from which we extract the peak frequency of oscillation. To ensure that we can measure the frequency of systems in the stable spiral regime, which, if unperturbed, do not show rhythmic activity, we add a single pulse perturbation with amplitude 1 to I ext at 3 seconds of simulation time. While it is possible to estimate analytically the frequency of the WC model around a fixed point, this approximation is only valid in the vicinity of the fixed point, and when the system settles in a periodic orbit, the oscillation frequency decreases as the system orbits away from the fixed point (S6 Appendix). For this reason, we analyze the effects of E-I homeostasis on oscillation frequency by solving the system numerically. Results Single and Multiple Modes of Homeostasis in the Wilson-Cowan Model In cortical networks, E-I homeostasis has the function of maintaining stable firing rates in the face of perturbations in activity levels of PY neurons [ 17 , 18 ]. In the experimental setting, such perturbations are often related to sensory deprivation [ 26 , 48 ], which perturbs the magnitude of sensory input reaching early sensory cortices. While E-I homeostasis returns population firing rates to their average values pre-perturbation, it is not clear what is the exact value of such target firing rates, although studies suggest they correspond to the activity of cortical networks at the edge of bifurcation [ 26 , 27 ]. For these reasons, we explore the response of the WC model relying on different modes of homeostasis, under different combinations of external input I ext and target firing rate ρ . To do that, we use the equations derived in the Methods section to estimate the steady state solution of different types of homeostasis as a function of I ext and ρ , together with the results of linear stability analysis around the fixed point ( Fig 4 ). Here, for simplicity, we show results for all values of ρ between 0 and 0.4, in steps of 0.005. However, we present model dynamics across values of ρ from 0 to 1 in S2 Fig . Download figure Open in new tab Fig 4. Behavior of the Wilson-Cowan Model under Single Modes of Homeostasis as a Function of I ext and target firing rate ρ . In all plots, dashed lines represent the transition from stable fixed point to stable spiral, while solid lines show the Andronov-Hopf bifurcation between a stable spiral and a limit cycle. The parameters were estimated using a timestep of 0.1 ms. For each mode of homeostasis, we include a diagram of the Wilson-Cowan model, indicating which model components are modulated by homeostatic plasticity. Blank areas relate to combinations of parameters for which the system has no solution (i.e. not possible to find a stable fixed point corresponding to ρ ). Our first conclusion is that single modes of homeostasis can be roughly divided into three groups, depending on how they shape the bifurcation between damped and sustained oscillations in the WC model. The first group includes G E homeostasis ( Fig 4a ) and the regulation of intrinsic excitability by modulation of µ E and σ E ( Fig 4d ). For both, homeostasis modulates the value of ρ where the Hopf bifurcation between damped and sustained oscillations is found. In this case, the bifurcation shifts toward higher values of ρ as the external input is increased. For example, in the case of µ E and σ E homeostasis, a target firing rate of ρ = 0.12 can correspond to limit-cycle oscillations, a stable spiral, or a stable fixed point. Therefore, models with higher I ext , under E-I homeostasis, are less likely to engage in sustained oscillations. The second group encompasses synaptic scaling of inhibition ( Fig 4b ), for which the bifurcation moves toward lower values of ρ as the external input is increased. Therefore, for two nodes with the same target firing rate and different inputs, the one with stronger input might be in the limit-cycle regime, while its counterpart displays damped oscillations. Therefore, in this case, homeostatically regulated models with higher external input are more likely to enter the limit cycle regime. Finally, for the third group, consisting of µ E homeostasis ( Fig 4c ), the bifurcation between the stable spiral and limit cycle regimes is only dependent on the target firing rate, occurring at ρ ≈ 0.12. In this case, model dynamics can be regulated solely through the target firing rate. Furthermore, our results also suggest that all four modes of homeostasis are robust to changes in I ext and ρ , being able to maintain activity close to the target firing rate, provided it is not too low (see, for example, blank space in Fig 4d ). However, we note that, in models relying on the homeostasis of inhibition ( Fig 4b ), the modulations in model parameters required to ensure that excitatory firing rates are maintained at ρ are orders of magnitude larger in comparison with the other modes of homeostasis. For example, for higher values of I ext and lower target firing rates ( ρ ), the value of c EI needs to be scaled to up to 250, which is two orders of magnitude larger than the default value of 2.5. As demonstrated in the following sections, this can have a strong impact on node dynamics, particularly the intrinsic oscillation frequency. Beyond the single modes of homeostasis, we investigate the dynamics of the WC model under multiple modes of plasticity, and how they are shaped by external inputs and the target firing rate ρ . The results for homeostasis of G E and c EI are presented in Fig 5a . Here, the first noteworthy result is homeostasis no longer requires substantial modulation of model parameters, particularly c EI , which can minimize the impact of E-I homeostasis in model dynamics that might be caused by substantial variation in the model parameters [ 43 ] ( S3 Fig ). On another note, the modulation of the bifurcation point by homeostasis is similar to what is observed for homeostasis of G E , shifting the bifurcation point toward higher values of ρ when the external input is increased. This indicates that, when combining homeostasis of parameters with multiplicative (i.e G E ) and subtractive (i.e. c EI ) effects on model activity, the effects on dynamics will be dominated by the multiplicative parameters. This is already visible for homeostasis of µ E and σ E ( Fig 4d ), where the parameter space is similarly dominated by the effect of changing σ E . Download figure Open in new tab Fig 5. Behavior of the Wilson-Cowan Model under Multiple Modes of Homeostasis as a Function of I ext and target firing rate ρ . In all plots, dashed lines represent the transition from stable fixed point to stable spiral, while solid lines show the Andronov-Hopf bifurcation between a stable spiral and a limit cycle. The parameters were estimated using a timestep of 0.1 ms. For each combination of modes, we include a diagram of the Wilson-Cowan model, indicating which model components are modulated by homeostatic plasticity. Blank areas relate to combinations of parameters for which the system has no solution (i.e. not possible to find a stable fixed point corresponding to ρ ). Regarding the homeostasis of G E , c EI , and the excitability threshold µ E ( Fig 5b ), we observe that the magnitude of variation of the model parameters is further reduced ( S3 Fig ). In addition, while the dominance of G E in shaping the parameter space is still present, the combined action of c EI and µ E is apparent in the extension of the limit-cycle region to higher values of I ext . Conversely, when plasticity of intrinsic excitability modulates µ E and σ E [ 23 , 24 ] ( Fig 5c ) , σ E has the opposite effect, extending the stable spiral region in comparison with the model with G E and c EI homeostasis only. Finally, it is worth noting that the combination of plasticity of µ E and σ E with other modes of homeostasis also makes it more robust to lower target firing rates (compare blank spaces in Figs 4d and 5c ). In conclusion, across all modes of homeostasis, the influence of multiplicative parameters such as G E and σ E dominates node dynamics, leading to a modulation of the bifurcation between the stable spiral and limit-cycle regimes, shifting it toward higher values of ρ when I ext increases. Conversely, for c EI homeostasis, the opposite effect occurs. Finally, for µ E homeostasis, our results suggest that the bifurcation only depends on ρ . Furthermore, when combining multiple modes of homeostasis, the effect of multiplicative parameters dominates, shifting the Hopf bifurcation to higher values of ρ as the external input is increased. On another note, the combination of multiple modes of homeostasis contributes to reducing the magnitude of parameter variations, which can be particularly strong for the plasticity of inhibitory synapses. Effect of E-I Homeostasis on the Natural Frequency of Oscillations The generation of rhythmic activity in the WC model depends on the push-pull between excitation and inhibition. Therefore, given that E-I homeostasis can modulate the strength of excitation and inhibition, it is relevant to explore the effect of homeostasis on the frequency of oscillation. Following the procedure introduced in the methods section, we analyze how the intrinsic frequency of oscillation of the WC model changes based on the combination of ρ and I ext , for different implementations of E-I homeostasis ( Fig 6 ). Download figure Open in new tab Fig 6. Peak frequency of oscillation of the Wilson-Cowan Model under Different Modes of Homeostasis as a Function of I ext and target firing rate ρ . We present the peak frequency of oscillation across the parameter space for all modes of homeostasis. The colormap is centered around 40 Hz, the frequency of oscillation of the default model. In addition, we also display the bifurcations between the stable fixed point and stable spiral regimes (solid line) and between the stable spiral and limit cycle (dashed line). Starting with G E homeostasis, we observe that, in the stable spiral regime, the frequency of oscillation is largely independent of I ext and increases up to ρ ∼ 0.27, where it reaches a peak close to 125 Hz, after which the frequency starts decreasing with ρ . As nodes enter the limit cycle, the frequency of oscillation also decreases with ρ , but to a lesser extent. In addition, for values of ρ close to 0.1, the frequency of oscillation corresponds to the low gamma band (30-60 Hz). Conversely, c EI homeostasis has a strong impact on the frequency of oscillations, which become particularly fast around the Hopf bifurcation. In this region of the parameter space, the model can display oscillations of up to 250 Hz due to the high magnitude of inhibition required to maintain activity at the target firing rate ( Fig 4b ). Since oscillations are generated through the reciprocal interactions between excitatory and inhibitory populations, increasing c EI has a pronounced impact on the frequency of oscillations. Accordingly, both forms of homeostasis of intrinsic excitability do not have such a strong effect on the frequency of oscillation. For example, for µ E homeostasis, the frequency only depends on ρ , peaking close to the Hopf bifurcation at ∼ 45 Hz. Conversely, coupled homeostasis of µ E and σ E impacts the frequency of oscillation, since σ E not only scales the external input but also the strength of the excitatory and inhibitory inputs to pyramidal populations. Therefore, this type of homeostasis modulates the frequency of oscillations more strongly than the adaption of µ E only. Nonetheless, since scaling σ E has the same effect on the inhibitory and excitatory couplings, the impact on the frequency is not as strong as synaptic scaling of either G E and c EI . In addition, similarly to c EI and µ E homeostasis, the frequency of oscillation peaks around the Hopf bifurcation. In addition to the impact of single modes of homeostasis, we further present how the conjugation of multiple modes shapes oscillations in the WC model ( Fig 6 ). Starting with G E and c EI homeostasis, since c EI does not vary as strongly as for c EI homeostasis (see Figs 4b and 5a ), the frequency remains closer to 40 Hz in a broader region of the parameter space, compared to c EI homeostasis, even though it still peaks at higher values. Here, the dominance of multiplicative parameters such as G E is further suggested by the similarity between the parameter spaces of G E and G E + c EI homeostasis. Furthermore, while the inclusion of µ E homeostasis has no substantial effects compared to the homeostasis of G E and c EI only, the implementation of plasticity of intrinsic excitability with µ E and σ E reduces the impact of homeostatic plasticity on the frequency of oscillation, which remains in the low-gamma range in a larger region of the explored parameter space. All in all, we demonstrate that, beyond the maintenance of target firing rates, different modes of homeostasis can have profound impacts on node dynamics such as the frequency of oscillations. Importantly, this result has implications for the oscillatory activity in cortical areas with different densities of excitatory and inhibitory receptors, which might arise due to the effect of E-I homeostasis. Therefore, we stress that this effect should be considered not only when analyzing models with E-I homeostasis but, more importantly, when interpreting results relating gamma oscillation frequency with the density of different receptors in the cortex. Effect of Integration Time Step (dt) on Node Dynamics The method derived to estimate the homeostatic values of model parameters as a function of the target firing rate ρ requires the numerical integration of the WC model. Therefore, from a methodological standpoint, it is relevant to investigate the impact of the integration time step ( dt ) on the results. More specifically, the rhythmic activity in the WC model is a consequence of the interplay between excitation and inhibition, leading to a constant push-pull that generates oscillations that might be damped or sustained. That said, given that the integration time step can influence how fast inhibition and excitation evolve, the behavior of the model can be affected, shifting the position of the bifurcation. For this reason, we analyze the position of the Andronov-Hopf bifurcation as a function of ρ and I ext and how it varies depending on the level of external input. Our results suggest that, for all modes of homeostasis and independently of I ext , the position Hopf-Bifurcation is affected by changing the integration dt ( Fig 7 ). Some models are minimally impacted, particularly when relying on either type of plasticity of intrinsic excitability, likely because these modes of homeostasis do not directly affect the E-I interactions in the model. Conversely, particularly in modes of homeostasis relying on the plasticity of c EI , the integration dt has a strong influence on the bifurcation, with slower dt values extending the limit-cycle regime to lower values of ρ . Download figure Open in new tab Fig 7. Hopf-Bifurcation as a function of the integration time step ( dt ). For each mode of homeostasis, dashed black lines represent the analytical bifurcation between the stable spiral and limit cycle regimes as a function of ρ and I ext . In addition, we present the bifurcation lines computed numerically using the Euler method with different integration time steps, ranging from 10 − 5 to 10 − 3 seconds. Importantly, when analyzing the analytical bifurcation portraits, it can be observed that c EI homeostasis can be more accurately grouped with µ E homeostasis as types of plasticity for which the bifurcation depends only on ρ . However, when integrating the system numerically, we observe the behavior reported in the previous section, where the bifurcation shifts toward lower values of ρ with increases in I ext . The existence of two groups of homeostatic mechanisms can be explained by the way each parameter shapes the external input in the WC model ( Eq 1 ). More specifically, as previously hinted at, c EI and µ E have a subtracting effect on the external input (direct in the case of µ E and indirect in the case of c EI ), while G E and σ E have a multiplicative relation with I ext . This results in the different dependence of dynamics on ρ and I ext . While here we focus on the dynamics of isolated WC neural masses, the implications of these differences on the network context should be subject to further study, especially regarding network stability. That said, since large-scale models are often solved numerically and given the non-negligible effect of integration dt , our results suggest that, when analyzing the bifurcation dynamics of the WC model, especially in the network context, the effect of dt should be accounted for and minimized. Furthermore, from the inspection of our results, values of dt lower than 0.2 milliseconds generally provide a sufficient approximation of the analytical behavior of the model while maintaining simulation times tractable. While this might not be the case for the homeostasis of c EI , where integration invariably affects the independence of the bifurcation on I ext , we suggest that the numerical bifurcation portrait, as opposed to the analytical one, should be considered when studying the dynamics of models relying only on c EI homeostasis. Bifurcation Control in Models with Fast and Slow Inhibition As stated in the previous section, rhythmic activity in cortical networks can be a consequence of the interplay between excitation and inhibition, leading to a constant push-pull that generates oscillations that might be damped or sustained. Our default WC model is based on the generation of gamma rhythms (∼40 Hz) by the PY-PV loop, due to the fast GABAergic synapses that fast-spiking PV interneurons establish in PY neurons [ 52 , 53 , 55 ]. However, other classes of inhibitory interneurons modulate pyramidal activity through GABA B receptors, which have slower dynamics [ 51 , 62 , 63 ]. Furthermore, this can lead to the resonance at rhythms slower than gamma (i.e. alpha/beta), as observed in the deeper layers of the cortex [ 64 – 67 ]. That said, E-I homeostasis may shape the dynamics of cortical networks differently depending on the time constants of inhibitory synapses. Therefore, we explore how the time constant of inhibition ( τ I ) interacts with the modulation of the Andronov-Hopf bifurcation by each mode of homeostasis ( Figure 8 ). Download figure Open in new tab Fig 8. Hopf-Bifurcation as a function of the time constant of the inhibitory population ( τ I ). For each mode of homeostasis, we present the analytical Hopf-bifurcation as a function of ρ and I ext , in models with τ I ranging from 5 to 100 ms. In all simulations, the excitatory time constant ( τ E ) was set to 2.5 ms. Our results indicate that, even though each mode of homeostasis displays a different bifurcation portrait, a common effect of increasing the inhibitory time constant can be discerned: as inhibition becomes slower, the bifurcation shifts toward lower values of ρ . Furthermore, while this effect is common across all modes of homeostasis, it is less strong in the model relying on the plasticity of intrinsic excitability through the threshold ( µ E ) and slope ( σ E ) of the input-output function (see µ E + σ E panel in Figure 8 ). The main implication of our result is that, generally and regardless of the implementation of E-I homeostasis, cortical networks with slower inhibition might be more likely to engage in sustained oscillations In summary, our results suggest that the critical point between damped and sustained oscillations is shaped not only by the mechanisms of E-I homeostasis but also by the time constant of inhibition. More specifically, for a given level of input I ext , slower inhibition shifts the bifurcation to lower target firing rates ( ρ ). The implications of this behavior in cortical networks are discussed in detail in the Discussion. Discussion In this work, we apply concepts from the field of dynamical systems, such as bifurcation theory and linear stability analysis, to explore the behavior of the Wilson-Cowan model under different modes of E-I homeostasis. To do this, we gather information from the literature on the mechanisms through which cortical networks regulate their balance to maintain stable firing rates [ 17 , 18 , 42 ] and translate them to equivalent adaptations in the parameters of the Wilson-Cowan model. This approach allows for the derivation of the model parameters that allow for systems to be poised at a given target firing rate ρ under different levels of incoming input. In addition, to account for the concurrent action of multiple homeostatic mechanisms and study their synergies, we derive a method to estimate steady-state parameters of models under multiple modes of homeostasis, under the assumption that all modes operate with the same timescale. One of the most prominent advantages of our approach is methodological and concerns the computational traceability of models with E-I homeostasis. Normally, the implementation of this form of plasticity in large-scale models requires that simulations are run for a substantial amount of simulation time, ensuring that model parameters reach a steady state (see, for example, [ 40 , 41 ]). Therefore, being able to estimate the steady-state solution of E-I homeostasis analytically may allow for a more extensive exploration of E-I balance in large-scale models, which would require massive computational resources otherwise. Furthermore, within the topic of multiple modes of homeostasis, it is worth pointing out that, even though we assumed equality of timescales between different types of homeostasis to reach the expressions for the homeostatic parameters [ 17 , 20 ], our mathematical framework allows for the estimation of steady-state parameters in models of different homeostatic timescales (S3 Appendix). Beyond the practical benefits of our approach, it also allows for a deeper understanding of the dynamics of the balanced Wilson-Cowan model and, consequently, of networks of coupled excitatory and inhibitory neurons in the cortex. For this reason, we explore the behavior of the model under homeostasis of excitatory-to-excitatory synapses, excitatory-to-inhibitory synapses, and intrinsic excitability, in addition to the most common implementation of homeostasis in inhibitory-to-excitatory connections [ 14 , 36 – 39 ]. In this context, our bifurcation analysis suggests that the modes of homeostasis can be divided into two groups. The first, including plasticity of c EI and µ E , relates to models where the bifurcation between a fixed point solution and a limit cycle is only dependent on the target firing rate ρ for most of the domain. In a network context, this would mean that all nodes could be poised in the same dynamical regime regardless of heterogeneities in external inputs. However, when simulating model activity by solving the system numerically, the integration dt has a strong effect on models with c EI homeostasis. More specifically, the bifurcation is shifted toward lower values of ρ as the external input is increased, meaning that nodes with higher inputs are more likely to be in the limit cycle regime. Conversely, in the other group, including plasticity of G E and coupled µ E and σ E , the bifurcation depends on both the external input and the target firing rate. Therefore, in the network context, nodes can maintain the same firing rate while being in different regimes relative to the bifurcation. Regarding the advantages of each type of plasticity, there is no clear picture at this level of analysis. On one hand, research suggests that cortical networks self-regulate toward quasi-criticality by regulating firing rates [ 18 , 26 ]. Under the assumption that cortical networks achieve this through the maintenance of a uniform set point firing rate, methods such as µ E homeostasis would be at an advantage. On the other hand, the study of [ 26 ] relates to sensory areas and it is unclear if all areas of the cortex have the same optimization goal point relative to criticality. In this context, some studies suggest that cortical dynamics are optimized when nodes have heterogeneous distances to the bifurcation [ 32 ]. This could be either achieved by having heterogeneous target firing rates or, in our model, through modes of homeostasis modulating multiplicative parameters (e.g. G E or σ E ), for which this heterogeneity emerges naturally from different levels of external input. Furthermore, recent results suggest that distance to criticality (equivalent to the distance to bifurcation) tracks the hierarchy of visual areas, in line with the idea of heterogeneity in distance to bifurcation across the cortex [ 35 ]. Within our framework, this could be a consequence of the interaction between different levels of input or heterogeneity in microcircuitry [ 7 ] and the effects of particular modes of homeostasis modulating the distance to the bifurcation. Therefore, future studies should focus on the comparison between the behavior of network models under the different modes of homeostasis, investigating how the different types of bifurcation portraits shape network dynamics. Nonetheless, our results agree with the edge-of-bifurcation theory of brain dynamics [ 27 , 30 , 31 , 68 , 69 ], demonstrating that bifurcation control can be implemented in cortical networks through homeostatic regulation of firing rates, following previous theoretical [ 27 ] and empirical studies [ 26 ]. More importantly, we offer a plausible explanation for the heterogeneity in the distance to bifurcation across cortical areas [ 32 , 35 ], by suggesting that this can result from the interaction between E-I homeostasis and local heterogeneities in incoming inputs or microcircuitry [ 7 ], without requiring heterogeneous target firing rates. Here, we also study how each mode of homeostasis shapes the frequency of oscillations generated by the model. Recent research shows that, in the WC model, the frequency is modulated by changes in the coupling parameters and external input [ 46 ]. This results from the fact that oscillations in the WC model emerge from the bidirectional interaction between E and I neurons. Shortly, the E population excites itself and the I population, which leads to an increase in E and I activity, until the latter reaches a point where it can compensate for the recurrent excitation and external input. At this point, both E and I activities start decreasing until the I activity is no longer able to compensate for the recurrent excitation, which leads to a new cycle. Therefore, it is to be expected that changes in G E or c EI , which regulate how quickly excitation and inhibition evolve, also lead to a modulation in the duration of cycles. For this reason, we study how homeostasis shapes the frequency of oscillations ( Fig 6 ). For example, our results demonstrate that c EI homeostasis, which is the most common method implemented in large-scale models [ 36 – 41 ], has a profound effect on the oscillatory behavior of the model. This is because changes in the external input require large modulations of c EI which strongly affects the frequency of oscillation. For this reason, even though the target firing rate might be common, nodes with different external inputs will oscillate at different frequencies. Conversely, for G E homeostasis, the frequency of oscillation has a smaller range of variation across different levels of input. While changing G E should modulate the frequency due to its effect on recurrent excitation, this is compensated by the modulation of external inputs by the same parameter. Therefore, through different types of homeostasis, nodes with varying inputs can have varying levels of heterogeneity in oscillation frequency. How does this relate to cortical dynamics? Results in the human cortex suggest that the frequency of gamma oscillators in the visual cortex correlates with the density of inhibitory receptors [ 70 ]. Furthermore, recent results show that this correlation extends to other cortical areas, possibly tracking the cortical hierarchy [ 71 ]. This is in line with our results regarding modes of homeostasis involving c EI , where the frequency of oscillation increases with the inhibitory coupling ( S4 Fig ). However, our results show that, when implementing homeostasis of inhibition only, the frequency can vary up to 200 Hz, which is outside the range of variation reported in [ 70 , 71 ] (∼35-100 Hz). Importantly, this issue can be solved through the integration of multiple modes of homeostasis. In this context, our results demonstrate that, since the magnitude of variation in parameters is smaller, there is a weaker impact of homeostasis on the frequency of oscillation when multiple modes are integrated. For example, under combined homeostasis of G E and c EI ( S4 Fig ), we find that the frequency still varies with the level of external input, and thus the strength of inhibition, but it remains within the range reported by [ 70 , 71 ]. This can be beneficial for long-range synchrony between cortical areas in the gamma range, which is essential for neural communication [ 72 ] and behavior [ 73 ]. That said, we argue that our results are in line with empirical data, showing a correlation between the strength of inhibition and the frequency of intrinsic gamma oscillations. More importantly, we provide a mechanistic explanation for this correlation, suggesting that it reflects the effects of E-I homeostasis maintaining stable firing rates in cortical areas that receive varying levels of external input or have differences in microcircuitry [ 7 ]. Importantly, our results also suggest that the timescale of inhibitory synapses can influence the bifurcation between damped and sustained oscillations. More specifically, given two systems with the same implementation of homeostasis, target firing rate ( ρ ), and external input ( I ext ), the one with slower inhibition is more likely to engage in sustained oscillations. This effect may be relevant for cortical dynamics, since cortical networks are known to generate slower (alpha/beta) oscillations in deep layers, while faster (gamma) rhythms are observed in superficial layers [ 64 – 67 ]. Therefore, our results on the time constant of inhibition should be interpreted in light of these laminar-specific rhythmic activity patterns. First, under the assumption that deep and superficial layers of the same local circuit have similar target firing rates and input levels, our results would suggest that deeper layers are more likely to engage in persistent oscillations. This effect could be counteracted, for example, by the maintenance of lower target firing rates in deep layers. However, extensive evidence suggests that the opposite is true, with higher spontaneous firing rates in deep versus superficial layers [ 74 – 78 ]. On the other hand, there is also the possibility that deeper layers receive stronger inputs which, when combined with the higher target firing rates, could set them in a similar regime relative to the bifurcation as superficial layers ( Figure 8 ). While models of the canonical inter-laminar microcircuit account for stronger external inputs to deeper layers [ 75 ], detailed anatomical studies suggest instead that the proportion of deep and superficial inputs received by cortical areas depends on their position on the cortical hierarchy [ 66 , 79 , 80 ]. However, in cortical networks, most excitation comes from local circuits [ 4 – 8 ]. In this context, inter-laminar interactions are dominated by excitatory projections from layer 3 to layer 5, while layer 5 neurons strongly project to layer 3 interneurons [ 66 , 75 ]. It is then reasonable to infer that pyramidal neurons in deeper layers might receive stronger excitatory inputs than their counterparts in superficial layers. Therefore, the stronger inputs and higher mean firing rates in deep layers could compensate for the effect of slower inhibition on the bifurcation of cortical dynamics, posing deep and superficial layers in a similar activity regime relative to the edge-of-bifurcation. This should, however, be studied with the appropriate level of detail. For that reason, we suggest that future studies focus not only on the laminar differences in mean firing rates and levels of input, but also investigate signatures of edge-of-bifurcation dynamics, such as long-range temporal correlations [ 81 , 82 ], in the slow and fast oscillations of deep and superficial layers. In addition, computational models accounting for the laminar patterns of fast and slow inhibition [ 66 , 67 ] should be employed to study how inter-laminar interactions, homeostatic plasticity and heterogeneity in external inputs and target firing rates affect the dynamics of local cortical circuits. Finally, it should be pointed out that SST interneurons, which are the main source of slow inhibition to PY neurons [ 51 ], show no signatures of participating in E-I homeostasis [ 12 , 21 ]. Nonetheless, the effect of inhibitory time constants on the bifurcation is evident even in models without homeostasis of inhibition (see G E or µ E panels in Figure 8 ). Therefore, we argue that it is still worth studying how E-I homeostasis shapes the dynamics of local circuits with fast and slow inhibition. On another note, although we derived the mathematical solution for c IE homeostasis, we did not delve deeply into the study of this type of plasticity. While some results point out that the excitatory synapses in PV neurons are also scaled in response to perturbations in PY firing rates [ 43 , 59 ], this mechanism may only come into play when activity is higher than the target [ 43 ]. Conversely, our analysis suggests that, at least for the WC model, it is not a sufficiently robust mode of homeostasis precisely when perturbations (i.e. higher I ext ) lead to higher firing rates. This is because homeostasis of c IE acts through the disynaptic modulation of E firing rates by increasing the excitation of the inhibitory population, which provides feedback inhibition. However, inhibitory activity has an upper limit due to the sigmoid non-linearity of the model, which means that inhibition cannot be increased beyond the maximum level and, therefore, this mode of homeostasis breaks down under higher levels of external input (see S4 Appendix). Nonetheless, although this possibility is not explored here, this mode of homeostasis might interact with other modes, contributing to the maintenance of E-I balance. A limitation of our approach is that we focus on the slow regulation of E-I balance and circuit dynamics through distinct mechanisms of E-I homeostasis [ 17 , 18 ]. However, it is known that E-I balance can be modulated in faster timescales. For example, during the formation of neural ensembles, fast E-I balance modulation can be recruited to regulate the speed-accuracy trade-off [ 16 ]. Similarly, recent results show that neuromodulatory signals can be used to adjust E-I balance in cortical networks and regulate their sensitivity to noise during decision-making tasks [ 83 ] or to control the gating of sensory inputs during learning [ 84 ]. Nonetheless, we argue that the slower homeostatic regulation of dynamics towards the edge of bifurcation is essential to support the fast modulation of E-I balance, by ensuring that cortical networks are poised in a regime where they can be quickly brought to either a state of high or low excitability depending on task demand. For this reason, we suggest that, while homeostasis does not interfere with faster mechanisms of E-I balance modulation, due to the separability of their timescales [ 17 ], it might still be essential to ensure that such mechanisms can be quickly recruited to shape circuit dynamics in behaviorally relevant timescales. To conclude, the framework developed and explored in this work benefits the study of the role of E-I balance, and its different homeostatic mechanisms, regarding the emergence of network dynamics and interactions in large-scale networks. In addition, our extensive study of node bifurcation dynamics allows for a deeper understanding of how nodes react in the network context, with implications for network stability and the emergence of collective patterns. More importantly, our results show that the regulation of firing rates in cortical PY neurons can also serve to control edge-of-bifurcation dynamics in cortical networks. Importantly, the bifurcation between stable and damped oscillations is shaped differently by distinct modes of homeostasis, offering a plausible explanation for the heterogeneity in distance to the bifurcation found across cortical areas, with likely implications for the behavior of large-scale networks. In addition, we show that E-I homeostasis, particularly when involving inhibitory synapses, has an impact on the frequency of gamma oscillations, which can explain the variation of gamma frequencies with the density of inhibitory receptors across cortical areas. Finally, our results evidence the robustness of relying on multiple modes of homeostasis, which can maintain stable firing rates with minimal perturbations to node dynamics, ensuring that gamma oscillations remain within a physiological range. All in all, we assert the pivotal role of E-I balance in cortical dynamics, posing E-I homeostasis as a robust mechanism for the control of not only firing rates but also edge-of-bifurcation dynamics in cortical networks. Supporting information Download figure Open in new tab S1 Fig Fig Activation Functions F ( x ) with different combinations of µ E and σ E = Kµ E . Note that the value of F (0) is always the same as long as σ E = Kµ E Download figure Open in new tab Download figure Open in new tab S2 Fig Fig Model parameters and dynamics for multiple values of ρ Download figure Open in new tab Download figure Open in new tab S3 Fig Fig Ranges of Variation of Model Parameters Under Different Modes of Homeostasis Download figure Open in new tab S4 Fig Fig Oscillation Frequency as Function of c EI for Different Modes of Homeostasis We present the relationship between the intrinsic frequency of oscillation and the strength of inhibition c EI in models with c EI (Top Left), G E + c EI (Top Right), G E + c EI + µ E (Bottom Left) and G E + c EI + µ E + σ E (Bottom Right) homeostasis. In addition, for each, we plot the relationship in models with different target firing rates ( ρ ), ranging from 0.1 to 0.14. S1 Appendix Derivation of Nullclines and Fixed Point Expression The expression for the r E nullcline can be obtained from Eqs 1 and 2 through the following chain of operations: Then, we isolate the r I term on one side of the equation to obtain the nullcline, expressing r I as a function of r E and the model parameters: where log denotes the natural logarithm and exp the exponential function. Similarly, we can derive the nullcline of r I by equating to 0, yielding: S2 Appendix Homeostasis of Excitatory Synapses Unto Inhibitory Neurons Although the evidence is scarce, some studies indicate that fast-spiking interneurons may also scale the strength of their excitatory synapses [ 43 ] ( Fig 2 ) and that this process depends mainly on pre-synaptic activity, with the function of maintaining PY firing rates [ 59 ]. This type of homeostasis can be implemented by scaling the parameter c IE , which quantifies the strength of the excitatory inputs to the inhibitory neural mass. By reordering Eq (5), we obtain: The most notable difference of c IE homeostasis relative to the two other types of synaptic scaling is the lack of robustness to higher levels of external input ( Figure 11 ). More specifically, if I ext ≳ 2.65, the system is not able to maintain any desired fixed point through homeostasis of excitatory synapses unto the inhibitory population (see blank space in Figure 11 ). This finding can be trivially explained given the dynamics of our model. For this type of homeostasis, when excitatory activity is too high, the model increases the excitation of inhibitory populations so that they inhibit the excitatory population more. However, because of the sigmoid non-linearity of the model, there is a cap on the activity of the inhibitory population at r I = 1. Therefore, regardless of how much the inhibitory population is excited, the maximum inhibitory input it can provide to the excitatory neural mass is equal to c EI × 1. For this reason, after the point where homeostasis requires the saturation of the inhibitory population to maintain , the system can no longer respond to increases in I ext and, thus, it is not able to maintain homeostasis of excitatory activity. For this reason, we do not consider this mode of homeostasis in our further exploration of E-I balance in the WC model. Download figure Open in new tab Fig 11. Behavior of the Wilson-Cowan Model Under Homeostasis of c IE . (Left) Homeostatic value of c IE and (Right) corresponding system behavior for different combinations of and I ext . S3 Appendix Computing steady-state parameters in systems with multiple modes of homeostasis To derive an analytical method to compute the steady-state parameters of models under multiple modes of homeostasis we start by defining equations describing the dynamics of homeostasis. Starting with plasticity of G E and c EI , we can write their dynamics as: Shortly, if the firing rate of the excitatory neural mass r E is larger than the target firing rate ρ, G E will decrease and c EI will increase, bringing the firing rate down due to reduced excitation and enhanced inhibition. While this implementation of c EI homeostasis differs from the common approach in large-scale models [ 36 , 37 , 39 , 85 ], in S4 Appendix we argue that it is more sensible for this version of the Wilson-Cowan model. Given these equations, one can write the value of G E at time t as follows: where is the initial value of G E at t = 0. Similary, we can write c EI ( t ) as: Considering that the term inside the integral is the same in Eqs 25 and 26, we can write: Therefore, in simpler terms, the magnitude of variation in both parameters between any two points in time is the same, provided that the timescales of homeostatic plasticity ( τ homeo ) are also the same. Similarly, for modes of plasticity involving the modulation of intrinsic excitability, and considering that, if firing rates are higher than ρ , the firing threshold µ E should increase, we can write the dynamics of µ E as: leading to: Then, we can write the following equality, relating µ E to G E and c EI : However, up until now, we have assumed that all forms of homeostasis operate in the same timescales. That said, it is still possible to solve the equations when the timescales are different. Let us assume that the timescales of c EI and µ E can be related to the timescale of G E ( τ homeo ) through the ratios R 1 and R 2 , respectively. In that case, we write: which can be rewritten as: Therefore, we can re-write the previous equality ( Eq. 30 ) as: This can be then employed to derive the new equations determining the steady-state values of the model parameters under multiple modes of homeostasis considering that the timescales of each mode can be different. That considered, for the plasticity of G E and c EI , we have: Similarly, for plasticity of G E , c EI and µ E , we obtain: Finally, for plasticity of G E , c EI , µ E and σ E , we reach the following expression: Where . In Figure 12 , we demonstrate that our analytical method can accurately predict the steady-state values of model parameters with different initial values and time constants. Download figure Open in new tab Fig 12. Evolution of model parameters under different modes of homeostasis. Solid lines represent the model parameters at any given point in time. Conversely, dashed lines represent the predicted values using Eqs 34 , 35 and 36. For all simulations, we used I ext = 2 and ρ = 0.1, and the equations were integrated using the Euler method with a time step of 0.2 ms. S4 Appendix Dynamics of c EI Homeostasis The most common implementation of homeostatic plasticity of inhibition relies on the fact that this plasticity can be disabled by abolishing firing in both the post-synaptic excitatory and pre-synaptic inhibitory neurons [ 17 ]. For this reason, [ 14 ] implemented this mode of homeostasis for feedforward inhibitory synapses in the context of networks of rate-based neurons with the following formulation: where ρ represents the target average firing rate. This implementation has been later translated into the context of large-scale modeling, particularly in approaches relying on the Wilson-Cowan model or similar neural-mass models [ 36 , 37 , 39 – 41 , 85 ]. However, here we argue that this is may not the most appropriate approach for the current implementation of the Wilson-Cowan model. Let us start by considering that the time constant of homeostatic plasticity is normally orders of magnitude slower than the timescale of population activity, ensuring the separation of timescales. Therefore, for the analysis of homeostasis, it can be considered that the dynamics of Equation 37 dominate the system, while r E and r I quickly reach their respective steady-states, where r I = F I ( c IE r E ). Therefore, we can write Equation 37 as: Now, let us consider that we have two different initial systems for which the steady-state firing rate and , with d > 0. Therefore, we have: Given that F I is sigmoid and, thus, a monotonously increasing function, we have that: which means that: Therefore, in two systems where the firing rate is equally distant from ρ , the one where the firing rate is lower than ρ will necessarily converge more slowly toward the target than the one where r E is above ρ . For this reason, under this implementation of homeostasis, balanced WC nodes recover more slowly from decreases in firing rates than from increases. However, there is no empirical evidence for this behavior, with results suggesting that homeostasis of inhibition operates in similar timescales regardless of the direction of perturbation [ 17 , 20 , 21 , 58 , 86 ]. We illustrate the influence of the r I term on the derivative of c EI in Figure 13 and compare it with the case where . Here, the aforementioned behavior can be observed more clearly, since the magnitude of the derivative of c EI is considerably lower when r E < ρ = 0.1. In addition, comparison with r E − ρ allows for the visualization of a further issue with this implementation. Since the sigmoid non-linearity is bound between 0 and 1, r I ( r E − ρ ) < ( r E − ρ ), this form of plasticity will always be slower than modes that are only dependent on the difference between r E and ρ , such as synaptic scaling of excitation [ 17 , 19 , 43 ], even when considering equal timescales. This effect is more pronounced for values of r E lower than ρ , since r I approaches 0, further slowing down the homeostasis of c EI . Download figure Open in new tab Fig 13. as a function of r E . We present two cases where (solid blue line) or (dashed black line). For illustration purposes, we consider that τ homeo = 1 s and ρ = 0.1. For these reasons, we consider that, in the current implementation of the WC model, the most reasonable expression for c EI homeostasis is: as opposed to the most common implementation [ 36 , 37 , 39 – 41 , 85 ]. Furthermore, we argue that it may still be in line with literature. More specifically, [ 17 ] states that “in contrast to synaptic scaling of excitatory synapses, homeostatic regulation of inhibition is a non-cell-autonomous process that either requires changes in both pre- and postsynaptic activity simultaneously or is triggered by global changes in network activity . ” . In our case, since we are working with mean-field models, which abstract the aggregate network activity into a single variable r E or r I , it can be considered that r E is sufficient to track network activity. Not only that, but, given the separation in timescales between neural dynamics and homeostasis, r I can be understood as F I ( c IE r E ), which means that perturbing r E will lead to a perturbation of r I in the same direction, scaled by F I . For these reasons, we consider that Equation 42 represents a reasonable implementation of homeostasis of inhibition in the Wilson-Cowan model, and possibly other mean-field approaches, particularly when activity is constrained by a sigmoid non-linearity. S5 Appendix Derivation of Jacobian Matrix for the Wilson-Cowan Model Let us consider that the equations of the Wilson-Cowan system can be written as: where: Then, the first-order part of the Taylor series of the Wilson-Cowan around a given point can be written as: Given that a is, by definition, a fixed point of the system f ( a ) = g ( a ) = 0. Therefore, we have a linear approximation of the system around a , from which we obtain the Jacobian matrix ( J ) containing the partial derivatives of f ( x ) and g ( x ). The partial derivatives of each function can be expressed as a function of the model parameters as follows. Furthermore, since the derivative of the canonical sigmoid function is F ( x )(1 − F ( x )), we can use the chain rule to derive the second part of as follows, where x = G E c EE r E − c EI r I + G E I ext : Similarly, we can derive the other components of the Jacobian matrix as: S6 Appendix Analytical Derivation of Oscillation Frequency Around a Fixed Point In a given fixed point of any 2-dimensional dynamical system, if the eigenvalues of the Jacobian have an imaginary component (i.e. Tr ( J ) 2 − 4 Det ( J ) 0 then the solution to the eigenvalue problem of the system yields: which corresponds to the following solution of the system, where v ± represents the eigenvectors of the Jacobian: Therefore, given that the e iωt term represents a sinusoidal wave in time, we can estimate the frequency of oscillation around the fixed point from the Jacobian matrix, given that ω = 2 πf , as: This approximation, however, is only valid around the fixed point. For example, in systems with a repelling fixed point, where dynamics settle on a limit cycle trajectory, the frequency of oscillation will decrease as the system orbits away from the fixed point ( Fig 14 ). Download figure Open in new tab Fig 14. Simulated r E and respective frequency of oscillation under c EI homeostasis for different combinations of I ext and fixed point r E . The black dashed line represents the analytical prediction for the frequency of oscillation around the fixed point Footnotes https://gitlab.com/francpsantos/wc_node_homeostasis References 1. ↵ Douglas RJ , Martin KAC , Whitteridge D. A Canonical Microcircuit for Neocortex . Neural Computation . 1989 ; 1 ( 4 ): 480 – 488 . doi: 10.1162/neco.1989.1.4.480 . OpenUrl CrossRef 2. Martin KA . The Wellcome Prize lecture. From single cells to simple circuits in the cerebral cortex . Quarterly Journal of Experimental Physiology (Cambridge, England) . 1988 ; 73 ( 5 ): 637 – 702 . doi: 10.1113/expphysiol.1988.sp003190 . OpenUrl CrossRef PubMed 3. ↵ Feldmeyer D , Qi G , Emmenegger V , Staiger JF . Inhibitory interneurons and their circuit motifs in the many layers of the barrel cortex . Neuroscience . 2018 ; 368 : 132 – 151 . doi: 10.1016/j.neuroscience.2017.05.027 . OpenUrl CrossRef PubMed 4. ↵ Douglas RJ , Martin KAC . Neuronal circuits of the neocortex . Annual Review of Neuroscience . 2004 ; 27 : 419 – 451 . doi: 10.1146/annurev.neuro.27.070203.144152 . OpenUrl CrossRef PubMed Web of Science 5. Benucci A , Verschure PFMJ , König P. High-order events in cortical networks: A lower bound . Physical Review E . 2004 ; 70 ( 5 ): 051909 . doi: 10.1103/PhysRevE.70.051909 . OpenUrl CrossRef 6. Roerig B , Chen B. Relationships of local inhibitory and excitatory circuits to orientation preference maps in ferret visual cortex . Cerebral Cortex (New York, NY: 1991) . 2002 ; 12 ( 2 ): 187 – 198 . doi: 10.1093/cercor/12.2.187 . OpenUrl CrossRef PubMed Web of Science 7. ↵ Wang XJ . Macroscopic gradients of synaptic excitation and inhibition in the neocortex . Nature Reviews Neuroscience . 2020 ; 21 ( 3 ): 169 – 178 . doi: 10.1038/s41583-020-0262-x . OpenUrl CrossRef 8. ↵ Barral J , D Reyes A . Synaptic scaling rule preserves excitatory-inhibitory balance and salient neuronal network dynamics . Nature Neuroscience . 2016 ; 19 ( 12 ): 1690 – 1696 . doi: 10.1038/nn.4415 . OpenUrl CrossRef PubMed 9. ↵ Douglas RJ , Martin KA . A functional microcircuit for cat visual cortex . The Journal of Physiology . 1991 ; 440 ( 1 ): 735 – 769 . doi: 10.1113/jphysiol.1991.sp018733 . OpenUrl CrossRef PubMed Web of Science 10. Wehr M , Zador AM . Balanced inhibition underlies tuning and sharpens spike timing in auditory cortex . Nature . 2003 ; 426 ( 6965 ): 442 – 446 . doi: 10.1038/nature02116 . OpenUrl CrossRef PubMed Web of Science 11. Okun M , Lampl I. Instantaneous correlation of excitation and inhibition during ongoing and sensory-evoked activities . Nature Neuroscience . 2008 ; 11 ( 5 ): 535 – 537 . doi: 10.1038/nn.2105 . OpenUrl CrossRef PubMed Web of Science 12. ↵ Xue M , Atallah BV , Scanziani M. Equalizing excitation-inhibition ratios across visual cortical neurons . Nature . 2014 ; 511 ( 7511 ): 596 – 600 . doi: 10.1038/nature13321 . OpenUrl CrossRef PubMed Web of Science 13. ↵ van Vreeswijk C , Sompolinsky H. Chaos in neuronal networks with balanced excitatory and inhibitory activity . Science (New York, NY) . 1996 ; 274 ( 5293 ): 1724 – 1726 . doi: 10.1126/science.274.5293.1724 . OpenUrl Abstract / FREE Full Text 14. ↵ Vogels TP , Sprekeler H , Zenke F , Clopath C , Gerstner W. Inhibitory Plasticity Balances Excitation and Inhibition in Sensory Pathways and Memory Networks . Science . 2011 ; 334 ( 6062 ): 1569 – 1573 . doi: 10.1126/science.1211095 . OpenUrl Abstract / FREE Full Text 15. Beggs J , Timme N. Being Critical of Criticality in the Brain . Frontiers in Physiology . 2012 ; 3 . 16. ↵ Sadeh S , Clopath C. Excitatory-inhibitory balance modulates the formation and dynamics of neuronal assemblies in cortical networks . Science Advances . 2021 ; 7 ( 45 ): eabg8411 . doi: 10.1126/sciadv.abg8411 . OpenUrl CrossRef PubMed 17. ↵ Turrigiano G. Too many cooks? Intrinsic and synaptic homeostatic mechanisms in cortical circuit refinement . Annual Review of Neuroscience . 2011 ; 34 : 89 – 103 . doi: 10.1146/annurev-neuro-060909-153238 . OpenUrl CrossRef PubMed Web of Science 18. ↵ Wen W , Turrigiano GG . Keeping Your Brain in Balance: Homeostatic Regulation of Network Function . Annual Review of Neuroscience . 2024 ; doi: 10.1146/annurev-neuro-092523-110001 . OpenUrl CrossRef 19. ↵ Turrigiano GG , Leslie KR , Desai NS , Rutherford LC , Nelson SB . sActivity-dependent scaling of quantal amplitude in neocortical neurons . Nature . 1998 ; 391 ( 6670 ): 892 – 896 . doi: 10.1038/36103 . OpenUrl CrossRef PubMed Web of Science 20. ↵ Maffei A , Turrigiano GG . Multiple Modes of Network Homeostasis in Visual Cortical Layer 2/3 . Journal of Neuroscience . 2008 ; 28 ( 17 ): 4377 – 4384 . doi: 10.1523/JNEUROSCI.5298-07.2008 . OpenUrl Abstract / FREE Full Text 21. ↵ Maffei A , Nelson SB , Turrigiano GG . Selective reconfiguration of layer 4 visual cortical circuitry by visual deprivation . Nature Neuroscience . 2004 ; 7 ( 12 ): 1353 – 1359 . doi: 10.1038/nn1351 . OpenUrl CrossRef PubMed Web of Science 22. ↵ Keck T , Keller GB , Jacobsen RI , Eysel UT , Bonhoeffer T , Hübener M. Synaptic scaling and homeostatic plasticity in the mouse visual cortex in vivo . Neuron . 2013 ; 80 ( 2 ): 327 – 334 . doi: 10.1016/j.neuron.2013.08.018 . OpenUrl CrossRef PubMed Web of Science 23. ↵ Desai NS , Rutherford LC , Turrigiano GG . Plasticity in the intrinsic excitability of cortical pyramidal neurons . Nature Neuroscience . 1999 ; 2 ( 6 ): 515 – 520 . doi: 10.1038/9165 . OpenUrl CrossRef PubMed Web of Science 24. ↵ Nataraj K , Le Roux N , Nahmani M , Lefort S , Turrigiano G. Visual deprivation suppresses L5 pyramidal neuron excitability by preventing the induction of intrinsic plasticity . Neuron . 2010 ; 68 ( 4 ): 750 – 762 . doi: 10.1016/j.neuron.2010.09.033 . OpenUrl CrossRef PubMed Web of Science 25. ↵ Wen W , Turrigiano GG . Developmental Regulation of Homeostatic Plasticity in Mouse Primary Visual Cortex . The Journal of Neuroscience . 2021 ; 41 ( 48 ): 9891 – 9905 . doi: 10.1523/JNEUROSCI.1200-21.2021 . OpenUrl Abstract / FREE Full Text 26. ↵ Ma Z , Turrigiano GG , Wessel R , Hengen KB . Cortical Circuit Dynamics Are Homeostatically Tuned to Criticality In Vivo . Neuron . 2019 ; 104 ( 4 ): 655 – 664.e4 . doi: 10.1016/j.neuron.2019.08.031 . OpenUrl CrossRef PubMed 27. ↵ Freeman WJ . A field-theoretic approach to understanding scale-free neocortical dynamics . Biological Cybernetics . 2005 ; 92 ( 6 ): 350 – 359 . doi: 10.1007/s00422-005-0563-1 . OpenUrl CrossRef PubMed Web of Science 28. ↵ Poil SS , van Ooyen A , Linkenkaer-Hansen K. Avalanche dynamics of human brain oscillations: Relation to critical branching processes and temporal correlations . Human Brain Mapping . 2008 ; 29 : 770 – 777 . doi: 10.1002/hbm.20590 . OpenUrl CrossRef PubMed Web of Science 29. ↵ Nanda A , Johnson GW , Mu Y , Ahrens MB , Chang C , Englot DJ , et al. Time-resolved correlation of distributed brain activity tracks E-I balance and accounts for diverse scale-free phenomena . Cell Reports . 2023 ; 42 ( 4 ): 112254 . doi: 10.1016/j.celrep.2023.112254 . OpenUrl CrossRef 30. ↵ Freeman WJ , Holmes MD . Metastability, instability, and state transition in neocortex . Neural Networks . 2005 ; 18 ( 5 ): 497 – 504 . doi: 10.1016/j.neunet.2005.06.014 . OpenUrl CrossRef PubMed Web of Science 31. ↵ Cabral J , Luckhoo H , Woolrich M , Joensson M , Mohseni H , Baker A , et al. Exploring mechanisms of spontaneous functional connectivity in MEG: how delayed network interactions lead to structured amplitude envelopes of band-pass filtered oscillations . NeuroImage . 2014 ; 90 : 423 – 435 . doi: 10.1016/j.neuroimage.2013.11.047 . OpenUrl CrossRef PubMed Web of Science 32. ↵ Deco G , Kringelbach ML , Jirsa VK , Ritter P. The dynamics of resting fluctuations in the brain: metastability and its dynamical cortical core . Scientific Reports . 2017 ; 7 ( 1 ): 3095 . doi: 10.1038/s41598-017-03073-5 . OpenUrl CrossRef PubMed 33. Lintas A , Rovetta S , Verschure PFMJ , Villa AEP Arsiwalla XD , Verschure P. Why the Brain Might Operate Near the Edge of Criticality . In: Lintas A , Rovetta S , Verschure PFMJ , Villa AEP , editors. Artificial Neural Networks and Machine Learning - ICANN 2017 . Cham : Springer International Publishing ; 2017 . p. 326 – 333 . 34. ↵ Turkheimer FE , Leech R , Expert P , Lord LD , Vernon AC . The brain’s code and its canonical computational motifs . From sensory cortex to the default mode network: A multi-scale model of brain function in health and disease. Neuroscience and Biobehavioral Reviews . 2015 ; 55 : 211 – 222 . doi: 10.1016/j.neubiorev.2015.04.014 . OpenUrl CrossRef PubMed 35. ↵ Harris B , Gollo LL , Fulcher BD . Tracking the distance to criticality in systems with unknown noise ; 2024 . Available from: http://arxiv.org/abs/2310.14791 . 36. ↵ Hellyer PJ , Jachs B , Clopath C , Leech R. Local inhibitory plasticity tunes macroscopic brain dynamics and allows the emergence of functional brain networks . NeuroImage . 2016 ; 124 (Pt A): 85 – 95 . doi: 10.1016/j.neuroimage.2015.08.069 . OpenUrl CrossRef 37. ↵ Naskar A , Vattikonda A , Deco G , Roy D , Banerjee A. Multiscale dynamic mean field (MDMF) model relates resting-state brain dynamics with local cortical excitatory-inhibitory neurotransmitter homeostasis . Network Neuroscience . 2021 ; 5 ( 3 ): 757 – 782 . doi: 10.1162/netna00197 . OpenUrl CrossRef 38. Deco G , Ponce-Alvarez A , Hagmann P , Romani GL , Mantini D , Corbetta M. How Local Excitation-Inhibition Ratio Impacts the Whole Brain Dynamics . Journal of Neuroscience . 2014 ; 34 ( 23 ): 7886 – 7898 . doi: 10.1523/JNEUROSCI.5068-13.2014 . OpenUrl Abstract / FREE Full Text 39. ↵ Abeysuriya RG , Hadida J , Sotiropoulos SN , Jbabdi S , Becker R , Hunt BAE , et al. A biophysical model of dynamic balancing of excitation and inhibition in fast oscillatory large-scale networks . PLoS computational biology . 2018 ; 14 ( 2 ): e1006007 . doi: 10.1371/journal.pcbi.1006007 . OpenUrl CrossRef 40. ↵ Castaldo F , Páscoa dos Santos F , Timms RC , Cabral J , Vohryzek J , Deco G , et al. Multi-modal and multi-model interrogation of large-scale functional brain networks . NeuroImage . 2023 ; 277 : 120236 . doi: 10.1016/j.neuroimage.2023.120236 . OpenUrl CrossRef 41. ↵ Páscoa dos Santos F , Vohryzek J , Verschure PFMJ. Multiscale effects of excitatory-inhibitory homeostasis in lesioned cortical networks: A computational study . PLOS Computational Biology . 2023 ; 19 ( 7 ): e1011279 . doi: 10.1371/journal.pcbi.1011279 . OpenUrl CrossRef 42. ↵ Gainey MA , Feldman DE . Multiple shared mechanisms for homeostatic plasticity in rodent somatosensory and visual cortex . Philosophical Transactions of the Royal Society of London Series B, Biological Sciences . 2017 ; 372 ( 1715 ): 20160157 . doi: 10.1098/rstb.2016.0157 . OpenUrl CrossRef PubMed 43. ↵ Rutherford LC , Nelson SB , Turrigiano GG . BDNF has opposite effects on the quantal amplitude of pyramidal neuron and interneuron excitatory synapses . Neuron . 1998 ; 21 ( 3 ): 521 – 530 . doi: 10.1016/s0896-6273(00)80563-2 . OpenUrl CrossRef PubMed Web of Science 44. ↵ Wilson HR , Cowan JD . Excitatory and inhibitory interactions in localized populations of model neurons . Biophysical Journal . 1972 ; 12 ( 1 ): 1 – 24 . doi: 10.1016/S0006-3495(72)86068-5 . OpenUrl CrossRef PubMed Web of Science 45. ↵ Cowan JD , Neuman J , van Drongelen W. Wilson-Cowan Equations for Neocortical Dynamics . The Journal of Mathematical Neuroscience . 2016 ; 6 ( 1 ): 1 . doi: 10.1186/s13408-015-0034-5 . OpenUrl CrossRef PubMed 46. ↵ Li X , Li Z , Yang W , Wu Z , Wang J. Bidirectionally Regulating Gamma Oscillations in Wilson-Cowan Model by Self-Feedback Loops: A Computational Study . Frontiers in Systems Neuroscience . 2022 ; 16 . 47. ↵ Ibata K , Sun Q , Turrigiano GG . Rapid synaptic scaling induced by changes in postsynaptic firing . Neuron . 2008 ; 57 ( 6 ): 819 – 826 . doi: 10.1016/j.neuron.2008.02.031 . OpenUrl CrossRef PubMed Web of Science 48. ↵ Hengen KB , Lambo ME , Van Hooser SD , Katz DB , Turrigiano GG . Firing rate homeostasis in visual cortex of freely behaving rodents . Neuron . 2013 ; 80 ( 2 ): 335 – 342 . doi: 10.1016/j.neuron.2013.08.038 . OpenUrl CrossRef PubMed Web of Science 49. ↵ Felleman DJ , Van Essen DC . Distributed hierarchical processing in the primate cerebral cortex . Cerebral Cortex (New York, NY: 1991) . 1991 ; 1 ( 1 ): 1 – 47 . doi: 10.1093/cercor/1.1.1-a . OpenUrl CrossRef PubMed Web of Science 50. Keller A. Intrinsic synaptic organization of the motor cortex . Cerebral Cortex (New York, NY: 1991) . 1993 ; 3 ( 5 ): 430 – 441 . doi: 10.1093/cercor/3.5.430 . OpenUrl CrossRef PubMed Web of Science 51. ↵ Tremblay R , Lee S , Rudy B. GABAergic Interneurons in the Neocortex: From Cellular Properties to Circuits . Neuron . 2016 ; 91 ( 2 ): 260 – 292 . doi: 10.1016/j.neuron.2016.06.033 . OpenUrl CrossRef PubMed 52. ↵ Wang XJ , Buzsaki G. Gamma oscillation by synaptic inhibition in a hippocampal interneuronal network model . The Journal of Neuroscience: The Official Journal of the Society for Neuroscience . 1996 ; 16 ( 20 ): 6402 – 6413 . doi: 10.1523/JNEUROSCI.16-20-06402.1996 . OpenUrl Abstract / FREE Full Text 53. ↵ Buzsaki G. Rhythms of the Brain . Oxford University Press ; 2006 . Available from: https://academic.oup.com/book/11166 . 54. Tiesinga P , Sejnowski TJ . Cortical enlightenment: are attentional gamma oscillations driven by ING or PING? Neuron . 2009 ; 63 ( 6 ): 727 – 732 . doi: 10.1016/j.neuron.2009.09.009 . OpenUrl CrossRef PubMed Web of Science 55. ↵ Buzsaki G , Wang XJ . Mechanisms of gamma oscillations . Annual Review of Neuroscience . 2012 ; 35 : 203 – 225 . doi: 10.1146/annurev-neuro-062111-150444 . OpenUrl CrossRef PubMed Web of Science 56. ↵ Deco G , Jirsa V , McIntosh AR , Sporns O , Kötter R. Key role of coupling, delay, and noise in resting brain fluctuations . Proceedings of the National Academy of Sciences . 2009 ; 106 ( 25 ): 10302 – 10307 . doi: 10.1073/pnas.0901831106 . OpenUrl Abstract / FREE Full Text 57. ↵ Cabral J , Castaldo F , Vohryzek J , Litvak V , Bick C , Lambiotte R , et al. Metastable oscillatory modes emerge from synchronization in the brain spacetime connectome . Communications Physics . 2022 ; 5 ( 1 ): 1 – 13 . doi: 10.1038/s42005-022-00950-y . OpenUrl CrossRef 58. ↵ Barnes SJ , Sammons RP , Jacobsen RI , Mackie J , Keller GB , Keck T. Subnetwork-Specific Homeostatic Plasticity in Mouse Visual Cortex In Vivo . Neuron . 2015 ; 86 ( 5 ): 1290 – 1303 . doi: 10.1016/j.neuron.2015.05.010 . OpenUrl CrossRef PubMed 59. ↵ Chang MC , Park JM , Pelkey KA , Grabenstatter HL , Xu D , Linden DJ , et al. Narp regulates homeostatic scaling of excitatory synapses on parvalbumin-expressing interneurons . Nature Neuroscience . 2010 ; 13 ( 9 ): 1090 – 1097 . doi: 10.1038/nn.2621 . OpenUrl CrossRef PubMed Web of Science 60. ↵ Gerstner W , Kistler WM , Naud R , Paninski L. Neuronal Dynamics: From Single Neurons to Networks and Models of Cognition . Cambridge : Cambridge University Press ; 2014 . Available from: https://www.cambridge.org/core/books/neuronal-dynamics/75375090046733765596191E23B2959D . 61. ↵ Izhikevich EM . Dynamical Systems in Neuroscience: The Geometry of Excitability and Bursting . The MIT Press ; 2006 . Available from: https://direct.mit.edu/books/book/2589/Dynamical-Systems-in-NeuroscienceThe-Geometry-of . 62. ↵ Benardo LS . Separate activation of fast and slow inhibitory postsynaptic potentials in rat neocortex in vitro . The Journal of Physiology . 1994 ; 476 ( 2 ): 203 – 215 . doi: 10.1113/jphysiol.1994.sp020124 . OpenUrl CrossRef PubMed Web of Science 63. ↵ Tamás G , Simon AL Anna , Szabadics J. Identified Sources and Targets of Slow Inhibition in the Neocortex . Science . 2003 ; 299 ( 5614 ): 1902 – 1905 . doi: 10.1126/science.1082053 . OpenUrl Abstract / FREE Full Text 64. ↵ Maier A , Adams G , Aura C , Leopold D. Distinct Superficial and Deep Laminar Domains of Activity in the Visual Cortex during Rest and Stimulation . Frontiers in Systems Neuroscience . 2010 ; 4 . 65. van Kerkoerle T , Self MW , Dagnino B , Gariel-Mathis MA , Poort J , van der Togt C , et al. Alpha and gamma oscillations characterize feedback and feedforward processing in monkey visual cortex . Proceedings of the National Academy of Sciences . 2014 ; 111 ( 40 ): 14332 – 14341 . doi: 10.1073/pnas.1402773111 . OpenUrl Abstract / FREE Full Text 66. ↵ Mejias JF , Murray JD , Kennedy H , Wang XJ . Feedforward and feedback frequency-dependent interactions in a large-scale laminar network of the primate cortex . Science Advances . 2016 ; 2 ( 11 ): e1601335 . doi: 10.1126/sciadv.1601335 . OpenUrl FREE Full Text 67. ↵ Sanchez-Todo R , Bastos AM , Lopez-Sola E , Mercadal B , Santarnecchi E , Miller EK , et al. A physical neural mass model framework for the analysis of oscillatory generators from laminar electrophysiological recordings . NeuroImage . 2023 ; 270 : 119938 . doi: 10.1016/j.neuroimage.2023.119938 . OpenUrl CrossRef 68. ↵ Deco G , Jirsa VK . Ongoing Cortical Activity at Rest: Criticality, Multistability, and Ghost Attractors . Journal of Neuroscience . 2012 ; 32 ( 10 ): 3366 – 3375 . doi: 10.1523/JNEUROSCI.2523-11.2012 . OpenUrl Abstract / FREE Full Text 69. ↵ Schirner M , Kong X , Yeo BTT , Deco G , Ritter P. Dynamic primitives of brain network interaction . NeuroImage . 2022 ; 250 : 118928 . doi: 10.1016/j.neuroimage.2022.118928 . OpenUrl CrossRef 70. ↵ Kujala J , Jung J , Bouvard S , Lecaignard F , Lothe A , Bouet R , et al. Gamma oscillations in V1 are correlated with GABA(A) receptor density: A multi-modal MEG and Flumazenil-PET study . Scientific Reports . 2015 ; 5 : 16347 . doi: 10.1038/srep16347 . OpenUrl CrossRef PubMed 71. ↵ Kujala J , Ciumas C , Jung J , Bouvard S , Lecaignard F , Lothe A , et al. GABAergic inhibition shapes behavior and neural dynamics in human visual working memory . Cerebral Cortex . 2024 ; 34 ( 2 ): bhad522 . doi: 10.1093/cercor/bhad522 . OpenUrl CrossRef 72. ↵ Arnulfo G , Wang SH , Myrov V , Toselli B , Hirvonen J , Fato MM , et al. Long-range phase synchronization of high-frequency oscillations in human cortex . Nature Communications . 2020 ; 11 ( 1 ): 5363 . doi: 10.1038/s41467-020-18975-8 . OpenUrl CrossRef PubMed 73. ↵ Harris AZ , Gordon JA . Long-Range Neural Synchrony in Behavior . Annual Review of Neuroscience . 2015 ; 38 (Volume 38, 2015): 171 – 194 . doi: 10.1146/annurev-neuro-071714-034111 . OpenUrl CrossRef PubMed 74. ↵ Sakata S , Harris KD . Laminar Structure of Spontaneous and Sensory-Evoked Population Activity in Auditory Cortex . Neuron . 2009 ; 64 ( 3 ): 404 – 418 . doi: 10.1016/j.neuron.2009.09.020 . OpenUrl CrossRef PubMed Web of Science 75. ↵ Potjans TC , Diesmann M. The cell-type specific cortical microcircuit: relating structure and activity in a full-scale spiking network model . Cerebral Cortex (New York, NY: 1991) . 2014 ; 24 ( 3 ): 785 – 806 . doi: 10.1093/cercor/bhs358 . OpenUrl CrossRef PubMed 76. Bastos AM , Loonis R , Kornblith S , Lundqvist M , Miller EK . Laminar recordings in frontal cortex suggest distinct layers for maintenance and control of working memory . Proceedings of the National Academy of Sciences . 2018 ; 115 ( 5 ): 1117 – 1122 . doi: 10.1073/pnas.1710323115 . OpenUrl Abstract / FREE Full Text 77. Senzai Y , Fernandez-Ruiz A , Buzsáki G. Layer-Specific Physiological Features and Interlaminar Interactions in the Primary Visual Cortex of the Mouse . Neuron . 2019 ; 101 ( 3 ): 500 – 513.e5 . doi: 10.1016/j.neuron.2018.12.009 . OpenUrl CrossRef PubMed 78. ↵ Self MW , van Kerkoerle T , Goebel R , Roelfsema PR . Benchmarking laminar fMRI: Neuronal spiking and synaptic activity during top-down and bottom-up processing in the different layers of cortex . NeuroImage . 2019 ; 197 : 806 – 817 . doi: 10.1016/j.neuroimage.2017.06.045 . OpenUrl CrossRef PubMed 79. ↵ Markov NT , Ercsey-Ravasz M , Van Essen DC , Knoblauch K , Toroczkai Z , Kennedy H. Cortical high-density counterstream architectures . Science (New York, NY) . 2013 ; 342 ( 6158 ): 1238406 . doi: 10.1126/science.1238406 . OpenUrl Abstract / FREE Full Text 80. ↵ Markov NT , Ercsey-Ravasz MM , Ribeiro Gomes AR , Lamy C , Magrou L , Vezoli J , et al. A Weighted and Directed Interareal Connectivity Matrix for Macaque Cerebral Cortex . Cerebral Cortex . 2014 ; 24 ( 1 ): 17 – 36 . doi: 10.1093/cercor/bhs270 . OpenUrl CrossRef PubMed Web of Science 81. ↵ Linkenkaer-Hansen K , Nikouline VV , Palva JM , Ilmoniemi RJ . Long-Range Temporal Correlations and Scaling Behavior in Human Brain Oscillations . Journal of Neuroscience . 2001 ; 21 ( 4 ): 1370 – 1377 . doi: 10.1523/JNEUROSCI.21-04-01370.2001 . OpenUrl Abstract / FREE Full Text 82. ↵ Hardstone R , Poil SS , Schiavone G , Jansen R , Nikulin V , Mansvelder H , et al. Detrended Fluctuation Analysis: A Scale-Free View on Neuronal Oscillations . Frontiers in Physiology . 2012 ; 3 . 83. ↵ Amil AF , Verschure PFMJ. Supercritical dynamics at the edge-of-chaos underlies optimal decision-making . Journal of Physics: Complexity . 2021 ; 2 ( 4 ): 045017 . doi: 10.1088/2632-072X/ac3ad2 . OpenUrl CrossRef 84. ↵ Puigbò JY , Maffei G , Herreros I , Ceresa M , González Ballester MA , Verschure PFMJ. Cholinergic Behavior State-Dependent Mechanisms of Neocortical Gain Control: a Neurocomputational Study . Molecular Neurobiology . 2018 ; 55 ( 1 ): 249 – 257 . doi: 10.1007/s12035-017-0737-6 . OpenUrl CrossRef 85. ↵ Hadida J , Sotiropoulos SN , Abeysuriya RG , Woolrich MW , Jbabdi S. Bayesian Optimisation of Large-Scale Biophysical Networks . NeuroImage . 2018 ; 174 : 219 – 236 . doi: 10.1016/j.neuroimage.2018.02.063 . OpenUrl CrossRef 86. ↵ Dorrn AL , Yuan K , Barker AJ , Schreiner CE , Froemke RC . Developmental sensory experience balances cortical excitation and inhibition . Nature . 2010 ; 465 ( 7300 ): 932 – 936 . doi: 10.1038/nature09119 . OpenUrl CrossRef PubMed Web of Science View the discussion thread. Back to top Previous Next Posted June 06, 2024. Download PDF Data/Code Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following Excitatory-Inhibitory Homeostasis and Bifurcation Control in the Wilson-Cowan Model of Cortical Dynamics Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share Excitatory-Inhibitory Homeostasis and Bifurcation Control in the Wilson-Cowan Model of Cortical Dynamics Francisco Páscoa dos Santos , Paul FMJ Verschure bioRxiv 2024.06.05.597312; doi: https://doi.org/10.1101/2024.06.05.597312 Share This Article: Copy Citation Tools Excitatory-Inhibitory Homeostasis and Bifurcation Control in the Wilson-Cowan Model of Cortical Dynamics Francisco Páscoa dos Santos , Paul FMJ Verschure bioRxiv 2024.06.05.597312; doi: https://doi.org/10.1101/2024.06.05.597312 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Neuroscience Subject Areas All Articles Animal Behavior and Cognition (7643) Biochemistry (17717) Bioengineering (13910) Bioinformatics (42018) Biophysics (21480) Cancer Biology (18629) Cell Biology (25537) Clinical Trials (138) Developmental Biology (13392) Ecology (19935) Epidemiology (2067) Evolutionary Biology (24356) Genetics (15617) Genomics (22531) Immunology (17755) Microbiology (40438) Molecular Biology (17200) Neuroscience (88706) Paleontology (667) Pathology (2840) Pharmacology and Toxicology (4832) Physiology (7657) Plant Biology (15171) Scientific Communication and Education (2046) Synthetic Biology (4304) Systems Biology (9828) Zoology (2272)