Synaptic Plasticity Can Produce and Enhance Direction Selectivity

The discrimination of the direction of movement of sensory images is critical to the control of many animal behaviors. We propose a parsimonious model of motion processing that generates direction selective responses using short-term synaptic depression and can reproduce salient features of direction selectivity found in a population of neurons in the midbrain of the weakly electric fish Eigenmannia virescens. The model achieves direction selectivity with an elementary Reichardt motion detector: information from spatially separated receptive fields converges onto a neuron via dynamically different pathways. In the model, these differences arise from convergence of information through distinct synapses that either exhibit or do not exhibit short-term synaptic depression—short-term depression produces phase-advances relative to nondepressing synapses. Short-term depression is modeled using two state-variables, a fast process with a time constant on the order of tens to hundreds of milliseconds, and a slow process with a time constant on the order of seconds to tens of seconds. These processes correspond to naturally occurring time constants observed at synapses that exhibit short-term depression. Inclusion of the fast process is sufficient for the generation of temporal disparities that are necessary for direction selectivity in the elementary Reichardt circuit. The addition of the slow process can enhance direction selectivity over time for stimuli that are sustained for periods of seconds or more. Transient (i.e., short-duration) stimuli do not evoke the slow process and therefore do not elicit enhanced direction selectivity. The addition of a sustained global, synchronous oscillation in the gamma frequency range can, however, drive the slow process and enhance direction selectivity to transient stimuli. This enhancement effect does not, however, occur for all combinations of model parameters. The ratio of depressing and nondepressing synapses determines the effects of the addition of the global synchronous oscillation on direction selectivity. These ingredients, short-term depression, spatial convergence, and gamma-band oscillations, are ubiquitous in sensory systems and may be used in Reichardt-style circuits for the generation and enhancement of a variety of biologically relevant spatiotemporal computations.


Introduction
Control of animal behavior often requires the discrimination of the direction of movement of sensory images. In some behaviors, including tracking [1][2][3][4][5] and postural balance [6][7][8][9], animals must determine the direction of sensory slip to generate appropriate compensatory movements to stabilize the sensory image on the receptor array. For other behaviors, including prey capture [10,11], animals must respond to the direction of motion of prey relative to the sensory background. A neural correlate of these functions was described over 40 years ago: Hubel and Wiesel characterized central neurons in the mammalian visual system that exhibited preferential responses to particular directions of movement of sensory images [12][13][14]. Following this discovery, directionselective response properties have been found in a diversity of animal species and across different sensory modalities, from the visual cortex of cats [12] and somatosensory cortex of monkeys [15] to the electrosensory midbrain of weakly electric fish [16,17].
Of particular interest are the direction selective responses of midbrain neurons observed in a species of weakly electric fish, Eigenmannia virescens. These neurons exhibit an unexpected enhancement of direction selectivity by the addition of concomitant naturally occurring sensory oscillations in the gamma frequency range [17]. These oscillations strongly induce short-term synaptic depression in these midbrain neurons, and measures of this depression correlate significantly with direction selectivity [17][18][19]. These correlations suggest that a depression-based mechanism might underlie both the generation of direction selectivity and its enhancement by the addition of gamma-band oscillations. Here we propose a parsimonious model to describe and explore the essential features of this mechanism.
The model is based on a conceptual framework for motion processing known as Reichardt detectors: information from two spatially separated channels with asymmetric temporal properties combine via a nonlinear operation on a downstream (post-synaptic) neuron to produce direction selective responses [20,21] (Figure 1). When supplied with a stimulus moving in the preferred direction, the temporal shift can compensate for the spatial separation, allowing the inputs from the two channels to interact constructively. The temporal phase shift is typically modeled as a pure delay or a low-pass filter.
For direction-selective neurons in V1, Chance et al. [22] propose that the dynamical differences between synapses that exhibit short-term synaptic depression and those that do not may provide a mechanism for generating both the asymmetrical temporal properties and the nonlinear operation required by an elementary Reichardt circuit (the electrosensory midbrain of Eigenmannia also exhibit these requisite ingredients for Reichardt-style selectivity based on shortterm depression [19,23].) For a depressing synapse, the magnitude of the response in the post-synaptic cell decreases during repetitive activation [24][25][26][27]. Short-term synaptic depression involves two dynamic processes with distinct time constants: the faster process, with a time-constant on the order of tens to hundreds of milliseconds, can be attributed to the depletion of the supply of readily releasable synaptic vesicles, while the slower process, with a time-constant of seconds to tens of seconds, can be attributed to the mechanisms for the replenishment of this supply [25,26].
Here we propose a parsimonious model that describes a mechanistic linkage between short-term synaptic depression and direction selectivity, based on a Reichardt-style circuit. We further test the possibility that the enhancement of direction selectivity by concomitant gamma-band oscillations may be mediated by short-term synaptic depression. Global synchronous oscillations in activity may arise endogenously, as occurs in cortical and other circuits [28], or exogenously, as occurs in weakly electric fish from the interaction of the electric fields of nearby conspecifics [29] and the jamming avoidance response [30]. In the model, these oscillations induce depression, which can lead to an enhancement of direction selectivity to moving objects. We systematically explore the effects of variations of biologically relevant parameters of the model and evaluate the results in relation to electrophysiological data from a population of motion-sensitive electrosensory neurons in the midbrain of weakly electric fish [17].

Results
Information from two spatially separated receptive fields converges onto a post-synaptic neuron via dynamically different synapses: one that exhibits short-term synaptic depression and the other that does not ( Figure 1). The spatial separation of the receptive fields combined with the differences in temporal dynamics of the synaptic inputs satisfies the requirements for an elementary Reichardt motion detector. In this model, only the depressing synapse contributes state to the model, which consists of one or two variables whose dynamic evolution is governed by uncoupled and identical (up to parameters) nonlinear ordinary differential equations.

Direction Selectivity Is Mediated by the Fast Process
We have found that the one-state model, which includes only the fast process of short-term synaptic depression, exhibits direction selectivity ( Figure 2A). Since short-term synaptic depression creates a phase advance in the synapse, a moving stimulus that first passes through the nondepressing area leads to a simultaneous arrival (a constructive combination) of signals from both synapses (Figure 2A, blue). Movement in the opposite direction leads to asynchronous arrival of information (Figure 2A, red). These results are similar to those reported previously [22].
The response to the sine-wave grating is nearly identical from cycle to cycle over time. In this case, the time constant of depression is fast relative to the period of the stimulation so that the depressing synapse has sufficient time to return to its initial state during the dark phase of each cycle.

Enhancement Is Mediated by the Slow Process
In the two-state model, which includes both the fast and slow processes associated with short-term synaptic depression, neurons exhibit direction selectivity that enhances from cycle to cycle of a sustained sine wave grating ( Figure 2B). In the first cycle, the response is nearly identical to the response of the one-state model. However, in each subsequent cycle there is a total reduction in the probability of firing and in the total number of spikes for both directions of motion. This reduction in probability of firing is asymptotic.
This overall reduction in firing nonetheless increases the direction index (reported in the caption to Figure 2 and defined in Model) by increasing the relative difference between the responses to the preferred and non-preferred directions of movement. Intracellularly, this enhancement effect will occur as long as the stimulus is maintained even if the depression limits the PSPs so that they do not reach the spiking threshold. In extracellular recordings, however, there is a possibility that the depression could lead to a complete elimination of spiking responses to the moving stimulus.

Responses to Intermittent Stimuli
The sine-wave gratings that we used are sustained stimulisuch stimuli that might arise during image-stabilization tasks. In contrast, many behaviors, such as prey capture, involve spatiotemporally localized, or intermittent stimuli. We have examined the performance of the model to this class of stimuli.

Author Summary
Short-term synaptic plasticity is ubiquitous in brain circuits, but its function in sensorimotor processing remains unclear. We propose a parsimonious model of motion processing using short-term depression to produce directionally selective responses. In the model circuit, information from two spatially separated receptive fields is combined after being asymmetrically processed by synapses that either exhibit short-term synaptic depression or do not. Motion in a preferred direction leads to a constructive interaction between the two channels; motion in the opposite direction does not. The model represents short-term synaptic depression as two processes with distinct time constants. The faster process alone suffices to generate direction selectivity in the circuit. The slow process, in contrast, can enhance direction selectivity to sustained stimuli. Therefore, the slow process mediates a form of attentional shift from alert, where the neuron responds more vigorously, to discriminating, where the neuron responds more selectively with fewer spikes. This explains a previously observed enhancement of direction selectivity in weakly electric fish in the presence of global synchronous gamma-band oscillations. These findings suggest a mechanistic connection between gamma-band oscillations and attention.
Our intermittent stimulus consists of the temporal sequence defined by Equation 2, which is a 1.5 cycle sine-wave pulse. Prior to the arrival of the stimulus, we initialized the system with at least 3 seconds of a spatially homogeneous stimulus of intermediate intensity, which we call 50% grey (see Model). At the arrival of the pulse, the model lies in approximately the same state as it does for the first cycle of the sine grating. As a result, the responses to the first cycle of the grating and to the intermittent stimulus are nearly identical (compare Figures 2B and 3A). For the same parameter values, the responses differ only because the stimuli are subtly different: the sine grating stimulus appears in both receptive fields at the same time, but at different phases, whereas the 1.5 cycle pulse first appears in one receptive field then moves to the other and disappears.
We also tested the model's response to an intermittent stimulus that was initialized not with a uniform background but rather with global synchronous gamma-band oscillations. These sorts of oscillations occur exogenously in groups of weakly electric fish [29] and endogenously in many CNS circuits [28]. In the model, the gamma-band oscillations drive activity simultaneously in both afferents which activates both the fast (0 , D(t) ( 1) and slow (0 , S(t) ( 1) processes in the depressing synapse (see Model).
The response of the model to the moving pulse after 3 seconds of global stimulation compares to its asymptotic response to a persistent sine grating (compare Figures 2B and  3B). The response in this condition is more ''sparse'' than in the grey-initialized condition-the responses are reduced due to the activation of the slow process associated with shortterm synaptic depression. The code is more sparse in that fewer spikes more reliably encode information-the direction of movement. Depending on the values of the parameters, this reduction in spiking can lead to an enhancement of direction selectivity ( Figure 3A versus 3B) or a reduction of the direction selectivity ( Figure 3C versus 3D).
We varied the contributions of the depressing and nondepressing synapses in the model and measured the response to the moving pulse in both the grey initialized and gammaband initialized conditions ( Figure 4A and 4B). Both plots show that direction selectivity reaches a maximum along a ray from the origin corresponding to an optimal ratio of depressing to nondepressing synapses.
To determine under which conditions the gamma-band initialization will lead to an enhancement of direction selectivity, we subtracted the surfaces in Figure 4A and 4B. The maximum enhancement was found to occur along a ray in which the depressing synapses make a greater contribution than the nondepressing synapses ( Figure 4C, magenta). In addition, we found a region in which the combinations of depressing and nondepressing synapses lead to a reduction in The blue channel exhibits short-term synaptic depression (D), which nonlinearly filters the stimulus, advancing the timing of the channel's peak response relative to the timing of the peak of the stimulus (blue dashed lines). The red channel exhibits no depression (ND), which linearly scales the stimulus, leaving the timing of the channel's peak response to coincide with the timing of the stimulus (red dashed lines). Differences in temporal processing interact with a spatial separation of the two channels' receptive fields to produce direction selectivity. A moving stimulus pulse first activating the red channel leads to the summation of coincident peaks in the post-synaptic cell. A stimulus pulse moving in the opposite direction leads to the summation of disparate peaks and hence a weaker response in the post-synaptic cell. (B) A schematic outline depicting the ascending electrosensory system in Apteronotus and Eigenmannia and other Gymnotiform genera. In short, electrosensory information from receptors in the skin (RF 1 and RF 2) project topographically onto the electrosensory lateral line lobe (ELL). Neurons in the ELL in turn project topographically onto neurons in the torus semicircularis (Torus) in the midbrain via the lateral lemniscus (LL). Midbrain afferents include both depressing and nondepressing synapses that converge on to individual neurons: this convergence of information meets the requirements for the proposed elementary Reichardt motion detector. We used three categories of stimuli: global stimuli (social signals that stimulate the entire sensory surface simultaneously), a localized moving bar (shown), and a larger moving sinewave grating. doi:10.1371/journal.pcbi.0040032.g001 direction selectivity with the addition of the gamma-band oscillations ( Figure 4C, green).

Relations between Measures of Direction Selectivity and Post-Synaptic Potential Depression
In intracellular recordings of midbrain neurons in Eigenmannia, the addition of an exogenous gamma-band oscillation resulted in an enhancement of direction selectivity to a moving bar stimulus [17]. In many neurons a correlate of the activity of inputs that experience short-term synaptic depression was observed: the amplitude of post-synaptic potentials (PSPs) declined on a cycle-by-cycle basis to a sustained gamma-band oscillation [17] ( Figure 5). This ''PSP depression'' has been shown likely to result from short-term synaptic depression and not other mechanisms [18,19]. We tested the model with identical stimuli and made an identical measurement of ''PSP depression'' [18]. PSP depression is the magnitude of the decline in amplitude of PSPs measured at or near the soma, and is therefore a sum of the synaptic activity, including both depressing and nondepressing synaptic inputs, to the neuron ( Figure 5).
In the model, PSP depression was strongest where the ratio of depressing synapses to nondepressing synapses was high, but the total number of synapses was low ( Figure 5). Surprisingly, adding more depressing synapses actually decreased the measure of short-term synaptic depression. This measure consisted of a ratio of the response (maximum depolarization above resting potential) to the first cycle of global synchronous stimulation to the average responses to later cycles, after the transient had decayed. As more depressing synapses were added, both the numerator and the denominator of this ratio increased, but they did so in a way such that the value of this ratio decreased.
In Eigenmannia, strong correlations were observed between the magnitude of PSP depression measured in each neuron and the magnitude of direction selectivity to the moving bar in both grey and gamma-band initialized conditions ( Figure  6A). We tested whether any simple set of parameters in the model could reproduce these relations.
We considered four hypothetical distributions, asking if each reproduced the qualitative observations drawn from the sample of neurons within the electrosensory midbrain of Eigenmannia virescens [17]. The qualitative observations were that the measures of PSP depression and direction selectivity were positively correlated in both grey and gamma-band initialized conditions, and that direction selectivity was increased by the addition of the gamma-band oscillation ( Figure 6A).
We tested hypothesized distributions supported on onedimensional restrictions of the parameter space in which 14 Cyan and magenta shows firing rate, F(t) (beneath V ' (t), and timing of action potentials determined from F(t) by a Poisson process (above). Blue and cyan (top subpanel in both A and B): preferred direction (r ¼ 1); red and magenta (bottom subpanel): non-preferred direction (r ¼À1). Small numbers are the expected action potentials in each cycle. Black dotted line indicates threshold for action potentials. In the one-state model, the direction index (defined in Model) is constant across all cycles at 0.48. In the two-state model, the direction index increases from 0.49 to 1-the two-state model becomes more direction selective with time. Other parameters of model are For convenience, we describe our parameter space restrictions in terms of numbers of synapses (with fixed synaptic weights, see Model). Distribution 1 assumes that the total number of synapses remains constant (80), but the ratio of depressing to nondepressing synapses varies. Distribution 2 assumes the number of depressing synapses remains constant (80) but the number of nondepressing synapses varies. Distribution 3 assumes the ratio of depressing to nondepressing synapses remains constant (5/3) but the total number varies. Finally, distribution 4 assumes the number of nondepressing synapses remains constant (12) but the number of depressing synapses varies. The restrictions of parameter space associated with these distributions are plotted in Figure 6B.
The resulting relationships between the measures of PSP depression and direction selectivity for each distribution are plotted in Figure 6C-6F. Distributions 6E and 6F match the three qualitative features seen in the population of neurons observed in the midbrain of Eigenmannia.

Discussion
We examined the roles of short-term synaptic depression in the generation and enhancement of direction-selective responses. A one-state model that includes the fast process (hundreds of milliseconds) associated with short-term synaptic depression can produce direction selectivity in an elementary Reichardt motion detector. The addition of a second state, the slow process of short-term depression (seconds to tens of seconds), can lead to an enhancement of direction selectivity. This enhancement is a form of sparsification: fewer spikes more accurately encode direction of movement.

Modulation of Spatiotemporal Processing via the Slow Process
In our model, the activation of the slower process of shortterm depression is quantified by the state variable S(t) that depends upon the stimulus history. If there has been little recent stimulation (recent with respect to the time constant, which in this model is set to 3 seconds), then the neuron resides in a state in which it responds vigorously to stimuli that are moving in any direction.
On the other hand, if there has been recent stimulation, the depressing synapses will be depressed, and as a result the neuron will respond less vigorously. As shown in Figure 4C, the ratio of the contributions of depressing and nondepressing synapses determines whether the neuron will be more or less directionally selective in the depressed state that results from recent stimulation. In a population of midbrain neurons in Eigenmannia, recent stimulation leads to an enhancement of direction selectivity [17].
Recent stimulation shifts depressing synapses from a highly responsive state to a more depressed state. The difference between these two states may correspond to vigilance and focus (at least in Eigenmannia) and can be seen as a form of attention. In this way, the current value of the slow variable S(t) corresponds to an attentional state associated with the neuron. Indeed, the state of the slow process could be critical to specific computations. Although we only explored the responses to moving stimuli, these results could be applied more broadly to any computation in the brain that involves the temporal comparison of information that converges from independent pathways.
Indeed, we have shown that any stimulus that activates the slow process can lead to a shift in the computational properties of elementary Reichardt circuits. As a result, any change in the activity patterns, whether they be stimulusdriven or endogenous, could affect neural computations through similar mechanisms. This feature may provide an opportunity for animals to use behavior to modulate computations in the brain.

Behavioral Modulation of Direction Selectivity
Animals may modulate the state of their synapses and hence the degree of direction selectivity in central neurons using behavior. This form of behavioral modulation requires that: 1) the behavior generates patterns of activity that elicit short-term depression in downstream neurons and 2) that these patterns of activity do not interfere with the motion processing.
Evidence for the behavioral modulation of directionselectivity is seen in the Jamming Avoidance Response (JAR) of weakly electric fish. In the JAR, the electric fields of fish that are within about a meter of each other interact to The responses in each plot are aligned to the peak of the PSP from the depressing synapse. This facilitates the comparison between the preferred and non-preferred directions of movement, but as a result, the onsets of the moving stimuli are not aligned. In the preferred direction, the object passes first through the nondepressing part of the receptive field and then enters the depressing part. Parameters of the model for this and all remaining figures are the same as listed in Figure 2 except for the synaptic factors, which vary. For A and B (magenta): g n ¼ 0.2, g d ¼ 8. For C and D (green): g n ¼ 0.6, g d ¼ 4. The stimulus parameters used throughout for the pulse stimulus are f g ¼ 20 Hz, A g ¼ 0.7, or A g ¼ 0 (respectively, with or without initializing global synchronous oscillations), t 0 ¼À1,000 msec, t s ¼ 0 msec, t f ¼ 5,500 msec, u ¼ 1358; r, p 0 , and p 1 differed depending on direction of movement: produce oscillations in electroreceptor activity across the entire body [30]. In this way, the stimulation leads to global, synchronous activity across the receptor array. In the wild, the frequencies of these oscillations are most commonly in the gamma frequency band, from 20 to 80 Hz [29]. Laboratory experiments have shown that lower frequency oscillations, below approximately 8 Hz, impair the perception of moving objects [31][32][33]. These gamma oscillations are encoded by electroreceptors and propagate through the ascending electrosensory system. In the midbrain, these oscillations match the stimulation frequencies that best elicit short-term depression [18,19]. As a result, the ongoing oscillations that occur in social situations dramatically modulate short-term depression, leading to an enhancement of direction selectivity to intermittent stimuli [17].
Another more general example of behavioral modulation of temporal processing in plasticity-based elementary Reichardt circuits could include movement-induced self-stimulation of sensory receptors. For example, were an electric fish to remain motionless in a tube, the slow process will be in its initial state for midbrain electrosensory neurons, whereas if the fish were to move back and forth within the tube, the slow process may be activated. The neurons would be more responsive and less selective while the animal remained motionless and would be less responsive but more selective when the animal was moving relative to nearby objects.

Intrinsic Modulation of Spatiotemporal Computations
If behavior can be used to generate patterns of brain activity that change the state of depressing synapses to alter spatiotemporal computations in the brain, then intrinsic activity in brain circuits could possibly have the same effect. The intrinsic activity would have to have two properties: 1) the patterns of activity must elicit short-term depression in downstream neurons and 2) these patterns of activity must not interfere with the spatiotemporal computation.
There are numerous examples of endogenous oscillations that occur at all levels of the CNS [28]. If the output of these endogenous synchronous oscillations converge on neurons that perform spatiotemporal computations, such as motion processing, through synapses that experience short-term synaptic depression, then the endogenous oscillations could have the same effects on computation that have been reported in Eigenmannia. The correlation between attentional processes and the emergence of gamma-band oscillations in cortical and other circuits may support this idea. Perhaps the gamma-band synchronous activity shifts elementary Reichardt circuits from a more responsive but less selective state to a less responsive but more selective state.

Parameter Distributions in Populations of Neurons
The population of neurons observed in the midbrain of Eigenmannia showed a positive correlation between a measure of short-term synaptic depression (PSP depression) and direction selectivity. We tested the hypothesis that this relation follows from the proposed elementary Reichardt circuit that uses short-term synaptic depression. The model, however, clearly demonstrated that this relation is but one possible outcome. Indeed, without assumptions constraining the distribution of parameter values, the model does not make any specific prediction about the relationship between PSP depression and direction selectivity in populations of neurons.
Thus, the relationship previously observed in the midbrain of weakly electric fish [17] is likely associated with functional constraints beyond the elementary Reichardt circuit. These functional constraints may be related to the control of specific behaviors. For example, in tracking behavior [4,5,34], fish make compensatory movements to stabilize an image on the sensory array. Future studies will determine, via neural system identification, how a population of direction-selective neurons may encode the sensorimotor transfer function inferred from behavioral performance [5].

Dynamic Receptive Field Structure
A key feature of the model is the convergence of information from spatially separated locations on the receptor array. The receptive field of a neuron that uses the mechanisms associated with this model should be composed of regions that differ in relation to short-term synaptic plasticity: the response of the neuron to stationary, highly localized stimuli at different locations should show differences in measures of short-term depression.
In the simplest case, the receptive field has two regions, e.g., caudal and rostral. If the neuron exhibits little or no shortterm depression when the stimulus is in the caudal region and strong depression to local stimulation in the rostral region, then the model predicts that the neuron will respond more strongly to caudal-to-rostral movement than to rostral-tocaudal movement [23]. Preliminary evidence obtained from midbrain neurons in weakly electric fish are consistent with this hypothesis: gamma-band stimulation in subregions of the receptive fields of midbrain neurons elicit different levels of short-term synaptic depression (personal observations). Nevertheless, one could envision far more complicated receptive field structure leading to direction selectivity using similar mechanisms.

Quantitative Validation of the Model
By posing a parsimonious model describing the transformation of the spatiotemporal stimulus into the neural response, we have taken the first step in a rigorous identification of the underlying neural system. The steps remaining include (for each neuron) validating or falsifying the model, estimating its parameters and comparing our model to alternative models. The field of system identification offers systematic and rigorous approaches to these remaining problems. For example, the stochastic model that includes action potential timing can be used to determine the likelihood that an experimentally observed train of action potentials was generated using mechanisms captured by the model. Systematic exploration of the parameters can be made to achieve the maximum likelihood estimates [35,36]. This procedure can be repeated for each neuron studied in the population to get more rigorous estimates of the distribution of parameter values.

Model
The structure of the model is generic: it does not incorporate any specific features from particular animal systems. Nevertheless, the basic structure of the model is inspired by a model for direction selectivity in V1 neurons by Chance et al. [22]. Further, stimuli (including moving objects and global synchronous gamma-band oscillations) approximate those used in electrophysiological experiments in weakly electric fishes, such as Apteronotus leptorhynchus and Eigenmannia virescens [37]. The parameters of the model, including the time constants for the fast and slow processes of depression, are similar to those observed in mid-brain neurons in Eigenmannia [18].
Electrosensory information from receptors in the skin project topographically onto the electrosensory lateral line lobe, which in turn projects onto the torus semicircularis in the midbrain ( Figure 1B). Midbrain afferents include both depressing and nondepressing synapses that converge onto individual midbrain neurons. We use three categories of stimuli global stimuli (social signals that stimulate the entire sensory surface simultaneously), a localized moving bar, and a larger moving sinewave grating.
Our computational model explains how direction selectivity arises from known features of midbrain neurons of the weakly electric fish Eigenmannia. Our model reproduces and explains the surprising experimental result that global synchronous electrosensory oscillations experienced prior to a local moving stimulus enhance direction selectivity of the local stimulus in weakly electric fish [17]. Moreover the model captures the diversity of the directionally selective neurons in the midbrain with respect to the observed correlation between direction selectivity and a measure of short-term synaptic depression.
Our model has been adapted, and substantially simplified, from a previously published model of direction selectivity [22]. The Chance et al. model incorporates a large number of dynamic variables and fixed parameters, many of which, we show, are not needed to reproduce the phenomena. We have significantly simplified the model, capturing the relevant features with a minimal number of dynamic variables (two), and substantially fewer fixed parameters (16), making the model amenable for the purpose of system identification [38]. We have provided our model code, written in MATLAB, in Protocol S1.
We model the response of a neuron to a spatiotemporal stimulus as a cascade of four elements: the afferents, the synapses, the synaptic conductance, and the cell membrane. In our simplest formulation of the model, only the depressing synapses have state, i.e., dynamics that depend on state variables that ''remember'' the stimulus history. State variables ''remember'' the stimulus history in the sense that differential equations integrating the element's input determine how they evolve. The other elements in the cascade are all memoryless, i.e., do not depend on any history-dependent quantities. For such elements, the input at every moment completely determines their output at that moment.
Our minimal instantiation utilizes a single depressing synapse (which has a memory in that it exhibits activitydependent reduction in efficacy) and a single nondepressing synapse (memoryless). Based on neurophysiological observations, the short-term synaptic depression dynamics involve two time constants, corresponding to faster and slower depression processes, each involving a single state variable. Chance et al. [22] attribute direction selectivity to the faster process of depression and propose that contrast adaptation might involve the slower process of depression. As evidence, they note that both contrast adaptation and the slower process of depression exhibit similar time scales. We use the more general term sparsification, in lieu of contrast adaptation, because the slow process is driven by any persistent stimulus that causes depression, not just to changes in contrast of the scene. We will build on Chance's observations by elucidating another role for activation of the slower process: enhancement of direction selectivity.
The input stimulus. The input to the system is the stimulus intensity as a function of time t and position x: I(t,x). We consider only four classes of stimuli or temporal sequences thereof: constant functions for all x and t; global synchronous oscillations in which I(t,x) varies sinusoidally with time, independent of x; moving sine gratings, used by Chance et al. to test their own model; and a moving sinusoidal pulse (with a 1.5 cycle period), analogous to a moving bar. Following Chance et al. [22] and others we normalize the intensity function so that I(t,x) 2 [À1,1] with I(t,x) ¼ 0 representing 50% grey.
The temporal sequence we use for a moving pulse stimulus is defined (for t 0 t t f ) as follows: Iðt; In Equation 1, k is the spatial wavelength that sets the unit for space: without loss of generality, we say k ¼ 3608. Likewise f is the temporal frequency that determines the speed of the motion of the pulse. In this case, f determines more than just the scaling of time, because f interacts with the time constants of the model. The parameter u determines the phase of the pulse in units of k, in our case degrees. Finally r ¼ 61 determines the direction of motion of the stimulus.
In Equation 2, the stimulus parameters p 0 and p 1 determine the moving pulse boundaries. Setting p 0 À p 1 ¼ 3p gives the pulse a period of 1.5 cycles. The sign r appears in the conditions on the right-hand side of Equation 2 to reverse the direction of inequalities, needed because p decreases with time for a positively directed stimulus (r ¼ 1) and increases with time for a stimulus in the opposite direction (r ¼À1). We choose p 0 and p 1 so that sin(p 0 ) ¼ sin(p 1 ) ¼ 0, cos(p 0 ) ¼ r, and cos(p 1 ) ¼ Àr to determine a time-symmetric pulse that darkens from 50% grey on its boundaries.
The parameters t 0 and t f are, respectively, the initial and final times of the stimulus. The stimulus parameter t s is a switching time between stimuli in the temporal sequence. At time t s the moving pulse appears, usually outside, but moving toward, the model cell's receptive field. Ahead of the pulse (and also appearing at time t s ) lies initializing input. Behind the pulse lies 50% grey. The initializing input consists of either global synchronous oscillations in the gamma band (A c . 0, 20 Hz f c 80 Hz) or 50% grey (A c ¼ 0).
To compare our simulated responses with experimental data we include a preinitialization phase with 50% grey (t 0 t , t s ), prior to the gamma oscillations. The preinitialization phase allows a meaningful quantification (independent of the model's arbitrary initial state) of short-term synaptic depression. This quantification compares the model's response to the first cycle of global synchronous oscillations following preinitialization to the response to a later cycle after the gamma oscillations have affected the state of the model. Equations 1 and 2 can, for the appropriate choice of stimulus parameters, also define all the other stimuli we employ (see figure captions). For example, the persistent sine grating (with a initialization before first cycle) can be defined by setting p 0 ¼ ' and p 1 ¼ À'. Finally, note that we do not count the stimulus parameters as part of the 16 fixed parameters of the model neuron, because they determine properties of the stimulus rather than the model cell.
The afferents. The afferents, the first element in our cascade, transform the input signal, I(t,x), into the firing rates, R d (t) and R n (t), of two pre-synaptic cells in the depressing and nondepressing channels, respectively. The receptive fields of these two cells are spatially separated, with centers at x d and x n . Following Chance et al. [22] receptive fields are separated by 908 relative to the stimulus wavelength k to maximize the model's response, however, the model is robust to the receptive field spacing relative to stimulus wavelength.
In the interest of model parsimony, we have replaced the kernel of the spatiotemporal integration from Chance et al. [22] with a delta function, rendering the afferents memoryless and confining the receptive fields to a single point. The input/ output transformation of the afferents is given (for I 2 fd,ng) by the rectified linear equation The parameters R b and R c are, respectively, the baseline firing rate and the contrast-dependent rate factor.
The synapses. These afferent firing rates (R d (t), R n (t)) serve as input to the second element, the synapses, consisting of parallel depressing and nondepressing channels. In the case of the nondepressing synapse, inputs are trivially passed to the output: R n (t) 7 ! R n (t). The depressing synapse, on the other hand, is the only component of our model with state. Its input/output relationship is given by the transformation R d (t) 7 ! (R d (t),D(t),S(t)) where (D(t),S(t)) is the synapse state which evolves according to Here s D and s S are the fast and slow time constants of depression and d and s are the fast and slow depression strengths, taken to be between 0 and 1. In Chance's model, D(t) and S(t) are reduced by the factors d and s, respectively, subsequent to an action potential arriving at the corresponding depressing synapse, as timed by a Poisson process with inhomogeneous rate R d (t). Equations 4 and 5 result from averaging (taking the expected values of) the resulting stochastic differential equations, eliminating the randomness in the model. In our model, each depressing synapse is governed by two uncoupled, deterministic, nonlinear differential equations. Setting the depression strength s (or d) to 1 removes the corresponding state (after some input-independent transient); setting both strengths to 1 renders the synapse memoryless and equivalent to the nondepressing synapse.
The synaptic conductance. The third element of the cascade is the model of the synaptic membrane conductance: In the corresponding element of the Chance model, G E (t) is a state variable whose evolution is given by a differential equation with a time constant s E ¼ 2 msec, two orders of magnitude faster than the time constants s D and s S of our model. In Equation 6 we replace Chance's dynamic state variable with its asymptotic value, removing the memory of the element. This manipulation has only a small effect on the overall dynamics, a difference that is not important for reproducing the phenomena considered here.
We call the coefficients c d and c n the synaptic factors, depressing and nondepressing respectively. The term ''synaptic factors'', in place of the more common ''synaptic weights'', suggests a second more biologically relevant interpretation of the model. In this interpretation, each synapse in fact represents a synapse class, encompassing the contribution of a population of similar synapses that project onto the same midbrain neuron. The synapses within each class are identical in the sense that they have identical properties, receive identical input, and maintain an identical state. Under this interpretation, the synaptic factor equals the product of the individual synaptic weights and the number of synapses in the respective class. By covarying these two new parameters, class population and synaptic factor, the model demonstrates the diversity in the population of toral neurons that was observed by Ramcharitar et al. [17].
The membrane: Membrane and action potentials. We consider three alternatives for modeling the cell membrane. The first candidate is the classical leaky-integrate-and-fire mechanism where V 0 is the resting potential, V E is the synaptic reversal potential, and s m is the membrane time constant: Equation 7 applies to the intervals between action potentials. We say the model fires an action potential when the membrane potential reaches a threshold potential V t . When an action potential occurs, the membrane potential discontinuously resets to V r , then again evolves according to Equation 7. By our choice of parameters V 0 , V r , V t , V E , the synaptic current is excitatory.
Notice that the above formulation requires the membrane potential to be a state of the system. However, the time constant of the membrane is roughly an order of magnitude less than the time constant for the faster depression process, so we may again approximate this third state algebraically. This alternative to the leaky-integrate-and-fire mechanism replaces the membrane potential with its instantaneous asymptotic value given the present state of the conductance. This calculation averages the reversal potentials for the leak current (i.e., the resting potential) and the synaptic current: This approximation eliminates the reset dynamics to avoid chattering between V ' (t) and V r when V ' (t) . V t , limiting its use to modeling membrane potential with action potentials blocked.
Whereas the second alternative eliminates action potentials from the model, the third alternative predicts the action potential firing rate during these intervals. This firing rate can be fed to an inhomogeneous Poisson process to predict the timing of the action potentials within the periods of activity, or, more importantly, to determine a statistical model of the timing, useful for identifying the system from extracellular data [35,36]. Alternatives two and three can be combined in a model without adding an additional state variable to the system.
We say that the instantaneous firing rate of the neuron, assuming V ' (t) . V t , equals the reciprocal of the time required for the membrane potential to reach the action potential threshold, V t , from the reset potential, V r , assuming that the value of the synaptic conductance remains fixed at its present value, G E (t). If V ' (t) , V t , such a traversal cannot happen and we say the firing rate is 0. Because the membrane dynamics, given by Equation 7, are linear, the firing rate, as we have defined it, is an exponential function that can be calculated in closed form. Nevertheless, we find it helpful to make one further simplifying approximation. We assume that the rate of change of the membrane potential during the interval between reset and firing is approximately constant and equal to its value at threshold V t . The closer V r lies to V t , the better this approximation. With this simplifying assumption, our firing rate calculation reduces to a rectified linear algebraic equation: where R s ¼ V E À V t s m ðV t À V r Þ ; ð10Þ Here R s and G f are constants. The quantity F 0 (t) represents the firing rate of the cell that is unbounded by a refractory period. Biological neurons have such bounds. For example, a refractory period of s R puts an upper limit of 1/s R on the firing rate of the cell. Without such a limit the firing rate of the model in response to certain stimuli can grow much larger than is biologically plausible. Fortunately by assuming the following saturation in firing rate, we can incorporate a refractory period into our model cell, without adding state to the membrane: where F(t) is the firing rate output of the cell, bounded by the refractory period s R . Note that if the firing rate F(t) determines the timing of action potentials by a Poisson process then it remains possible, though not particularly likely, that the cell will fire two (or more) action potentials within any given time interval of duration of s R . This unlikely possibility will not prevent identifying the system through fitting parameters of the model to experimental data.
Measures of the response. Finally, we consider several measures to quantify the response of the model to various stimuli. For the sine wave grating we use the direction index as a function of stimulus cycle number. Our calculation of direction index involves the expected number of action potentials for the jth cycle in the preferred (P j ) and nonpreferred (N j ) directions, calculated by integrating the firing rate over the corresponding period. If we assume 0 N j P j and P j 6 ¼ 0 (i.e., that the preferred direction has been correctly identified) then the following equation gives the direction index: Note that the quantities P j and N j can be defined in other ways, but as long as 0 Nj Pj and Pj 6 ¼ 0, the direction index lies between 0 and 1, reaches its maximum when N j ¼ 0, and its minimum when N j ¼ P j .
While the direction index remains the most common way to quantify direction selectivity, other quantifications can be useful. Ramcharitar et al. [17] used an unbounded measure of direction selectivity (referred to as magnitude of direction selectivity) to demonstrate a nearly linear correlation with a similar measure of short-term synaptic depression (called magnitude of PSP depression). To compare our results with these data, and because their experimental paradigm corresponds to the moving pulse stimulus for our model, we use this second measure to quantify the response of our model to the pulse stimulus. This second measure is the ratio (converted to dB) of the height of the depolarization of the membrane above the resting potential in the preferred direction to the same height in the non-preferred direction. A similar calculation quantifies short-term synaptic depression with the ratio of the response to the first cycle of global synchronous oscillations to a later cycle after the transient has decayed (PSP depression).