| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |
|
|
|
|||||||
|
|||||||||
1Toronto Western Research Institute, University Health Network; 2Departments of Medicine (Neurology), Physiology, and 3Institute of Biomaterials and Biomedical Engineering, University of Toronto, Toronto; and 4Department of Applied Mathematics, University of Waterloo, Waterloo, Ontario, Canada
Submitted 18 June 2004; accepted in final form 15 November 2004
 |
ABSTRACT |
|---|
 TOP ABSTRACT INTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
|
INTRODUCTION |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
10–20% of the total neuronal population. They are known to be very diverse in terms of their biochemical content, morphology, electrophysiological characteristics, and neuromodulator sensitivities (McBain and Fisahn 2001
). Importantly, networks of these inhibitory cells have been found to be responsible for the generation and control of rhythmic brain activities, molding the appropriate temporal output of pyramidal cells (Buzsáki and Chrobak 1995; Freund and Buzsáki 1996). For example, oscillatory activity in the hippocampus in the theta (8–12 Hz) and gamma (20–80 Hz) frequency bands occur during memory consolidation and spatial navigation (Bragin et al. 1995; Buzsáki 2002). Several experimental studies indicate that population (network) rhythms arise as a result of coherent activities in interneurons (e.g., Traub et al. 1999; Whittington et al. 1995; Wu et al. 2002). Extensive work on interneurons in recent years is giving rise to well-defined characteristics of interneurons with potential functional significance. For example, parvalbumin-containing basket cells seem to make up a relatively unmodifiable inhibitory network population responsible for generating gamma and theta rhythms (Freund 2003). To understand the generation and control of population rhythms, we need to know what are the important controlling factor(s) (i.e., biophysical constraints) or cellular mechanism(s) underlying interneuron coherence.
It was over a decade ago that modeling studies showed that it is possible to obtain synchronous output from purely inhibitory networks (Wang and Rinzel 1992). Since then, several modeling and theoretical studies of inhibitory networks have been performed (e.g., Skinner et al. 1994; vanVreeswijk et al. 1994; Wang and Rinzel 1993). The inclusion of heterogeneous inputs to inhibitory networks has been considered in later studies (Bartos et al. 2001, 2002; Maex and De Schutter 2003; Neltner et al. 2000; Tiesinga and José 2000; Wang and Buzsáki 1996; White et al. 1998), and from them, it is clear that heterogeneity is a significant factor that strongly affects the ability of inhibitory networks to synchronize. However, an exploration of how much heterogeneity can be tolerated by inhibitory networks under a range of synaptic conductance and time constant conditions has not been done. This is an important consideration since recent work by Bartos and colleagues (2001, 2002), in which inhibitory synaptic characteristics were directly measured, showed that the decay time constants were smaller (<5 ms) than what is typically used in inhibitory network models. Re-doing simulations similar to Wang and Buzsáki (1996), they found that their coherent network oscillations exhibited higher frequencies and were more robust to heterogeneities. Because the modeling studies (Bartos et al. 2001, 2002; Wang and Buzsáki 1996) were performed using large (100-cell) networks along with the consideration of other issues (such as amount of connectivity and electrical coupling), it is difficult to determine the underlying mechanisms that give rise to the different network characteristics (frequency and coherence). However, as shown by White and colleagues (1998), using smaller (two-cell) networks is helpful in understanding larger network dynamics. Specifically, they identified two different ways in which coherence could be lost. This occurred in two different regimes, termed phasic and tonic, as identified by different values of the ratio of time constant and network period. Their work shed light on why optimal synchronization at particular frequencies might be obtained in the larger networks. However, White et al.'s (1998) studies were restricted to regions of mild heterogeneity (<5% difference in intrinsic frequencies).
In this paper, we explore how the amount of input heterogeneity to two-cell inhibitory networks affects their dynamics. We find that heterogeneity can be used to determine an optimal set of model parameters for coherent rhythms. Further, we show that coherence properties are preserved in larger (10-cell) networks suggesting that a strategy of "embedding" small network dynamics in larger networks might be a useful way to understand the contribution of biophysically derived parameters to population dynamics in large networks.
|
METHODS |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
to represent the intrinsic properties of a hippocampal interneuron cell. The equations for each cell are given by
 | (1) |
 | (2) |
 | (3) |
is the steady-state activation of the sodium current and is given by: m = 
m/(m + 
m) where m(V) = –0.1(V + 35)/{exp[–0.1(V + 35)] – 1}, m(V) = 4exp[–(V + 60)/18]. h(V) = 0.07exp[–(V + 58)/20], h(V) = 1/{exp[–0.1(V + 28)] + 1}, n(V) = –0.01exp{–(V + 34)/exp[–0.1(V + 34)] – 1}, n(V) = 0.125exp[–(V + 44)/80], 
= 5. Maximal sodium, gNa, potassium, gK, and leak, gL, conductances are: 35, 9, and 0.1 mS/cm2, respectively. Reversal potentials, VNa, VK, VL, are 55, –90, and –65 mV, respectively, and the capacitance, C, is 1 µF/cm2.
We consider neuronal networks formed due to coupling via inhibitory synapses that are described by first order kinetics. The inhibitory synaptic current, Isyn, which would be added to the current balance equation of the postsynaptic cell, is given by
 | (4) |
 | (5) |
 | (6) |
= 6.25 ms–1 is the rate constant of the synaptic activation taken from (Bartos et al. 2001), and 
syn is the synaptic decay time constant. syn values as low as 1 ms were measured between hippocampal basket cells (Bartos et al. 2002), and so we explore a syn range of 1–10 ms in our simulations, with a 1-ms resolution. gsyn values of 0.05, 0.15, 0.25, 0.35, and 0.5 mS/cm2 are used because they encompass physiological values as measured in (Bartos et al. 2001, 2002). Iapp represents the applied or external drive to the cell, and we use this parameter to introduce heterogeneity into the system. We consider two-cell networks that are reciprocally coupled and 10-cell all-to-all coupled networks.
For the two-cell networks, the external drive to cell 1 or 2 is
 |
 |
. We define the percent heterogeneity, %Het, as
 | (7) |
syn, Iµ, and are examined. For the ten-cell networks, the Iapps are randomly chosen in the interval (–%Het,+%Het). All-to-all coupled networks are simulated and gsyn values are normalized so that each presynaptic cell "delivers" a maximal inhibition of gsyn/9.
A measure based on (White et al. 1998) is used to determine network "coherence" or the level of synchrony present in the 10-cell networks. This measure approximates the amount of overlap that exists between two spiking cells where the amount of overlap used is 20% of the period of the faster firing cell. Briefly, the spike trains of each pair of neurons are approximated by a series of square pulses of unit height and fixed width of 20% of the period of the faster firing cell. The shared area of the square pulses from each train that overlap in time is calculated for a given duration. This is equivalent to taking the cross-correlation at zero time lag. The average between all pairs is taken, and the last 2 s of the 10-cell simulations (which were 5 s in total duration) are used as the duration to calculate the coherence measures presented in the paper. Gaussian noise is included in some of the 10-cell simulations by adding it to the current balance equation of each cell (i.e., Eq. 1). Zero mean and SDs from half to twice the amount of heterogeneity in the system is used.
Network simulations are performed using a modified version of our in-house software, NNET (Murray 2004; Skinner and Liu 2003), which integrates the system of differential equations using CVODE (Cohen and Hindmarsh 1996). Unless otherwise stated, initial conditions (ICs) for the two-cell network simulations are: V1 = –58.7249, V2 = –55.0456, h1 = h2 = 0.9379, n1 = n2 = 0.1224, s1 = s2 = 0.1386. In some cases, these values were varied, but this was not done systematically.
A dynamical systems approach for studying nonlinear ordinary differential equations, such as Eqs. 1–3, involves bifurcation analyses, in which the dependence of solutions to the differential equations on various parameters is examined. Solutions include steady states and oscillations. A well-known feature of nonlinear systems is its ability to express multiple stable solutions, i.e., multistability. For a fixed set of parameters, one can determine equilibrium points for the system. For example, V values that do not change with time, dV/dt = 0 (see Eq. 1) and so on. It is the stability of these equilibrium points that determines the particular solution(s) that the system expresses (steady states, oscillations etc.). A dynamical system is said to undergo a bifurcation if the qualitative dynamics of the system are different above and below a particular value of a parameter. The value where the change occurs is called a bifurcation value. For example, the WB model (Wang and Buzsáki 1996) has an Iapp bifurcation value close to 0.2 µA/cm2 in that below this Iapp value, the model system expresses a steady-state voltage output, and above this value, the model system produces an oscillatory output, i.e., repetitive firing. Bifurcations are considered local or global depending on details regarding the number of equilibrium points and their stability characteristics. Common local bifurcations are Hopf bifurcations and saddle-node bifurcations, and common global ones are homoclinic bifurcations and saddle node of limit cycles bifurcations. Knowing the specific type of bifurcation that a system expresses as particular parameters are varied can allow one to predict the presence of multistable patterns. For example, a subcritical Hopf bifurcation arising in a bifurcation analysis using a particular parameter means that bistability exists for a range of values of the particular parameter. So, for a given parameter value in that range, the system expresses two stable patterns, say oscillations and steady states. A bifurcation diagram is a way to illustrate the stability of equilibrium points and periodic solutions as a particular parameter is varied (e.g., see Fig. 6). Because it is hardly ever possible to perform these calculations analytically, numerical techniques must be employed. A numerical continuation refers to using numerical approaches to find equilibrium points and periodic solutions and to examine how their stability changes as parameters are varied. Although challenging to set up and run, numerical continuations have the advantage over numerical simulations of easily finding multistability and finding highly accurate bifurcation values. Numerical solutions from simulations can only produce stable solutions. If multistability is present, the solution produced depends on the initial conditions.
|
) in the XPPAUT software package (Ermentrout 2002). Depending on the particular set of parameters, step sizes are adjusted to allow continuation of solutions along the various branches. However, the maximum stepsize used to obtain the results shown in Tables 3 and 4 is 0.001 (thus determining the resolution).
|
|
|
RESULTS |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
In homogeneous, two-cell inhibitory networks, in-phase (synchronous) and anti-phase oscillations are possible depending on the particular parameters in operation (Lewis and Rinzel 2003; Skinner et al. 1994; Wang and Rinzel 1992). When mild heterogeneity is introduced into the two-cell system, harmonic locking and asynchronous states are also obtained (White et al. 1998). For the ranges of parameters explored (see METHODS), we obtained five distinct patterns in our two-cell network simulations. These observed patterns are 1) near-synchronous, 2) near-antiphase, 3) varied phase-locking with small phase lags (i.e., 1-to-1, but period varies, and phase lags are <10%), 4) suppressed states, and 5) harmonic locking. Examples of each of these are shown in Fig. 1, A–E. In addition, we obtained dynamic states that did not express any clear pattern, referred to as asynchronous behavior. An example is illustrated in Fig. 1F. We considered a pattern to be asynchronous if we did not observe any of patterns after at least one second of simulation time when transients should no longer be apparent (given the time constants involved in the model). We examined 
25 spikes after transients before concluding that there was no clear pattern. Multistable patterns were uncovered in some parameter ranges when different initial conditions were explored. We also obtained patterns that are mixtures of the above patterns. Because these mixture variants did not occur in any obvious fashion, we did not examine them further here. Assuming that closely aligned firing of cells in inhibitory networks are necessary for generating population rhythms, one expects that patterns 1 (Fig. 1A) and 3 (Fig. 1C) are the most relevant.
|
syn and %Het are varied for given gsyn and Iµ values. As gsyn was increased, the amount of the %Het-syn plane that is represented by suppressed states increased. For smaller gsyn values, harmonic locking states dominated the plane. We also found that for smaller Iµ values more suppressed states occurred over the %Het-syn plane (not shown). The expression of suppressed states in the plane can be somewhat understood by results obtained in previous modeling studies by White et al. (1998). Although they did not specifically discuss how their network dynamics changed with variations in the amount of heterogeneity in the system, they did find that there are two different mechanisms by which networks could lose coherence in the presence of mild heterogeneity (<5% difference in intrinsic frequencies). For one of the mechanisms, in their so-called phasic regime, syn/T has smaller values relative to the tonic regime (where T refers to the network period). Coherence occurs in this regime, and is lost via suppressed states for smaller Iapp and larger gsyn values (see Fig. 2, bottom).
|
Because we are assuming that closely aligned firing of the cells in the inhibitory networks are important for generating population rhythms, we define the robustness of a two-cell network as its ability to express near-synchronous network oscillations (i.e., "coherence") in the face of heterogeneity of external inputs received (i.e., 
0). A stronger robustness means that the network can express coherence for larger amounts of heterogeneity. In Fig. 2 (top), we see that the robustness of the two-cell system increased as syn was increased up to just under 7%Het for syn =10 ms. We can consider a measure of the robustness of an inhibitory network system as the maximal %Het expressed for a given set of parameters. Therefore in Fig. 2 (top), the robustness or maximal %Het for syn = 1 is 4%. Interestingly, our simulations indicate that there might be an optimal set of parameters that maximizes the network robustness. In particular, Fig. 2 (bottom) shows that for syn = 5–6 ms (with Iµ = 3 µA/cm2,gsyn = 0.25 mS/cm2), the system still expresses near-synchronous oscillations in the face of >11%Het—for larger and smaller values of syn, the robustness is weaker.
In Table 1, we show parameter sets (Iµ, syn, gsyn) for which robust oscillations were obtained in our simulations. We only list cases for which the robustness or maximal %Het exceeded 7%, and this maximal %Het for each parameter set is shown. In addition, network frequency and phase information is indicated. For all the cases shown in Table 1, the network pattern at the maximal %Het is the near-synchronous described pattern 1 (see Fig. 1A for illustration). With heterogeneity beyond this maximal %Het, harmonic locking or suppressed states occurred, and asynchronous states (as defined in the preceding text) sometimes preceded the suppressed states as the heterogeneity was increased. For most of the parameter sets, %Het values less than the maximal %Het shown in Table 1 also resulted in near-synchronous oscillations. However, near-antiphase, varied phase-locking and asynchronous patterns sometimes occurred at lower %Het values (see Fig. 3), as well as bistable patterns (see Fig. 4). There would be an additional three cases to include in Table 1 from our simulations if we also encompassed maximally robust oscillations with varied phase-locking patterns (3; see Fig. 1B for illustration). In those cases, a further increase in heterogeneity gave way to suppressed behaviors at larger %Het values. These observations are naturally limited by the resolution of the simulations performed (see METHODS).
|
|
|
; White et al. 1998), the network frequency increased with decreasing syn. Any two rows (with different syn) in Table 1 where the %Het is the same can be compared to see that this is the case.
Because it was found that syn/T values separated tonic and phasic regimes [i.e., in reference to different ways in which coherence is lost (White et al. 1998)], we also calculated syn/T for our network oscillations. This is shown in Table 1 for the robust oscillations with >7%Het. We found syn/T to have "small" values in accordance with (White et al. 1998), who found that syn/T was smaller (<1, for their cellular model) for the phasic (as opposed to the tonic, >2) regime where coherence occurs. Together with the observation from our simulations that there seems to be an optimal set of parameters (Iµ, gsyn, syn) that maximizes the network robustness, one can define a window of syn/T values for which maximal robustness occurs.
As derived from several simulations, a nonmonotonic dependence of robustness on each of the parameters gsyn, Iµ, and syn occurs. This is shown in Fig. 5. The exact values for maximal %Het and syn/T of Fig. 5 are given by the boldfaced values in Table 1 for maximal %Het values or robustness measures exceeding 7%. These data indicate that it might be possible to determine a set of model parameters that give rise to maximally robust network oscillations. However, our simulations also show that the dynamics in these two-cell networks express much complexity. For example, Fig. 4 shows an example of bistability, and Fig. 3 shows how several dynamic patterns can be expressed for small changes in heterogeneity. The patterns shown in Fig. 3 were obtained for a chosen set of initial conditions (see METHODS). However, given the nonlinearity of the system, other stable patterns likely exist for these same parameter values. For the set of initial conditions used, Fig. 3 shows that increasing the heterogeneity can be helpful in bringing about "coherence" or near-synchronous solutions.
|
To establish the validity of our simulation observations regarding the existence of maximal robustness, we embarked on a bifurcation analysis (see description in METHODS). Our observations might depend on the nonlinear relationship between Iapp and the intrinsic frequency, but it is not obvious that a nonmonotonic dependence of robustness on parameters as shown in Fig. 5 should exist. We performed a series of numerical continuations of the two-cell heterogeneous network system. Our simulations led us to focus mainly on gsyn = 0.15, 0.25, 0.5 mS/cm2 values as where the most robust network oscillations were present (see Table 1). Continuations using either Iµ or as the bifurcation parameter (see Tables 2 and 3) for a range of syn values were done. We were able to determine the existence and stability of both the equilibria and periodic solutions of the system as the Iµ or parameter was varied while the other parameters remained fixed. This allowed us to obtain the precise dependence of coherence in the model network system on variations in specific parameter values.
|
syn (as in Fig. 6, A and B), the system went to a suppression oscillation when the near-synchronous oscillations lost stability. For smaller values of syn, the region of stability of the suppression oscillations decreased, and when the near-synchronous oscillations lost stability, the system exhibited either asynchronous or harmonically locked oscillations. We summarize the results from several numerical continuations in Table 2. All these particular continuations were performed using = 0.05 µA/cm2. Note that these continuations were performed with Iµ as the bifurcation parameter, and so the other parameters (including ) were fixed at particular values as indicated in Table 2. For each case, the minimum Iµ for which near-synchronous behavior is present is indicated (e.g., see Fig. 6B example). This Iµ together with determines the %Het for the system and this is also shown in Table 2. Note that decreasing Iµ corresponds to increasing the %Het in the system because of the nonlinear relationship between Iapp and the intrinsic frequency of the isolated WB model cell (Wang and Buzsáki 1996). For this choice of and with gsyn = 0.25 mS/cm2, we can see that as syn decreases, one obtains a %Het that increases until it reaches a maximum and then decreases. In other words, for a fixed difference in external input (as given by ), the minimal required amount of external drive to give rise to near-synchronous oscillations can both increase and decrease with increasing syn.
Near-synchronous and suppression solutions were found and followed in our numerical continuation studies. Other solutions such as varied phase-locking and antiphase oscillations are present in various parameter regimes (as found in the 2-cell simulations) but were not systematically followed in the studies performed here. A typical bifurcation diagram with as the bifurcation parameter is shown in Fig. 6C. As (and hence %Het) was increased, the stable near-synchronous oscillations were lost in one of two ways. They either ceased to exist via a saddle node of limit cycles (as shown in Fig. 6C) or they lost stability in a period doubling bifurcation before the saddle node of limit cycles. For large syn, when the stable near-synchronous oscillations were lost, the system exhibited suppression oscillations. As syn was decreased, the region of existence of the suppression oscillations was diminished. For syn small enough (e.g., see Fig. 6C), this region no longer overlapped with the region of stability of the near-synchronous oscillations. The system exhibited either harmonically-locked or asynchronous oscillations when the stable near-synchronous oscillations were lost, and suppression oscillations at larger values of or %Het as shown in Fig. 6C. In Table 3, we indicate the largest value of for which stable near-synchronous oscillations exist for each given set of parameter values as determined from several numerical continuations. As syn decreases, it is clear that a maximal %Het exists. Specifically, for gsyn = 0.25 mS/cm2, Iµ = 3 mA/cm2, and syn = 5.7 ms, there is a maximal %Het of 11.75%. We conclude that the nonmonotonic dependence of robustness on different parameters as observed in the simulations (see Fig. 5) is a well-defined phenomenon.
From 2- to 10-cell heterogeneous networks
Now that we have established that there exists a set of parameters that give rise to maximally robust oscillations, we would like to determine whether coherence properties observed in the two-cell networks are borne out in larger networks. That is, are the parameter regimes for robust oscillations in two-cell networks similar to parameter regimes in larger networks? Consider the simplified view that if one randomly distributes the amount of heterogeneity up to a maximal amount between a set of cells, then the possible patterns that the network expresses would include those observed in the two-cell network for the given heterogeneity. Therefore parameter values that give rise to coherent oscillations in larger (>2-cell) networks could be predicted from what is known about coherent oscillations in the two-cell networks.
In Table 4 we show coherence measure values obtained from several 10-cell network simulations. For the top part of Table 4, the simulations used parameter values of gsyn = 0.25 mS/cm2, Iµ = 3 µA/cm2, and syn = 1 or 5 ms. Using these same parameters, two-cell network simulations produced coherent oscillations 
5.8%Het with syn =1 ms, and 11.4%Het with syn = 5 ms. Furthermore, near-synchronous solutions were obtained when smaller %Het values were used in the two-cell simulations (using initial conditions given in METHODS). With syn = 1 ms, the 10-cell network simulations produced coherent oscillations with 3%Het, but not with 8%Het as given by their coherence values shown in Table 4. When syn was increased to 5 ms, the 10-cell networks with 3%Het were still very coherent, but now the ten-cell networks with 8%Het were also coherent. That is, using 8%Het, the coherence measure value obtained with syn = 5 ms is larger than with syn = 1 ms. However, as with syn = 1 ms, the coherence value obtained with 8%Het is smaller than with 3%Het but to a lesser extent. Given that the two-cell networks show "coherence" for 5.8%Het with syn = 1, but for 11.4%Het with syn = 5, it would appear that the coherence properties of two-cell heterogeneous networks are "preserved" in the larger networks. The 10-cell simulations were re-done with the addition of noise (see METHODS). We used Gaussian noise with SDs that were half of, equal to, and twice that of the level of heterogeneity in the system. In all cases, the difference in coherence, as "predicted" from the two-cell simulations was maintained (see Table 4).
As shown in Figs. 4 and 6B, bistable patterns can occur in the two-cell heterogeneous network system. The bistable pattern illustrated in Fig. 4 was obtained using parameters ofgsyn = 0.25 mS/cm2, Iµ = 1 µA/cm2, and syn = 1 ms. The initial conditions were such that the initial voltages were either different by 4 mV, as given in METHODS, or the same for the two cells, V1 = V2 = –59.5567. The initial values for the other variables (h, n, s) of each model cell were the same for the two cases and are given in METHODS (see Fig. 4 legend). With increasing heterogeneity, for the first case of initial conditions, we obtained near-antiphase patterns (4%Het), then varied phase-locking, and then near-synchronous patterns up to an 8%Het. Harmonic locking patterns arose with further increases in heterogeneity. For the second case of initial conditions, we obtained varied phase-locking and then near-synchronous patterns with increasing heterogeneity 8%Het, and then harmonic locking as in the first case of initial conditions. (The particular patterns for the case of 4%Het are shown in Fig. 4.) Using these parameter values in a 10-cell network configuration, we find again that the two-cell dynamics "predict" the coherence properties of the 10-cell networks. Consider two different initial conditions for the 10-cell networks that are derived from the two cases of initial conditions for the 2-cell networks. They are either the initial voltages for the 10 cells were different (but within 5 mV of each other) or the initial voltages were all the same value. As in the two-cell network simulations, the initial values for the other variables (h, n, s) of each model cell were the same for the two cases. We will refer to the first and second case of initial conditions (ICs) for the 10-cell networks as ICsA and ICsB, respectively, and they are given with Table 4. The bottom part of Table 4 provides coherence values obtained from 10-cell network simulations with ICsA and ICsB. Because the heterogeneity is spread randomly among the 10 cells, and if the dynamic behaviors seen in the 2-cell networks occur in the larger networks, then one might expect to see larger coherence values with ICsB than with ICsA. We found this to be the case with 5%Het. At and below this level of heterogeneity, the 2-cell network dynamics mostly included near-antiphase (i.e., noncoherent) patterns with initial conditions given by the first case, and so the vast difference in coherence values observed in the 10-cell simulations is supportive of the predictive nature of the 2-cell network dynamics for larger networks. However, when we increased the heterogeneity further to 8%Het in the 10-cell simulations, this large difference in coherence values using the two different ICs was no longer apparent. Because at and below this amount of heterogeneity, both sets of ICs now mostly involve near-synchronous patterns in two-cell networks, this does not negate the predictive aspect of the two-cell dynamics for larger (10-cell) network coherence. It suggests that noncoherent patterns (such as antiphase) need to dominate the two-cell network dynamics (at the particular level of heterogeneity) to manifest their effect in larger networks. Interestingly, when noise is added to the simulations (e.g., see the Table 4, last column, which is boldfaced), the dependence on the exact level of heterogeneity is removed. For both 5%Het and 8%Het, the 10-cell networks have a much lower coherence value when ICsA are used as near-antiphase patterns occur in the two-cell networks using the first case of initial conditions described in the preceding text. The preceding results support the possibility that coherence properties in larger networks can be predicted from dynamic patterns obtained in two-cell networks. However, one should caution that such observations would depend in some fashion on the particular choices of initial conditions, noise levels, and the dynamic constraints on different stable patterns expressed by nonlinear systems, as has been examined for homogeneous inhibitory networks (Golomb and Rinzel 1994). We examine this further in Skinner et al. 2005.
|
DISCUSSION |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
syn = 5.7 ms, where robust network oscillations with 12%Het occur (see Fig. 2 and Table 3). With these values, the two-cell heterogeneous network exhibits oscillations of 90 Hz frequency and with about a 10% phase lag between the two cells. Although a direct comparison with the 100-cell network simulation results of (Bartos et al. 2001, 2002; Wang and Buzsáki 1996) cannot be made, our results go some way in explaining the more robust and higher frequency oscillations observed by Bartos et al. (2001, 2002) in their simulations as compared with simulations performed by Wang and Buzsáki (1996). In addition, our simulations show that 2-cell coherence properties are maintained in larger 10-cell networks. Together with other studies by White et al. (1998) showing similar qualitative behaviors in 2, 10, and 100-cell networks, our work suggests that using the amount of heterogeneity (i.e., ) as a parameter in (2-cell) network studies allows the determination of optimal parameter values, which in turn are predictive of coherent regimes in larger networks. Biophysical detail and external input
In developing mathematical models of neuronal networks that include some sort of biophysical characteristics, it is useful to have constraints on a chosen model's parameters. Questions regarding model detail naturally arise: how much and what? For inhibitory network models targeted toward hippocampal cortex, there are network model explorations using cellular inhibitory models that are: detailed multi-compartment representations (e.g., Maex and De Schutter 2003; Traub and Bibbig 2000), single-compartment representations (e.g., Baker et al. 2002; Jalil et al. 2004; Tiesinga and José 2000; Wang 2002; Wang and Buzsáki 1996; White et al. 1998), or simpler integrate-and-fire units (e.g., Brunel and Wang 2003; Neltner et al. 2000). To a large extent, the amount of biological detail that makes sense to include in a neuronal model depends on the question being asked. While analyzing and determining underlying cellular mechanisms in network models give rise to constraints on some parameters for given dynamic outputs, one is still often faced with few or no constraints on other model parameters. A problem faced by all modelers is how best to both incorporate biophysical detail and do mathematical analyses to understand and predict the model output. The intrinsic cellular properties of a model always make a difference to the network output, but it may not be a critical difference. While an examination of time constant values and active voltage ranges for particular ion currents is suggestive of the importance of particular intrinsic properties, this is a qualitative estimation that is subject to nonintuitive, nonlinear outcomes. We suggest that our criterion of maximal heterogeneity could be used as a basis for evaluating how much of a "difference" various intrinsic details make to network coherence. That is, would the optimal Iµ, gsyn, syn values change with, say, using three (rather than 2) types of potassium currents?
A critical parameter in network models is the external or applied current, Iapp. However, there is no well-defined way to measure this parameter value as it represents the summation of synaptic and background current (noise and heterogeneities) in the dendritic tree, which in turn depends on the integration properties of the particular cell and the network in which it resides. Depending on whether the individual cell model includes a dendritic component and on what network size is being considered, the input could be described by a tonic current or Poisson statistics with or without noise, heterogeneities, and correlations (Brunel and Wang 2003; Destexhe et al. 2003; Maex and De Schutter 2003; McMillen and Kopell 2003). Interestingly, Maex and DeSchutter (2003) show that the difference between using tonic or afferent fiber input (Poisson statistics) is not important in extracting a resonant synchronization in their inhibitory networks. A modification of our suggestion above would be to determine the Iµ (i.e., Iapp) value that gives the most robust oscillations for the given cellular model. In this way, the choice of Iapp would be less arbitrary than simply using values that give rise to certain network frequencies.
Network modeling strategies and related works
The best strategy to use if one wants to understand network dynamics in which the biophysical details are not ignored is unclear. However, strategies for particular problems have been developed (Ermentrout and Chow 2002). For example, for problems involving long-distance synchronization, a mapping strategy was introduced by (Ermentrout and Kopell 1998) in which the biophysical equations can be reduced to a one-dimensional map in terms of spike timings of different neurons. The subsequent map is then easily analyzed and its predictions verified with simulations (e.g., Acker et al. 2003). A geometric approach, in which several trajectories (each corresponding to 1 cell) move around in a lower dimension phase space has been used to understand thalamic oscillations (Rubin and Terman 2000). Otherwise, one can test predictions from simpler (integrate-and-fire) models using more detailed biophysical models (e.g., see Wilson et al. 2004). A more general strategy using the theory of weakly coupled oscillators (Kuramoto 1995) allows one to mathematically reduce the network system to a set of equations only involving phases. These equations are then straightforward to analyze. While this approach is powerful and has been used to analyze networks of coupled cells (e.g., Chow 1998; Lewis and Rinzel 2003; Neltner et al. 2000), it is limited to neuronal units that exhibit stable limit cycles and that have sufficiently weak coupling.
The ability to continue solutions when doing bifurcation analyses is constrained by the complexity of the set of nonlinear equations. As shown in our work here, this is feasible using two-cell networks but is clearly impractical for large networks. However, because it seems that the coherent properties in two-cell networks can be preserved in larger networks, it might be a useful small-to-large network modeling strategy to determine optimal external input and/or intrinsic biophysical parameters from two-cell network analyses via the maximal heterogeneity criteria suggested above. Depending on what experimental data are available, it might also be reasonable to determine some optimal synaptic parameters in this way. Then, in performing larger network simulations, one would already have intuition from the two-cell network dynamics so that it may be possible to fundamentally understand any observed changes in coherence patterns of the large networks. In addition, one could attempt to mathematically reduce the intrinsic cell model with its optimal parameters to subsequently use in larger network simulations.
Our results are in accordance with other modeling studies of heterogeneous, inhibitory networks. In the 100-cell simulation studies by Bartos et al. (2001, 2002) and Wang and Buzsáki (1996), a Gaussian distribution of heterogeneities was used and different amounts were explored for specific parameter values (e.g., see Fig. 5 in Wang and Buzsáki 1996). Neltner et al. (2000) defined robustness of synchrony as the critical disorder at which the asynchronous state becomes linearly stable. With this definition, they obtained a maximal robustness using WB neurons when synaptic conductance was varied, for a given time constant and mean frequency. In another modeling study, Tiesinga and José (2000) suggest that stochastic weak synchronization (as opposed to strong synchronization), in which the particular cycle in which each cell fires is random, might underlie robust oscillations in inhibitory networks.
Basket cell networks in hippocampus
This study was partly motivated out of concern for how to constrain cellular models used in inhibitory networks. Developing models of individual interneurons in hippocampus that encompass the richness of their intrinsic properties is a major undertaking (e.g., Saraga et al. 2003). Furthermore, there is usually not enough experimental data for particular interneurons to warrant developing a very detailed model. However, if one then uses a simplified model, it is unclear how best to constrain the chosen parameters. Our study here suggests how optimal parameters might be determined. It is interesting to note that Bartos et al. (2002) found differences in synaptic time constants and peak synaptic conductances in DG, CA3, and CA1 regions of hippocampus in inhibitory basket cells. One could envisage using the particular syn and gsyn values to select out an Iµ value via our maximal heterogeneity criterion. This Iµ value in turn would give rise to a certain network frequency. One might speculate that particular hippocampal regions are "tuned" to particular frequencies. Basket cell networks are also connected by gap junctions, which adds another layer of complexity (Freund 2003; Fukuda and Kosaka 2000). This further emphasizes the importance of trying to understand local inhibitory circuits in the brain from the perspective of biophysical parameter values. Our work shows how this could be approached.
|
GRANTS |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
|
FOOTNOTES |
|---|
Address for reprint requests and other correspondence: F. K. Skinner, Toronto Western Research Institute, University Health Network, 399 Bathurst St., MP13-317, Toronto, Ontario M5T 2S8, Canada (E-mail:fskinner{at}uhnres.utoronto.ca)
|
REFERENCES |
|---|
|
TOPABSTRACTINTRODUCTIONMETHODSRESULTSDISCUSSIONGRANTSREFERENCES |
|---|
Baker PM, Pennefather PS, Orser BA, and Skinner FK. Disruption of coherent oscillations in inhibitory networks with anesthetics: the role of GABAA receptor desensitization. J Neurophysiol 88: 2821–2833, 2002.
Bartos M, Vida I, Frotscher M, Geiger JRP, and Jonas P. Rapid signaling at inhibitory synapses in a dentate gyrus interneuron network. J Neurosci 21: 2687–2698, 2001.
Bartos M, Vida I, Frotscher M, Meyer A, Monyer H, Geiger JRP, and Jonas P. Fast synaptic inhibition promotes synchronized gamma oscillations in hippocampal interneuron networks. Proc Natl Acad Sci USA 99: 13222–13227, 2002.
Bragin A, Jandó, G, Nádasdy Z, Hetke J, Wise K, and Buzsáki G. Gamma 40–100 Hz. oscillation in the hippocampus of the behaving rat. J Neurosci 15: 47–60, 1995.[Abstract]
Brunel N and Wang X-J. What determines the frequency of fast network oscillations with irregular neural discharges? I. Synaptic dynamics and excitation-inhibition balance. J Neurophysiol 90: 415–430, 2003.
Buzsáki G. Theta oscillations in the hippocampus. Neuron 33: 325–340, 2002.[CrossRef][ISI][Medline]
Buzsáki G and Chrobak J. Temporal structure in spatially organized neuronal ensembles: a role for interneuronal networks. Curr Opin Neurobiol 5: 504–510, 1995.[CrossRef][ISI][Medline]
Chow CC. Phase-locking in weakly heterogeneous neuronal networks. Physica D 118: 343–370, 1998.[CrossRef]
Cohen SD, and Hindmarsh AC. CVODE, a stiff/nonstiff ODE solver in C. Comput Phys 10: 138–143, 1996.
Destexhe A, Rudolph M, and Paré, D. The high-conductance state of neocortical neurons in vivo. Nat Rev Neurosci 4: 739–751, 2003.[CrossRef][ISI][Medline]
Doedel EJ. AUTO: A program for the automatic bifurcation analysis of autonomous systems. Congr Numer 30: 265–284, 1981.
Ermentrout GB. Simulating, Analyzing, and Animating Dynamical Systems: A Guide to XPPAUT for Researchers and Students SIAM, Philadelphia, http://www.math.pitt.edu/bard/xpp/xpp.html.2002.
Ermentrout GB, and Chow CC. Modeling neural oscillations. Physiol Behav 77: 629–633, 2002.[CrossRef][Medline]
Ermentrout GB, and Kopell N. Fine structure of neural spiking and synchronization in the presence of conduction delays. Proc Natl Acad Sci USA 95: 1259–1264, 1998.
Freund TF. Interneuron Diversity series: Rhythm and mood in perisomatic inhibition. Trends Neurosci 26: 489–495, 2003.[CrossRef][ISI][Medline]
Freund TF, and Buzsáki G. Interneurons of the hippocampus. Hippocampus 6: 347–470, 1996.[CrossRef][ISI][Medline]
Fukuda T and Kosaka T. Gap junctions linking the dendritic network of GABAergic interneurons in the hippocampus. J Neurosci 20: 1519–1528, 2000.
Golomb D and Rinzel J. Clustering in globally coupled inhibitory neurons. Physica D 71: 259–282, 1994.
Jalil S, Grigull J, and Skinner FK. Novel bursting patterns emerging from model inhibitory networks with synaptic depression. J Comput Neurosci 17: 23–37, 2004.
Kuramoto Y. Collective behavior of coupled oscillators. In: The Handbook of Brain Theory and Neural Networks, edited by Arbib M. Cambridge MA: MIT Press, 1995, p. 203–206.
Lewis TJ and Rinzel J. Dynamics of spiking neurons connected by both inhibitory and electrical coupling. J Comput Neurosci 14: 283–309, 2003.[CrossRef][ISI][Medline]
Maex R and De Schutter E. Resonant synchronization in heterogeneous networks of inhibitory neurons. J Neurosci 23: 10503–10514, 2003.
McBain C and Fisahn A. Interneurons unbound. Nat Rev Neurosci 2: 11–23, 2001.[ISI][Medline]
McMillen D and Kopell N. Noise-stabililized long-distance synchronization in populations of model neurons. J Comput Neurosci 15: 143–157, 2003.[CrossRef][ISI][Medline]
3, i.e., the open circle. C: 
bifurcation parameter results
, the maximal %Het (for which coherence occurs) as given by the left-hand ordinate axis; 
,