Functional Phase Response Curves: A Method for Understanding Synchronization of Adapting Neurons

Experimental

MODEL REFERENCE CURRENT INJECTION (MRCI) DYNAMIC CLAMP.

We recorded the membrane potential (V) of the spiking neurons and used dynamic clamp (Sharp et al. 1993a,b) to generate the synaptic current input to the spiking neurons on abdominal ganglion of Aplysia californica. The excitatory or inhibitory stimulus current was generated from an alpha-shaped synaptic conductance waveform by the MRCI dynamic-clamp system (Butera et al. 2001; Raikov et al. 2004). The equations to generate the alpha shaped conductance and the injected current are given by

(1)

Here, I_syn represents the synaptic current injected into the neuron; α(t) represents an alpha-shaped waveform with a time constant τ, normalized by the term e/τ to a height of one, which provides a conductance of unit amplitude regardless of the value of τ; g_syn represents the maximum conductance; i_trig is set as a square pulse with amplitude 1 and duration of 1 ms; and E_syn is set at 0 mV for excitatory input and −70 mV for inhibitory input.

Use of this differential form to generate an alpha-shaped conductance can accommodate any number of inputs arriving at any time without the need for any additional bookkeeping of input history.

Theoretical

MODEL NEURON NETWORKS.

A network consisting of two single-compartment WB model neurons was used. We chose this model for two reasons: 1) it is a relatively simple but sufficient model to produce quite similarly shaped PRCs to the Aplysia neurons that we record from and 2) it is easily modified to include adaptation. It has been shown that these Aplysia spiking neurons exhibit type I PRCs (Preyer and Butera 2005) for a wide range of excitatory stimulus magnitudes. The WB model produced type I PRCs not just of similar shape but of similar changes in shape as the conductance stimulus amplitude is increased, as revealed previously (Preyer 2007). The methods presented herein are general and not dependent on any particular model or on any particular PRC shape.

The model neurons were reciprocally coupled by identical inhibitory or excitatory synapses. The equations for each neuron are given by

where the capacitance C is 1 μF/cm²; V is the cell membrane voltage in mV; t is time in unit of ms; h is the inactivation of the sodium current; n is the activation of the potassium current; I_z is the adaptation current; and m_∞ is the steady-state activation of sodium current. Reversal potentials, E_Na, E_K, and E_L, are 55, −90, and −60 mV, respectively. The maximum conductances for sodium g_Na, potassium g_K, and leak g_L are 35, 9, 0.1 mS/cm², respectively. The scale factor γ is set to 5.

The adaptation current I_z, based on the formalism described by Ermentrout (1998), is given by

where a = 0.005 ms⁻¹, b = −5 mV, k = 2.0 mV, and g_z = 1 mS/cm². Because of the steep threshold for activation at −5 mV, this type of adaptation is induced only during spikes, similar to a calcium-dependent afterhyperpolarizing current (I_AHP) but not to a muscarinic current (I_M), which allows adaptation to be produced at subthreshold voltages (Ermentrout 1998; Prescott et al. 2006).

The synaptic current I_syn is given by

where g_syn is the maximum synaptic conductance; E_syn is the reversal potential, which is set to −75 mV for inhibitory perturbation and 0 mV for excitatory perturbations; s is the gating variable; V_pre is the voltage of the presynaptic cell; k_s = 6.25 ms⁻¹ is the rate constant of the synaptic activation; τ_syn is the synaptic decay time constant; and I_app represents the applied current to the model cell to control the spiking frequency. In our model networks, the two cells have identical values of I_app. The differential equations were integrated using a driver program written in C that called the Fortran subroutine radau5 (Hairer and Wanner 1991), using a variable step implicit fifth-order Runge–Kutte method on a Beowulf cluster running CentOS.

Experimental

When the fPRC protocol was applied as in Fig. 1B, the most dramatic changes in ISI length occurred when the train of fixed-delay stimuli was turned on or off. The ISI between the first and second red dots (after the pulse train was turned on) is much smaller than the intrinsic free-running period, as measured between the last blue dot and the first red dot; this change is primarily due to the effects of the phase resetting that produces an advance and therefore a faster frequency. On the other hand, the ISI between the last red dot and the first green dot (when the pulse train was turned off) is longer than the intrinsic period, presumably due to the increase in activation of an adaptive current caused by the increase in frequency. The corresponding time sequence of ISIs in the process of measuring this particular fixed delay (a single point in the fPRC) is illustrated in Fig. 1C. The ISI before the repetitive stimuli were turned on (blue dots) was almost steady at 420 ms; the ISI dropped to 220 ms (first red dot) once the pulse train was switched on at a fixed-delay time of 160 ms (indicated by the cyan star) and gradually reached a steady value around 280 ms (red dots); the ISI increased to 500 ms when the repetitive perturbations were turned off (green dots) and gradually reached the intrinsic period. The fixed-delay time was swept through values from zero to that of a single intrinsic period in a random order to minimize any systematic error in the experiment.

The time sequences of the ISIs reordered by the fixed-delay time for inhibitory and excitatory perturbations are illustrated in Fig. 2, A and B, respectively. The blue, green, and red dots are color-coded the same as in Fig. 1C. For excitation, the dramatic changes in the ISI, when the train of fixed-delay perturbations was turned on and off, are more obvious at smaller phases than at larger phases, whereas the converse is true for inhibitory perturbations. The duration of the transient period prior to converging to a stable ISI was longer at earlier phases for excitation and later phases for inhibition. This was true both for the transient at the onset of pulse stimulation and for the transient where the stimulation was turned off.

Before measuring the fPRC with a repetitive perturbation, we also measured a single-pulse PRC with the same experimental conditions. For a larger magnitude of perturbation [Fig. 3, A1 (inhibition) and B1 (excitation)] fPRCs (stars) deviated dramatically from the first-order PRCs (dots). For a smaller magnitude of perturbation [Fig. 3, A2 (inhibition) and B2 (excitation)] fPRCs are close to the first-order PRCs because both the second-order resetting and any adaptation are negligible. For excitatory cases, we found the results of Fig. 3B to be representative of all neurons measured. Specifically, in over 20 neurons measured, the first-order single-pulse PRC (f₁) exhibited advances at all phases; the second-order single-pulse PRC (f₂) exhibited delays at all phases; the fPRC was similar to the summation of f₁ and f₂, (f₁ + f₂), for sufficiently weak inputs; and the fPRC was distinct from f₁ + f₂ for strong inputs.

Computational and theoretical

Motivated by the experimental demonstration of adaptation and its effect on the functional PRC, we then attempted to reproduce similar results in a model, to analyze the stability of the protocol, to use the data to predict the synchronization behavior of a circuit comprised of two adapting neurons, and to test these predictions.

The simulated measurements of fPRCs for inhibitory and excitatory perturbations are shown in Fig. 4, A and B, respectively, which capture a trend similar to that of the experimental observations as shown in Fig. 2, A and B. The simulations also capture the longer convergence times for excitation for shorter delays and for inhibition for longer delays. If the currents are adapting, the protocol of measuring fPRCs allows the currents to adapt completely. The net result is a periodic steady state in which the neuron receives periodic inputs and is engaged in a 1:1 phase-locked periodic rhythm. Thus the fPRC accounts for the degree of adaptation that normally occurs during a periodic rhythm in a network context.

The fPRC obtained in Fig. 4, A and B is plotted as the dot-dashed curve in Fig. 4, C and D. Note that the fPRC does not exist at late phases due to stability problems described in the following text. Figure 4C shows that the first-order resetting (f₁, solid line) for inhibition consists of delays, whereas the second-order (f₂, dotted curve) and the third-order (f₃, long dashes) phase responses consist of advances. The sum of the three resetting components (f₁ + f₂ + f₃, short dashes) approaches the functional PRC and might approach it more closely if higher-order terms were included, but would require a fine resolution and the analysis of many cycles. On the other hand, Fig. 4D shows that the first-order resetting for excitation consists of advances, whereas the second- and third-order resettings consist of delays. Both the first-order resetting and the fPRC are limited by causality. In other words, the amount of advance produced is limited by the time remaining until the next spike would have occurred anyway. The causality limit falls along the line P_i(1 − f₁). For both the simulated and experimental cases of excitatory input (Figs. 4D and 3B1), the fPRC differs from the first-order PRC at early phase, but converges on it at later phases when the causality limit is reached.

For weak perturbations, the assumption that the resetting in subsequent cycles adds linearly is valid, so the sum of the first-, second-, and higher-order PRCs converges to fPRCs. The fPRCs are very close to the first-order PRCs because second- (and higher-) order resetting is small (Fig. 3, A2 and B2). The fPRCs are different from the first-order PRCs when the perturbation is strong in both the excitation and inhibition cases, as shown in Figs. 3, A1 and B1 and 4, C and D. For strong perturbations, first-order PRCs advance (delay) phase, whereas second- and third-order PRCs delay (advance) phase for excitation (inhibition). However, the second- and third-order resettings are not necessarily manifested under repetitive stimulation for strong perturbations because the neuron does not return to the original limit cycle between perturbations. In other words, the assumption of pulsatile coupling is violated. Therefore including the second- and third-order resettings can produce a worse match to the functional PRC at most phases than ignoring it (Fig. 4D). For less severe perturbations (experimental Fig. 3), adding second-order resetting can slightly improve the correspondence with the fPRCs.

Adaptation can be interpreted as changing the intrinsic period of the oscillator due to that activation of the adaptation current. The adaptation current used in the model summates temporally only during a spike and decays in between spikes. Therefore the greater the change in number of spikes per unit time, the greater the change in the level of activation of the adaptation current. An estimate of the amount of adaptation that has occurred can be obtained by observing the length of the ISI immediately after the fixed-delay protocol is terminated. We refer to this estimate as the modified intrinsic period (P_m in Fig. 4, A and B) and to the steady period achieved during the fixed-delay protocol as the network period (P_n in Fig. 4, A and B). The actual phase at which an input is received cannot be determined by using the delay normalized by P_i if P_m differs significantly from P_i; instead one should normalize by P_m. Furthermore, the single-pulse PRC measured in the absence of adaptation may be quite different from that measured in its presence, if that were possible. Therefore single-pulse PRCs are not strictly applicable in the presence of significant adaption, as we will show.

Note that the inhibitory PRC (Fig. 4C, solid line) is larger at late phases, whereas the excitatory PRC (Fig. 4D, solid line) is larger at early phases. Slow convergence to a stable phase relationship results because phase resetting changes the period, which causes adaptation that changes the effective intrinsic period (P_m) in the opposite direction. For inhibition, resetting increases the period, which decreases the activation of the adaptation current, which decreases the effective intrinsic period. Decreasing the effective intrinsic period increases the effective phase, which generally further increases the observed network period, slowing and opposing the effect of adaptation until the modified intrinsic period stops changing (adaptation is complete). A complementary process occurs for excitatory inputs. The interplay between phase resetting and adaptation is complex. Each pulse has an effect on the timing of the next spike due to resetting the phase of the ongoing oscillation, but it also affects the timing of the next spike by altering the amount of time available for the adaptation current to decay—and these two effects usually oppose each other. It is not easy to predict whether the opposing effects will cancel each other out or whether one will dominate.

Existence and stability criteria for two-neuron network

The theoretical method for predicting the existence and stability of phase-locked modes using single-pulse PRCs for networks of two bursting neurons has been previously developed (Oprisan et al. 2004; Sieling et al. 2009). We describe it here to show how the new methods based on the fPRC are an improvement. The assumed firing pattern in a two-cell network with a 1:1 locking is given in Fig. 5 A1. Neuron 1 (2) sends input to neuron 2 (1) when it fires. The recovery interval (tr) is defined as the time elapsed between when a neuron receives an input and when it next fires. The stimulus interval (ts) in a network is the time between when a neuron fires and when it receives an input. By definition of 1:1 phase locking, tr₁[n − 1] = ts₂[n] and tr₂[n] = ts₁[n].

The assumptions of the single-pulse PRC method can be summarized as follows: 1) each neuron in the network can be represented as a limit cycle oscillator, 2) the trajectory of each neuron returns to its limit cycle between inputs (this is the assumption of pulsatile coupling), and 3) the input received by each neuron has the same effect in the closed-loop circuit as in the open-loop circuit used to generate the phase response curves. The key idea is that under the assumption of pulsatile coupling, the first-order PRC (f₁) and the second-order PRC (f₂) can be used to calculate the values of the intervals as follows: ts = P_i[φ + f₂(φ)], tr = P_i[1 − φ + f₁(φ)], such that tr is a function of ts, tr = g(ts). When phase-locking occurs, tr₁ = ts₂ and tr₂ = ts₁, so the intersection of the tr₁ = g(ts₁) and ts₂ = g⁻¹(tr₂) curves yields the fixed points. The stability is calculated using the slopes of the first- and second-order phase resetting curves as in Oprisan et al. (2004).

In contrast, the fPRC requires only the third assumption just cited, that the fPRC is generated with an input similar to the one it will receive in the network. The key idea of the functional PRC is to characterize the response of each oscillator to inputs arriving periodically with a fixed delay to characterize the response of the steady phase-locked mode when the delay times within the network are perturbed away from the steady value. The functional PRC protocol (Fig. 5A2) allows the tabulation of the steady-state recovery interval (tr) in a given neuron as a function of the delay (δ) between when the neuron fires and when it next receives an input, such that tr = g(δ). Figure 5A3 illustrates the assumption that in a closed-loop circuit the recovery intervals for each neuron can be calculated as g(δ) by assuming that at steady state the stimulus interval is the fixed delay δ. Thus at steady state in a periodic locking δ₂ = g₁(δ₁) and δ₁ = g₂[g₁(δ₁)]. The modes that satisfy the existence criteria can be located graphically (Fig. 5B) at the intersection of tr₂[∞] = g₁(ts₁[∞]) and ts₂[∞] = g₂⁻¹(tr₂[∞]), where the minus 1 indicates the inverse function and tr₁[∞] = ts₂[∞] = δ₂ and tr₂[∞] = ts₁[∞] = δ₁. Although the delays are not fixed in the closed-loop circuit, the stability of the locking can be approximated by assuming a perturbation of two neurons that are actually coupled by fixed delays. Linearizing the map δ₁[n] = g₂{g₁(δ₁[n − 1])} provides a stability criterion that −1 < g*₂(ts*₂)g*₁(ts*₁) < 1, where the asterisk indicates the value at the fixed point and the prime indicates the slope at the fixed point. Note that all intervals are defined in terms of time rather than phase. The modes that satisfy the existence criteria are located graphically by finding the intersection of tr₁ = g₁(ts₁) and ts₂ = g₂⁻¹(tr₂), where tr₁ = ts₂ and tr₂ = ts_1. Then the stability criterion is applied to determine whether the mode is stable, that is, whether the inevitable perturbations away from the presumed fixed-delay times in the network will grow or shrink. Figure 5B shows the intersection of the tr₁ = g₁(ts₁) and ts₂ = g₂⁻¹(tr₂) curves for the WB model with adaptation and mutual excitation. The reciprocal fixed-delay paradigm described here does not completely characterize the circuit because it neglects the feedback interaction between oscillators, although the results in Fig. 5, C and D described in the following text confirm that this is a reasonable approximation. The idea is that if each neuron can lock to an input given at a fixed delay, when they are reciprocally coupled with the same effective (but not fixed) delays, the resultant locking should be stable.

The intersection method given earlier does not find exact synchrony because it presumes an alternating firing order. A separate check is run to determine whether the recovery intervals at a phase of zero are measurable, indicating that the fixed-delay protocol produces a stable locking and that they are equal. Period two modes that satisfy the mapping, tr₂[n] = g₂〈g₁{g₂[g₁(tr₂[n − 2])]}〉, could theoretically be found using the zeros of the firing time map Δ = g₂〈g₁{g₂[g₁(tr₂[n − 2])]}〉 − tr₂[n] (but see next section).

For perturbations from exact synchrony, the criterion −1 < g*₂(ts*₂)g*₁(ts*₁) < 1 should be applied at ts*₁ = 0 and ts*₂ = g(0). In other words, there are two possible perturbations from synchrony, one in which a given neuron slightly leads it partner and one in which it lags. When it leads, it will receive an input at a phase slightly greater than zero (to the right of zero). The network period is given by g(0); therefore when the neuron lags its partner it will receive an input at a time just less than g(0) [to the left of g(0)]. In some cases a stable locking could not be obtained for the fPRC protocol at a phase of g(0). In that case, there is no way to measure this slope from the left; thus in practice the criterion cannot be applied. However, the absence of a stable locking implies that synchrony in the two-neuron network would be vulnerable to perturbations and probably unstable. The inability to rigorously apply the criterion is of concern in only a narrow region on the border of stability in which the attraction of the stable locking at zero delay compensates for the unstable locking at a delay of one network period.

The choice of a type I PRC for the simulated example was purely phenomenological to match the experimental neuron. Neither the conventional PRC nor the functional PRC methods depend on the shape of the PRC for their implementation, although the shape of the PRC or fPRC determines which modes can exist as well as their stability.

Testing the predictions with networks of adapting model neurons

To determine any advantage of using the fixed-delay protocol, we compared the predictions made with multiple-pulse fPRCs and single-pulse PRCs. For a network consisting of two adapting WB model neurons reciprocally coupled by inhibitory and excitatory synapses, the predictions based on the fPRC protocol (green circles) gives a better agreement with the observations (red crosses) than the prediction by standard PRCs (blue circles) for both inhibitory (Fig. 5C) and excitatory (Fig. 5D) coupling. As indicated in Fig. 5C, in an inhibitory network, antiphase synchrony was perfectly predicted by the fPRC method (green circles) for all the parameters tested. In an excitatory network (Fig. 5D), the nearly synchronous mode, the nearly antiphase mode, and the perfect antiphase mode are correctly predicted by the fPRC method. The standard PRC method correctly predicted for only very limited ranges of small conductance for both inhibition and excitation networks.

With mutual inhibition (Fig. 5C) the antiphase mode is always observed and predicted for the range of values of synaptic conductance. In addition, in-phase synchrony (not shown) is observed in the range of conductance values 0.01 to 0.06 mS/cm², but is predicted from only 0.01 to 0.04 mS/cm² with multiple-pulse fPRCs. This is because at conductance values of 0.05 to 0.06 mS/cm² no stable locking fixed-delay protocol existed at a delay of one network period corresponding to the network period at a delay of zero (see above-cited methods). Therefore synchrony was predicted to be unstable because there was no way to calculate the required slope for a perturbation from synchrony. At conductance values from 0.07 to 0.18 mS/cm², near-leapfrog modes in which the order of firing changes every cycle (Maran and Canavier 2008; Oh and Matveev 2008) were observed to be bistable (Skinner et al. 2005) with the antiphase mode. These leapfrog modes are not shown because they could not be predicted as the fixed points of the map described in the preceding session. The stable lockings under the fixed-delay protocol at the appropriate delays required to measure the g functions used in the maps did not exist. On the other hand, the method based on single-pulse PRCs performed much more poorly, predicting stable synchrony at conductance values from 0.01 to 0.16 mS/cm².

Figure 5D compares the predictions made by fPRC (green circles) and standard PRC (blue circles) for excitation networks. The predictions generated using standard PRC protocol are qualitatively incorrect beyond about 0.08 mS/cm². At smaller values of the synaptic conductance, both methods predict a nearly synchronous mode in which the firing intervals between the two neurons alternate between short and long. Above this value of conductance, both the observed modes and the modes predicted using the fPRC protocol gradually change from near synchrony to near antiphase and finally to perfect antiphase at g_syn = 0.16 mS/cm² in which all firing intervals are equal (see also Fig. 5B). The intrinsic period of the uncoupled oscillators is 60 ms, so the phase resetting greatly shortens the period. For example, the network period is 21.75 ms at g_syn = 0.16 mS/cm². As the conductance is increased beyond 0.3 mS/cm², the network period is reduced to less than a sixth of its uncoupled value and the multiple-pulse protocol no longer accurately predicts the antiphase mode. This is because the gaps in the curves similar to those plotted in Fig. 5B prevent an intersection. Therefore the individual fixed-delay lockings that constitute the observed antiphase mode are unstable and are stabilized by the coupling in a manner that we cannot at this time explain or predict.

History of fixed-delay phase resetting experiments

The phase resetting effects of a single stimulus have been studied extensively (Glass and Mackey 1988; Winfree 1980), but the effects of repeated stimulation at a fixed delay following a marker event such as a spike have been studied much less. An early study (Pinsker and Kandel 1977) used the term “contingent” to describe pulses applied at a fixed delay after a marker event. Stimulation of a presynaptic cell was contingent on the beginning of a burst in the regularly bursting of left upper quadrant bursters in the abdominal ganglion of Aplysia. For some early and late phases, they found an adaptation effect for multiple-pulse PRCs. Early fixed-delay studies using repetitive burst perturbations were performed in the hearts of anesthetized dogs (Levy et al. 1972; Michaels et al. 1984). The cardiac pacemaker was driven to a synchronized frequency with the stimuli. The study of fixed-delay stimulation on spontaneous beating of chick heart cell aggregates (Kunysz et al. 1997), an experimental model of reentrant tachycardia, showed that the complex patterns of bursting behaviors arose from the interaction of phase resetting (thus the excitability) and overdrive suppression, which altered with the amplitude of the stimuli. The authors confirmed the finding with the theoretical studies on a nonlinear model and an ionic model. A study considering respiratory rhythm applications assumed a radial isochronic clock model in which a fixed delay implied that the input was received at a fixed phase, but at variable amplitude such that complicated dynamics ensued (Lewis et al. 1987). The trajectory was not constrained to return to the original limit cycle during the stimulation protocol.

More recently, single-pulse phase resetting curves (Foss and Milton 2000; Foss et al. 1997) were used to determine the existence and stability criteria for phase locking to a repetitive input at a fixed delay that could be much longer than the oscillator period; inputs at fixed delays shorter than the oscillator period were assumed to occur at a fixed phase and to be stable because the phase was not assumed to be able to change during the fixed-delay interval (second-order resetting was ignored, as was adaptation). In no case were the results of the fixed-delay protocol used to predict 1:1 phase-locking network activity.

Relationship of adaptation and PRC types

Previous studies, such as that of Mancilla et al. (2007) and Tsubo et al. (2007), focused on changing the frequency with an applied current, allowing adaptation to stabilize, then measuring a single-pulse PRC. This method is more suited to weak coupling in which the network frequency is assumed not to be very different from the intrinsic frequency of the component neurons. There are various biophysical mechanisms that underlie spike frequency adaptation. The most commonly studied are the calcium-dependent potassium current (I_AHP) and the slow voltage-dependent potassium current (I_M) and, less frequently, sodium-activated potassium currents. As mentioned previously, there are two basic PRC types, type I and type II; they have different shapes and are associated with different synchronization properties (Ermentrout 1996; Rinzel and Ermentrout 1998). Increasing I_M can switch the PRC of hippocampal pyramidal cells (Prescott et al. 2008) from type I to type II. This result was predicted (Prescott and Sejnowski 2008; Prescott et al. 2006) by dynamic analysis of a Morris–Lecar model with both adaptation currents of I_AHP and I_M. Conversely, a decrease in the amount of I_M is predicted to change the shape of the PRC from type II to type I. This switch was demonstrated in cortical pyramidal neurons (Stiefel et al. 2008, 2009) by using a cholinergic neuromodulator to down-regulate I_M. Modeling studies (Ermentrout et al. 2001; Jeong and Gutkin 2007) suggest that in type I neurons I_AHP flattens the PRC at the beginning of the firing cycle. Our methods do not require us to determine how the shape of the single-pulse PRC changes as a function of adaptation or of frequency, but rather depend on the simple and direct measurement of the network period on the time at which an input is received.

In our experiments, we found that the Aplysia neurons that we studied exhibit type I standard PRC under all the values of synaptic conductance that were tested (see Fig. 2 of Preyer and Butera 2005). In the simulations, we were able to match the adaptation observed in Aplysia neurons by using an I_AHP-like current that temporally summates only during spikes and decays in between. The theoretical method based on fPRCs successfully predicted the synchronization of excitatory and inhibitory networks consisting of two type I neurons with an I_AHP-like adaptation current. Our method, however, is general and completely independent of the underlying PRC types.

Limitations

A derivation for stability of a fixed-delay locking in the presence of adaptation was not obtained and we can obtain the recovery interval only as a function of the delay for stable lockings. It seems likely that if both neurons in the network can stably entrain in isolation then the network consisting of those intervals will be stable, but it is also conceivable that the feedback within the network could stabilize some delay values that are not stable in isolation. By analogy, Dror et al. (1999) showed a stable mutual locking of neurons that could not be periodically driven at the same phase they receive input in the mutual locking.

Summary

In our experimental and numerical studies, functional PRCs were generated by delivering a repetitive-pulse perturbation at a fixed-delay time triggered by spike onset at delays ranging from zero to the intrinsic period. Our experimental results on spiking neurons of Aplysia have shown that the fPRCs are distinct from the standard PRC for relatively strong perturbations that induce spike frequency adaptation. The simulated fPRCs for an adapting WB model neuron were in good agreement with experimentally generated fPRCs for strong inhibitory and excitatory perturbations. The predictions of 1:1 phase-locking network activity generated using fPRCs were in good agreement with the observations of phase locking in networks consisting of two adapting neurons coupled with excitatory or inhibitory synapses.

Previously the application of PRC methods has been limited to limit cycle oscillators. However, some components of central pattern generators, such as postinhibitory rebound neurons, do not oscillate when uncoupled from the network. Others may oscillate chaotically (Elson et al. 1999). If these neurons entrain to a periodic input, it is still possible to describe their recovery times as a function of the stimulus time and analyze the network using the methods developed in this study. Therefore the significance is that this method can correctly account for the effect of adaptation on synchronization and has the potential to be applied even to neurons that are not endogenous oscillators but participate in a periodic rhythm.

A limitation of the fPRC method is that a given neuron is assumed to receive only a single input per cycle; thus it cannot be easily generalized to account for a triphasic rhythm in three neurons, for example. However, as currently formulated, the fPRC method can easily be extended to the analysis of global synchrony in an N neuron network using a method analogous to that shown for the conventional PRC (Achuthan and Canavier 2009).

The fPRC method can also be extended to any 1:1 phase locking between two clusters of neurons, again by analogy with published PRC methods (Achuthan and Canavier 2009). Locking between clusters of similar or disparate populations may be widespread in the nervous system. It has been suggested that a synchronous population of cortical neurons locks with a population of thalamic neurons in absence seizures (Perez Velazquez et al. 2007). Similarly, populations of inspiratory neurons in the pre-Botzinger complex and expiratory neurons in the parafacial respiratory group (Janczewski and Feldman 2006; Joseph and Butera 2005; Mellen et al. 2003) have been proposed to phase lock during respiratory rhythmogenesis. Two clusters of entorhinal stellate cells (Pervouchine et al. 2007), with each cluster firing at theta frequency in antiphase with the other cluster, have been invoked to explain a beta peak seen by Cunningham et al. (2004) in slices of entorhinal cortex.