Fractional differentiation by neocortical pyramidal neurons

Neural systems adapt to changes in stimulus statistics. However, it is not known how stimuli with complex temporal dynamics drive the dynamics of adaptation and the resulting firing rate. For single neurons, it has often been assumed that adaptation has a single time scale. Here, we show that single rat neocortical pyramidal neurons adapt with a time scale that depends on the time scale of changes in stimulus statistics. This multiple time scale adaptation is consistent with fractional order differentiation, such that the neuron’s firing rate is a fractional derivative of slowly varying stimulus parameters. Biophysically, even though neuronal fractional differentiation effectively yields adaptation with many time scales, we find that its implementation requires only a few, properly balanced known adaptive mechanisms. Fractional differentiation provides single neurons with a fundamental and general computation that can contribute to efficient information processing, stimulus anticipation, and frequency independent phase shifts of oscillatory neuronal firing.


Introduction
In response to the sudden change of a stimulus, the firing rate of many neurons undergoes a quick, large change followed by slower adaptation to steady state. First described over 75 years ago 1, this adaptation causes the spike rate to accentuate changes in the input rather than constant values, and in this sense the neuronal response resembles a derivative of the input. In general, neural systems adapt to a change in stimulus statistics, which may serve to maintain maximal information transmission 2-4. Recently many neural systems have been shown to respond adaptively to changes not just in mean input, but also in input variance or contrast 3-10,11,12, higher-order moments 13, and correlation structure 14. In several systems it has also been shown that the time scale of firing rate adaptation is not static but depends on the stimulus, such as on the duration of stimulus steps 4,15,16. It is often unclear to what extent such adaptation of neural systems results from a property of the neural circuit or whether it is a property of single neurons. Neocortical neurons have adaptive processes that modify mean firing rate over time scales ranging from tens of milliseconds to tens of seconds 17-20, and have been implicated in single neuron adaptation to stimulus statistics 7,11. But for single neurons, the time course of firing rate adaptation is often taken to be fixed, or at least independent of the history of change in stimulus statistics.
The mean and variance of fluctuating inputs to cortical neurons can be modulated by changes in input synchrony 21,22, by changes in network excitability due to changes in brain state 23,24, or by stimulus properties. Here, we show that the firing rate in single rat neocortical pyramidal neurons adapts to such changes with multiple time scales, as revealed by corresponding time scales of change in stimulus statistics. We find that these multiple time scales are a signature of a more general property: the neuron's firing rate encodes slowly varying stimulus statistics through fractional order differentiation. This gentle form of differentiation emphasizes change in the stimulus while preserving low frequency information and produces a power law response to sudden stimulus transitions. The fractional differentiation model captures the firing rate response to a stimulus with complex time-varying statistics. Surprisingly, this operation can be biophysically implemented through a balance of only a few known adaptive mechanisms 17-20. Fractional differentiation provides single neurons with a form of adaptation in which no single time scale is preferred; these dynamics may tune the effective time scale of adaptation to the time scales of the stimulus. This operation expands our understanding of the basic processing capabilities of single neurons to the encoding of changes in stimulus statistics over behaviorally relevant time scales in the range of seconds to tens of seconds.
( Figure 1c) and the noisy stimuli with a step change in variance (Figure 1d), thus showing history-dependent adaptation as a property of a single neocortical neuron.
The close to linear stimulus dependence of the adaptation time scale suggests that the observed rate adaptation was generated by a scale-invariant (no preferred time scale), multiple time-scale process rather than a single exponential mechanism. This is consistent with power-law dynamics, as has been observed in some primary sensory receptors 25,26 and references therein. Power-law response to a step is one example of an input-output relationship for a system implementing fractional order differentiation. We tested whether the adaptive mechanisms of single neurons allow them to function as fractional differentiators of the somatic input current.
Fractional differentiation 27, see also supplemental Materials is a linear operation that can be written in the time (t) or frequency (f) domain as: (1) where r(t) is the firing rate response to a time-varying input x(t), h(t) is the fractional differentiating filter in the time domain, and H(f) is the Fourier-transformed filter in the frequency domain ( Figure 2). The constants k and r 0 account for the filter gain and the overall mean firing rate; these constants can be fixed by the mean and standard deviation of the observed firing rate. In the frequency domain, the filter H(f) is (i2πf) α , where α is the order of fractional differentiation. If a neuron functions as a fractional differentiator, the gain of the frequency response will be proportional to (2πf) α and each frequency component of the rate response will have a frequency-independent phase lead of απ/2 with respect to the stimulus (see Appendix and Supplementary Equations). For comparison, a single exponential filter leads to a frequency-dependent gain of 1/(1+4π 2 f 2 τ 2 ) and phase of arctan(−2πfτ). If these conditions hold, the order α of fractional differentiation can be determined from the gain-frequency relationship and phase leads. For example, when x(t) = sin(2πft), then the fractional derivative with order α is: d α x/dt α = (2πf) α sin(2πft + (απ)/2). Finally, fractional differentiation differs from integer order differentiation in that it is nonlocal. While integer order differentiation depends instantaneously on the input function, the result of fractional differentiation depends upon the history of the stimulus, as seen in the slow (power-law) decay of h(t) with time as in Figure 2a.

Single neurons as fractional differentiators
We tested this functional model of the neuronal response through linear analysis. We recorded firing rate responses of L2-3 neurons to injected sine wave (noiseless) currents and sine-modulated noise (constant mean) (Figures 3a and 3b). While previous work has focused on sudden stimulus changes and response characterization for high input frequencies (~1-1000 Hz)28,29, here we focused on low frequency responses (~0.03-1 Hz), for which stimulus periods are large compared to typical inter-spike intervals.
The average gain over all neurons in response to sine wave current and sine wave modulated noise decreased with increasing T in a power-law-like manner (Figure 3c). The exponent α can be determined from the slope of the gain-period relationship on a log-log curve. In response to both sine wave current and sine wave noise, the phase leads were close to frequency-independent (Figure 3d), and the value of α extracted from the phase shift was consistent with that estimated from the gain. The increased phase leads in response to sine wave currents for T = 1-2 sec may reflect a voltage-dependent activation of a conductance with a short time scale, because the mean membrane potential was relatively depolarized relative to that observed during the sine wave noise.
Our data show that mean firing rate responses to stimuli that change slowly relative to the inter-spike intervals are consistent with a model of fractional differentiation. In these experiments, we showed that slowly varying, noiseless currents and modulated noise with constant mean are transformed into firing rate according to the fractional differentiation model. In subsequent experiments, we focused on the responses to sinusoidally modulated noise currents with constant mean.
Because fractional differentiation, Eq. (1), is linear, we tested whether the neuronal response to single sine waves predicts the response to signals constructed as a sum of sine waves. Neurons were given noisy stimuli modulated by an envelope consisting of a sum of sine waves ( Figure 4a). These responses were compared to those observed when each period of sine wave noise was presented individually (Figure 4b). We found that the gain curves and phase leads derived from responses to individual sine waves or a sum of sine waves were not significantly different, suggesting that the firing rate of these neurons depends linearly on the input within the tested ranges of stimulus parameters.
A further test of linearity is given by the responses to square waves that we have already shown. While these curves were fit with single exponentials in order to extract a representative time scale for each curve (Figure 1), the existence of multiple time scales of adaptation implies that this is not an optimal procedure: both faster and slower time scales are apparent during some individual stimulus steps. Therefore, spike responses to square waves were fitted with the multiple time-scale response predicted by fractional differentiation, Eq. (1).
Fractionally differentiating a step function yields a symmetric power-law increase and power-law decrease governed by the parameter α (see Appendix and Supplementary Equations). In Figure 5a, we show the response (solid line) of a neuron to square wave current stimuli, for three different periods, overlaid with the predicted response (dashed line) derived from the linear filter with amplitude and exponent α derived from the sinusoidal experiments as in Figure 3c. In Figure 5b, we show the response to a more complex stimulus in which the standard deviation changed periodically between three different values. This induced upward adaptation when the stimulus stepped from the high to the medium variance, and downward adaptation when the step was from the low to the medium variance. Given appropriate scale factors, these dynamics were well captured by the fractional differentiator model without fitting (see Methods).
Thus, the dynamics of neuronal responses to complex stimuli were well predicted by the differentiating filter. We then least squares fit the square wave responses from Figures 1c and 1d with the response predicted by Eq. (1) using the parameter α. In this case, the best α was found for each neuron in response to square wave current or square noise. For square wave current and noise, α = 0.163 (SD = 0.034) and α = 0.137 (SD = 0.048), respectively (Figure 5c). The difference in error between fitting each curve with a different exponent and fitting all with a single fractional differentiating filter was small. For square wave current, the difference in the mean of the absolute error was 0.18 Hz/bin, while for square wave noise the difference was 0.06 Hz/bin, where errors were slightly smaller with exponential fitting. Differences in squared errors were significant for square wave noise (two-way ANOVA, F(1,3938) = 7.94, P = 0.005), but not for square wave current (two-way ANOVA, F(1,2864) = 2.92, P = 0.087). Thus, scale-invariant fractional differentiation provides an accurate description of the data, with fewer parameters, compared to fitting separate exponential time scales.
We can now compare the value of α as estimated in eight ways, Figure 5c. The overall mean was found to be α = 0.15 (SD = 0.06). These data show that the neuronal firing rate response to time-varying stimulus statistics is compactly approximated as fractional differentiation of order 0.15.

Implementation of fractional differentiation
It remains unclear how power-law dynamics are implemented biologically, although previous theoretical work suggests several possibilities 25,26,30-32. For single neurons, several mechanisms are known to underlie rate adaptation, including slow sodium channel inactivation18 and after-hyperpolarization (AHP) currents17,19,20, which decay with a mixture of time scales33. Because slow AHP (sAHP) currents are distinct from the spikegenerating mechanism, in contrast to slow sodium inactivation, it is relatively straightforward to manipulate these currents experimentally. Here, we experimentally manipulated slow AHP (sAHP) currents, as we expected that these currents contribute to the adaptation time scales as measured by the phase lead; increasing sAHP currents should increase phase leads, while decreasing sAHP currents should decrease phase leads.
First, a sAHP conductance with a decay τ = 3 sec was artificially added to patch-clamped neurons (L2-3 and L5) using spike-triggered dynamic current clamp ( Figure 6a; Supplementary Materials). Without the artificial conductance, phase leads were frequencyindependent, but with the artificial conductance, phase leads increased for T = 4-16 sec with a frequency dependence consistent with the analytical work of Benda and Herz 34. Thus, the addition of a sAHP conductance disturbed the balance of underlying adaptive mechanisms.
Next, again using sharp electrodes (L2-3 neurons), we pharmacologically partially blocked the early sAHP current (τ ~ 1 sec) using the 5-HT 2 receptor agonist α-methyl-5-HT6,7, reducing it by 63% (SD = 13%) (see Figure 6 legend and Supplementary Figs S5 and S6). This manipulation decreased phase leads for T = 4-16 sec (Figure 6b). However, the magnitude of the pharmacologically induced reduction of the sAHP was only weakly correlated with the magnitude of reduction in phase leads (Pearson correlation coefficients = [0.51, 0.21, 0.33] with P = [0.13, 0.55, 0.36] for T = [4,8,16] sec). Based on this and previous data showing that other mechanisms also contribute to spike frequency adaptation 18,19, we do not believe that sAHP currents are solely responsible for the observed phase leads. Rather, it is likely that a combination of mechanisms gives rise to fractional differentiation in these neurons. Here we have shown that manipulating the sAHP current affects phase leads, as expected if this current at least partially underlies multiple time-scale adaptation.
We examined whether a few AHP currents with different time constants are sufficient to generate the observed multiple time-scale adaptation. To examine this, we added to a standard, single-compartment Hodgkin-Huxley model neuron with no slow adaptive processes, first, 2-3 AHP currents and, second, two time scales of slow sodium inactivation (Figure 6a; see Supplemental Materials). AHP currents were added with time constants that roughly span the range of observed sAHP currents in neocortical pyramidal neurons 20. The amplitudes of these currents were precisely balanced with respect to one another to approximate fractional differentiation. For slow sodium inactivation, sodium current dynamics were augmented with slow inactivation gates having two time scales and with kinetics as measured in neocortical neurons18.
These models displayed multiple time scales of adaptation ( Figure 7a) and approximately frequency-independent phase leads (Figure 7b), as was observed in the real neurons ( Figures  1 and 3). The observation that only a few adaptive processes can underlie rate adaptation over a wide range of time scales is consistent with theoretical analysis showing that a single adaptive process affects the neuronal response over a wide range of frequencies34. Together, these results demonstrate that fractional differentiation can be implemented biophysically with known adaptive mechanisms.

Discussion
While single neurons respond to current fluctuations on very short time scales35, the mean firing rate can be considered to encode slowly varying statistical properties of the stimulus 4 or the stimulus envelope 36, acting as an additional channel of information 37. Here, we have shown that the dynamics underlying spike rate adaptation to stimulus steps in single cortical neurons functionally approximate fractional differentiation. This provides a general model for the firing rate response to time-varying stimulus statistics or envelope encoding. Far from simply accommodating to a new state, these adaptive dynamics convey detailed information about stimulus components. Further, this rate response is a linear function of the stimulus mean or envelope. Although extracting a stimulus envelope requires a nonlinear operation36, here the dynamics of spike generation and adaptation provide this operation, while retaining the neuron's ability to represent high frequency fluctuations. We considered responses to variations in mean and variance separately. Further work is required to evaluate the encoding of simultaneous and potentially correlated change in both statistics.
We have focused on the slower time scales that are relevant for time-varying firing rates. One could combine this approach with methods of capturing fast time scale dynamics e.g. 38,39 to give a model that captures both specific stimulus components causing single spikes and slower-varying components relevant for the firing rate. The simplest possibility is that the rate dynamics r(t) act as a multiplicative gain on the fast dynamics of spike generation. Consider for example the case of a slowly time-varying standard deviation σ(t). For a given σ, the occurrence of single spikes is determined by a gain function g σ , which acts on the input after convolution with a filter f σ Decomposing the input into an envelope and a fastvarying component as σ(t)η(t), this multiplicative model for the spike probability is where the rate r is determined by the envelope σ(t), as we have explored here. Such a description may hold when there is a separation of time scales between the short time scale action of f σ and g σ and the slowly time-varying stimulus parameter σ, and allows for changes in the instantaneous coding of inputs due to statistical properties of the stimulus, such as gain normalization4.
Previous work has shown that slow adaptation currents in neocortical neurons have multiple time scales 17,18,20. La Camera et al. 19 found that in vitro data from rat L5 pyramidal neurons could be modeled by a leaky integrate-and-fire (LIF) model neuron modified to include two spike-dependent adaptation currents and one facilitation current. The average fitted time constants for these currents were 48 ms, 5.8 s, and 580 ms, respectively. Due to the facilitating current, this modified LIF model 19 does not generate results consistent with fractional differentiation. However, a modified LIF model with an early sAHP current (τ ~ 1 s)20, instead of a facilitating current, can give results similar to Figure 7. Since the LIF and HH models can give similar results, it is likely that the fractional differentiation model does not depend on the details of spike generation. Other mechanisms of adaptation, such as slow sodium inactivation, cannot be readily implemented in an LIF model.
A fractional differentiator has a power-law response to stimulus steps. Power law responses have been observed in a range of neural systems, from single ion channels 15 to cognitive behavior 32,40. The presence of power law adaptation to a step and fractional differentiation are not synonymous; it is possible to exhibit power law response without the frequencyindependent phase property of the fractional differentiator. Fractional order dynamics have been observed in the vestibular-ocular system 30,41,42 and in the fly motion sensitive neuron H1 43. However, these results were from neural systems rather than single neurons. A range of mechanisms may contribute to fractional order dynamics and power laws, including circuit30 and synaptic 40 mechanisms, geometrical properties of cells 25 and dendrites 42, and the multiple inactivation states of sodium channels 31. Our results show that this computation can be carried out by single cortical neurons and can be implemented by known adaptation mechanisms in a spatially restricted region of the neuron.
Unlike a full derivative, fractional differentiation enhances responses to stimulus change but does not entirely remove information about very low frequency stimulus fluctuations. Retaining information about both the signal magnitude and the signal rate-of-change is a property of fractional order derivatives that leads to dependence on stimulus history. Unlike integer order derivatives, fractional orders are non-local. We found the fractional order to be α ~ 0.15. This gentle differentiation may reflect a requirement for the neural system to encode only a modest amount of rate-of-change information by a single neuron. Increasing the effective value of α to greater than 0.15 could be accomplished by a neural circuit through the sequential adaptation effects of multiple neurons and intervening synaptic dynamics 30,44. We find that the fractional differentiation model for spike rate changes holds over time scales from ~1-30 sec, a range that is behaviorally relevant for fluctuations in input statistics. It has been suggested that fractional order differentiation may be an important property of motor control systems that compensates for fractional order integration dynamics of muscle and tissue, resulting from their viscoelastic properties42. Here, we find that single neurons in the sensory-motor cortex have this property. This kind of rate encoding may also allow neurons to process information efficiently by matching response time scales to input time scales 16,45 or to temporally decorrelate, or whiten, the stimulus envelope46. The power law form of the gain may be related to the typical power-law spectra of natural stimuli 47,48, which would facilitate stimulus whitening. The constant phase shift allows phase to be modulated independently of frequency, which may be significant for the interaction of neocortical neurons with slow oscillatory frequencies in the brain 49. Finally, representing an input's derivative could allow neural circuits to predict future stimuli, similar to a Taylor approximation44.
Our results suggest that fractional differentiation is a fundamental and elementary computation of L2-3 neocortical pyramidal neurons, and provides a general model of the response of adapting neocortical neurons to time-varying stimuli. These intrinsic dynamics provide a substrate for a form of short-term memory, retaining stimulus information over the intermediate time scales of seconds to tens of seconds. While we have focused on L2-3 neurons distributed throughout the sensorimotor cortex, we found similar results among a small sample of L5 neurons. Fractional differentiation may be a general property of adapting neurons in the neocortex. A subset of L5 neurons show minimal slow firing rate adaptation 7, and these cells may not show fractional differentiation. While it may seem surprising that these dynamics have not been discovered previously, it is likely that the experiments that reveal these dynamics were not performed at the slow time scales investigated in the present study. The appearance of these dynamics in single neurons and our ability to disrupt them suggest that adaptation currents and slow sodium inactivation, among other adaptive mechanisms, may be finely tuned to contribute these multiple time scales. This tuning may achieve a balance between providing rate-of-change information and preserving a continuous response in the face of negative feedback.

Preparation of cortical slices
Five-to twelve-week-old Sprague Dawley rats were deeply anesthetized in a chamber filled with 4% isoflurane in oxygen and quickly decapitated. The rostral, caudal, and ventral portions of the brain were removed, and the remaining block was attached to the stage of a Vibratome tissue slicer (TPI, St. Louis, MO) using cyanoacrylate glue (Loctite 404; Loctite, Rocky Hill, CT), and immersed in ice-cold cutting solution containing (in mM) 220 sucrose, 3 KCl, 1 CaCl 2 , 5 MgCl 2 , 26 NaHCO 3 , 1.25 NaH 2 PO 4 , and 10 D-glucose. The cutting solution was bubbled with 95% O 2 /5% CO 2 to maintain pH at 7.4. Coronal slices (300 μm thick) were cut (~ +1 mm to −1 mm relative to Bregma) and transferred to a holding chamber filled with artificial CSF (ACSF) containing (in mM) 125 NaCl, 3 KCl, 2 CaCl 2 , 2 MgCl 2 , 26 NaHCO 3 , 1.25 NaH 2 PO 4 and 20 D-glucose, bubbled with 95% O 2 /5% CO 2 . The temperature of the holding chamber was maintained at 34°C for ~ 60 min and then allowed to cool to room temperature.

Current stimuli
Single neurons were stimulated with time-varying current stimuli in current clamp. Square waves were generated to switch between values of 1 and 2 with a range of oscillation periods. Similarly, standard sine waves were generated to vary between 1 and 2 with a range of periods. A positive holding current was applied to depolarize each neuron to within 10% of rheobase, near the threshold for firing. The basic periodic waveforms were scaled by a factor between 50 and 250 pA chosen individually for each cell such that the neuron fired throughout the stimulus block. In trials using noise stimuli, the simple periodic waveforms were multiplied by exponentially filtered (1 msec time constant) Gaussian white noise to give periodically modulated noise, with a standard deviation given by the square or sine wave. To produce current stimuli, these noise sequences were again DC shifted to near rheobase and scaled as described above. During the noise stimuli, this scaling produced membrane potential fluctuations with SD = ~ 2-8 mV. These stimulus sets were presented in blocks with periods T = [1,2,4,8,16,32] sec, where each period was presented for 96 sec, giving a total stimulus length of 576 sec, following which periods were presented in reverse order. On every individual cycle repetition, the particular noise sequence differed. Finally, modulated noise stimuli were created where the modulation envelope was generated as a sum of sinusoids. Four sine waves with periods T = [4,8,16,32] sec and phase shifts φ = [0, 1, 2, 3] rad were summed, and the result was used to multiply a Gaussian white noise sequence as above. In this case, the block length was 192 sec. Control blocks consisting of individual sine waves with T = [4,8,16,32] sec and block length 96 sec were presented to the same cells. Stimuli were scaled prior to injection, as described above.

Cell acceptance criteria and firing rates
Neurons were considered healthy and stable if their resting membrane potential was < −70 mV, spikes were > 65 mV in height, and resting membrane potential and rheobase (the current at which the neuron first began to spike) did not fluctuate by more than 10% between stimulus blocks (except during pharmacological or dynamic clamp manipulations). Spike times were recorded when the membrane potential crossed −10 mV and reached its peak, and subsequent upward crossings within 2 msec were not considered. Time-varying firing rates were found by taking a histogram of spike times.

Fitting exponentials to firing rates
Firing rate time courses, as in Figures 1A and 1B, were least-squares fit with the twoparameter exponential A exp(−t/τ), where the steady state was defined to be the last three bins of each half-period and 30 bins were used for each period. Cells for the square wave (n = 8) and square noise (n = 11) experiments were included if firing was >1 Hz throughout and firing rates did not show upward adaptation after a step increase. Upward adaptation in the high condition was observed in ~ 10% of neurons, and was usually seen in cells with unusually high firing rates, signs of poor cell health, or unstable recordings. Thus, experiments showing this type of response were excluded from the present analyses. Exponential τ's were not considered if the amplitude of the exponential relaxation, A, was less than 10% of the steady-state firing rate. This criterion excluded 14 of the total 228 measurements of τ.

Fitting sine waves to firing rates
The amplitude and phase lead of the response with respect to the stimulus were found by least-squares fitting the mean firing rate response, determined as previously by taking a histogram of spike times, of each neuron with A sin (2πt/T+φ), with the two parameters A and φ. Identical results were obtained via Fourier transform of the firing rate response, represented as a vector of zeros and ones designating bins in which there was not or was a spike; the dominant frequency amplitude and corresponding phase were easily identified. The gain for each T was found by dividing A by the stimulus scaling factor used for a given neuron, which ranged from 50-250 pA. For the sum of sines experiment, a vector of 0's and 1's to signify no spike or spike, respectively, was created with a sampling rate of 500 Hz; amplitudes and phases were found from the Fourier transform of this vector, with results nearly identical to those from least-squares fitting of the mean firing rate histogram, or from the Fourier transform of the histogram. Neurons were included if at least four blocks (12.8 min) of data in response to the sum of sines stimulus were obtained, in order that phases could be determined with sufficient accuracy. Neurons were excluded if they did not fire throughout T.

Fitting square wave data with fractional differentiator responses
The shapes of the firing rate curves were fit with the response that a fractional differentiator with order α has to a square wave stimulus. The fractional-differentiating filter was found as follows: In the frequency domain the filter was defined as A (i2πf) α for frequencies −960/T to 960/T with Δf=0.01/T. In the time domain, this filter was convolved with step functions of period T with a sample rate of 1920/T, where the non-zero-padded result was taken. Neighboring bins of this finely sampled response were averaged to match the bin number of the firing rate responses. As we simply wished to evaluate the exponent α for the population, we fixed the amplitude A by normalizing each response such that input and output had the same mean and standard deviation. For each experiment, which consisted of six periods, the single best α was found by least squares fitting to the firing rate response. Thus, fitting responses from a single experiment required 3 parameters (α plus normalization of mean and variance) instead of 36 parameters (six parameters for two exponential fits per T) as previously.

Using dynamic clamp to add sAHP currents
When the membrane potential crossed −10 mV, a pulse generator injected a 1 msec negative current pulse into an RC circuit, the voltage V RC of which relaxed with a 3 sec time constant. The RC circuit provided input to a custom-made analog dynamic-clamp circuit that produced an output G asAHP (V membrane − E K ), where G asAHP = kV RC = 0.05 nS/spike and E K =−100 mV. This was summed with the stimulus current.

Biophysical modeling
The single-compartment, conductance- given these kinetics, τ s is voltage-dependent with a peak of ~2300 msec at −50 mV, tapering off to ~500 msec by 0 mV. To approximate AHP currents with time constants of 0.3 sec and 6 sec as in Figure 7, k was set equal to 2/0.3 and 2/6, respectively. Equations were solved numerically using fourth-order Runge-Kutta integration with a fixed time step of 0.05 msec, and spike times were identified as the upward crossing of the voltage trace at −10 mV (resting potential = −65 mV) separated by more than 2 msec. For these data, the mean of the exponentially-filtered (τ = 1 msec) Gaussian white noise stimulus was 5.5 μA/cm 2 , and the SD of the square and sine waves varied between 20 and 32 μA/cm 2 .

Magnitude and phase of a fractional differentiating filter
Because complex numbers can be written in terms of amplitudes and phases, the response R(f) to a signal X(f) = C f exp(iθ f ) through a fractional-differentiating filter may be written as: where it can be seen that the magnitude of R(f) is f α X f and the phase of R(f) has a frequencyindependent phase shift of απ/2. If the signal x(t) is a step function, the resulting response after the step will be a power law decay, as the Fourier transform of a power law is a power law.

Supplementary Material
Refer to Web version on PubMed Central for supplementary material. In response to many repeats of (a) square wave current (top panel) or (b) square wave noise (top panel) with a period T, injected into the soma, the neuron fires action potentials (bottom panels), (c, d) The neuron's mean firing rate adapts both upward and downward (solid lines, bottom panels; 30 bins per period, see also Supplementary Materials). The upward and downward adaptation time courses were each least-squares fit by an exponential of the form A exp (−t/τ) (dashed lines). The upward and downward responses were separately fitted with a time constant τ for each stimulus period T. As T increased, τ similarly increased for both the (e) noiseless current (n = 8) and (f) constant mean, noisy (n = 11) stimuli. Data in (a) and (c) are from one example neuron while data in (b) and (d) are from another.  Sine wave current (a) or sinusoidally modulated noise (b) was injected, and elicited a sinusoidal response (60 bins per period) with a phase lead with respect to the stimulus. Data are from two example neurons, (c) The average gain of the response decreased with period T (n = [8,11] for sine wave currents and sine modulated noise, respectively), similar to a power law with exponent α, as seen from the log-log gain curves; the slope of the best-fit line is −α. (d) In response to sine wave current, the phase lead was frequency-independent Lundstrom  for T = 4-32 sec (n = 8, one-way ANOVA, F(3,28) = 1.63, P = 0.204) with mean 11.9° (α = 0.13), and in response to sinusoidally modulated noise was frequency-independent for T = 1-32 sec (n = 11, one-way ANOVA, F(5,60) = 0.51, P = 0.766) with mean 12.3° (α = 0.14). Error bars represent standard error. Neuronal response to sine wave noise envelopes is similar when periods T = [4,8,16,32] sec are presented simultaneously (with phases ϕ = [0, 1, 2, 3] rad) or individually. (a) The stimulus noise envelope is a sum of four sine waves; the amplitude of the neuronal response can be obtained from the Fourier transform of spike responses. Data is from an example neuron. (b) Gain curves were similar (n = 12, two-way ANOVA, F(1,88) = 0.02, P = 0.902) as were (c) phase leads (F(1,88) = 1.37, P = 0.245). For individually presented sine wave noise, the mean phase lead was 13.1° (α = 0.14), and for sum of sine wave noise the mean was 11.9° (α = 0.13). Error bars represent standard error. In (a), the gain amplitude and α for the predicted response from the fractional differentiator was determined for this neuron using individual sine waves as in Figure 3c; the mean of the predicted response was matched to the firing rate response. In (b) the predicted response was normalized to the mean and SD of the firing rate response. Average firing rates are shown in response to (a) square wave current with T = [4,8,16] sec, and (b) square wave noise with T = 24 sec and three levels of SD = [3,2,1,2] presented for [8,4,8,4] sec, respectively. We also show the response obtained from α = 0.05 and α = 0.25. (c) Bar graph shows the mean (and standard error) for α calculated in 8 ways in response to square and sine wave current (3 data sets) as well as to square and sine wave noise (5 data sets), where the method of finding α is in parentheses. n is the number of α's that contribute to the mean, and is equivalent to the number of contributing neurons. The overall mean was α = 0.15 (SD = 0.06), where none of the obtained values were significantly different (one-sample t-test, P = 0.05). Sine wave current data from T = [ 1 2] were excluded, as these phases were frequency-dependent (Figure 3d). Artificially increasing or pharmacologically blocking sAHP currents led to opposite changes in phase leads, (a) sAHP currents were added to patch-clamped neurons (n = 6) using dynamic clamp (1 msec pulse, τ = 3 sec, ΔG = 0.05 nS/spike, E r = −100 mV). Without the artificial sAHP conductance, phase leads were frequency-independent (one-way ANOVA, F(5,30) = 1.60, P = 0.190). With the artificial sAHP conductance, phase leads were strongly frequency-dependent (F(5,30) = 7.48, P = 0.000). The two conditions were significantly different (two-way ANOVA, F(1,60) = 27.9, P = 0.000), and phase leads increased for T = 4-16 s. (b) Application of α-methyl-5-HT reduced sAHP currents with τ ≈ 1 sec by 63% (SD = 13%), as measured 450-550 msec after a train of 30 spikes at 50 Hz from a pre-train membrane potential set to −65 mV by injection of positive holding current. Drug application (n = 10) altered phase leads (two-way ANOVA, F(1,108) = 9.39, P = 0.003), by reducing them for T = 4-16 sec. Error bars indicate standard error, and asterisks indicate significant differences (paired t-tests, P<0.05).  1c and 1d; as period T is further increased, the time scales become constant (insets), (b) Adding a single, exponential adaptive mechanism can increase phase leads over a wide range of time scales, and two such mechanisms yield a roughly constant phase lead over approximately a 50-fold change in T when the amplitudes of the adaptation mechanisms are properly balanced. Additional currents can lead to responses that vary less over T.