# Spatiotemporal dynamics of neuronal population response in the primary visual cortex

^{a}Department of Mathematics, Ministry of Education-Key Laboratory of Scientific and Engineering Computing and^{b}Institute of Natural Sciences, Shanghai Jiao Tong University, Shanghai 200240, China; and^{c}Courant Institute of Mathematical Sciences and^{d}Center for Neural Science, New York University, New York, NY 10012

See allHide authors and affiliations

Contributed by David W. McLaughlin, May 2, 2013 (sent for review August 11, 2012)

## Abstract

One of the fundamental questions in system neuroscience is how the brain encodes external stimuli in the early sensory cortex. It has been found in experiments that even some simple sensory stimuli can activate large populations of neurons. It is believed that information can be encoded in the spatiotemporal profile of these collective neuronal responses. Here, we use a large-scale computational model of the primary visual cortex (V1) to study the population responses in V1 as observed in experiments in which monkeys performed visual detection tasks. We show that our model can capture very well spatiotemporal activities measured by voltage-sensitive-dye-based optical imaging in V1 of the awake state. In our model, the properties of horizontal long-range connections with NMDA conductance play an important role in the correlated population responses and have strong implications for spatiotemporal coding of neuronal populations. Our computational modeling approach allows us to reveal intrinsic cortical dynamics, separating them from those statistical effects arising from averaging procedures in experiment. For example, in experiments, it was shown that there was a spatially antagonistic center-surround structure in optimal weights in signal detection theory, which was believed to underlie the efficiency of population coding. However, our study shows that this feature is an artifact of data processing.

- lateral long-range connection
- optimal detection theory
- pixel size
- population dynamics
- spatiotemporal patterns

Sensory processing, motor coordination, and higher brain functions are carried out by elaborate networks that are made up of a large pool of neurons. In some primary sensory areas, such as the primary visual cortex (V1), there is a one-to-one mapping of the external visual world onto the cortical sheet and it is represented by a pattern of neuronal population activity (1). Then, visual signals contained in these patterns are decoded further downstream, mediating behaviors (2). Population coding is very common in the brain (3). To understand these early sensory population encoding and decoding processes, it is important to measure the spatiotemporal responses of neuronal populations. Optical imaging techniques have been developed and successfully applied to both anesthetized and awake animals during the last few decades (4). For example, the voltage-sensitive dye (VSD) technique allows for the measurement of target-evoked cortical responses in V1 at high spatial and temporal resolutions (5, 6). This real-time optical imaging tool, combined with single-unit recordings, has opened up new opportunities for the exploration of fundamental mechanisms that underlie cortical development, function, and plasticity in the brain.

Using real-time optical imaging based on VSD, a series of experiments has been carried out to measure the spatiotemporal responses of V1 when monkeys performed visual detection tasks (7⇓⇓–10), and signal detection theory has been applied to analyze optimal neuronal decoding rules in comparison with monkeys’ behavioral detection sensitivity. It was found that (*i*) a small Gabor target of different contrast can elicit a population response described by a Gaussian profile in space, (*ii*) the statistics of responses at each pixel site are Gaussian, (*iii*) long-range (∼4 mm) spatial correlation of responses is nearly independent of contrast, (*iv*) the standard deviation (SD) of responses is nearly constant in time and is also independent of contrast, (*v*) the response onset to the Gabor stimulus is dependent on contrast and the responses begin to rise approximately simultaneously over the entire evoked V1 region but reach their peaks more rapidly at the center of the evoked region, and (*vi*) when the Gabor stimulus is turned off, the related offset dynamics are independent of contrast and they fall simultaneously with the same rate at all locations of the evoked region. In particular, surprisingly, it was found (7) that (*vii*) there is a center-surround structure in the optimal set of weights and the optimal whitening filter derived from the spatiotemporal activity using signal detection theory. Some of the above phenomena, for example the onset/offset dynamics, have been examined through a population-gain control model (10).

However, there still remain some important questions. How could all of the observed phenomena be related dynamically? How do they reflect the underlying properties of V1 cortical dynamics? What intrinsic synaptic and architectural mechanisms could be responsible for these observed spatiotemporal activities? To answer the above questions so as to deepen the understanding of possible functions implied by these phenomena, such as the above antagonistic center-surround structures proposed as properties of an efficient population decoding strategy (7), we use a large-scale, realistic computational neuronal network model of V1 to study the observed spatiotemporal V1 dynamics by VSD optical imaging in awake monkeys. Our model cortical network (MCN) incorporates many physiological facts (e.g., neural coupling architectures) derived from experiments (*Methods*). Here, we model cortical dynamics of an awake state, which differs from the previous study for the anesthetized state (11, 12) in that the cortical network here tends to be more inhibited. (The current work uses neurophysiological parameters for the monkey; see ref 11.) We demonstrate that our model can capture very well all of the above-listed experimental observations in refs. 7⇓⇓–10. Then, using our computational model, we further show that the spatial and temporal scales of NMDA conductance play an important role in the observed highly correlated responses and in the onset–offset dynamical asymmetry between the rising and falling time courses. Therefore, our work suggests a possible functional role of NMDA conductance, that is, as the synaptic mechanism responsible for the spatiotemporal activity of V1, this long-range, long-lasting synaptic coupling may greatly impact population coding as embedded in the spatiotemporal patterns of V1. We emphasize that our computational approach allows us to address how to separate effects mainly caused by experimental methods from those of intrinsic cortical dynamics in the model. For example, as our study shows, some observed phenomena (e.g., contrast independence of variance) could be a consequence of statistical averaging effects arising from data measured with a large bin size (7, 9). In addition, the center-surround structures in the optimal weights and whitening filter are consequences of processing experimental data with such a large bin size, and they are likely not related to the intrinsic functional organization of cortical dynamics.

## Results

To understand possible decoding and encoding strategies embedded in those properties of responses in V1 experiments, we present the spatiotemporal dynamics of our MCN (the architecture is shown in Fig. 1*A*, and *Methods* gives details) in comparison with experimental observations by following precisely the same data processing procedures as used in experiments.

### Spatial Dynamics.

We first discuss the spatial dynamics of voltage responses in both our MCN and experiments. Because the amplitude and spatial spread of the responses characterize neuronal populations that can contribute to the visual detection, we investigate this statistical property first by averaging the response in each site across all target contrasts and over a time interval (∼200 ms) after the onset of a small Gabor target (*Methods*); “site” here refers to an imaging pixel in the VSD experiment (7) and its corresponding group of neurons in MCN. As shown in Fig.1*B*, the size and the spatial structure of the target-evoked responses in our MCN are similar to those in experiments as seen in Fig. 1*C*. Although the retinotopic representation of the small Gabor stimulus on our model cortex through the cortical magnification factor of 3 mm per degree is only ∼0.5mm (8, 10), the evoked responses in both our MCN and V1 extend over an area of several square millimeters, indicating that a large population of V1 neurons carries target-related signals that could be used for visual detection (7, 9, 10). In addition, as in the experiment (Fig. 1*E*), the spatial distribution of our model responses can also be well fitted by a 2D Gaussian, as shown in Fig.1*D*.

To investigate the signal-to-noise characteristics and the reliability of population responses, we measure the variability of responses at each site. As seen in Fig. 1*F*, our model displays the following interesting phenomena similar to those found in experiments as shown in figures 3 a and b in the supplementary text of ref. 7 and figure 4b in ref. 9: For the SD at the site with the maximum signal-to-noise ratio, (*i*) there is no significant difference in the SD between the target-present and target-absent cases and (*ii*) the SD is nearly constant in time in both cases, namely, the variability is largely stimulus-independent. This property of population responses is surprising, in contrast to the variance of the spike count being proportional to the mean for single neurons (13, 14). However, we point out that the stimulus-independent property of the SD in the responses (as shown in Fig. 1*F*) may be related to the choice of the pixel size, which determines the total number of neurons for averaging in each pixel used in optical imaging VSD experiments (7, 9). To see that, we use different pixel sizes (note that smaller size contains fewer neurons and the response is averaged over neurons within a pixel). Then we reexamine the time course of SD at the site with the maximal signal-to-noise ratio. As shown in Fig. 1*G*, the SD of the response is largely stimulus-independent if the pixel size is sufficiently large [e.g., ¼ mm as used in experiments (7)] whereas the SD of the response increases with contrast if the pixel size is small (e.g., 1/32 mm). When the pixel size is large (e.g., ¼ mm), the response averaging is taken over thousands of neurons containing both the stimulus-evoked and spontaneous neuronal responses, as indicated in the experiment (10, 15); the fluctuations of spontaneous activities of neurons dominate in this case, giving rise to a nearly stimulus-independent SD. However, when the pixel size is small (e.g., 1/32 mm), the response is averaged over a small number (∼20) of nearby neurons. In such a case, the neurons used in average for each site tend to share similar properties and are strongly correlated with stimulus, and thus the SD becomes stimulus-dependent. This is also consistent with the previous experimental observation that for single neurons the variance of the spike count during a short interval is proportional to the mean (13, 14), which can be viewed as the extreme limit of our small pixel size.

Following the experiment, we also investigate the distribution of the Z-transformed responses across trials for all pixels. Fig. 1*H* shows that the probability density distribution of the voltage responses in our MCN is (*i*) approximately Gaussian and (*ii*) stimulus-independent. These two features are also observed in experiments as shown in figure 3 c and d in the supplementary text of ref. 7. However, we point out that this Gaussianity may simply arise from the central limit theorem because the response at each site comes from a summation of voltage responses around thousands of neurons, and is further averaged over a time interval (∼200 ms). It seems to be the case because this Gaussianity can be easily obtained in our MCN with no need of fine tuning at all.

Because the spatial correlation of responses could have a large impact on the gain of population decoding (16, 17), we measured the spatial correlation in our MCN, obtained by averaging over trials between the voltage responses at pairs of sites as a function of their distances. The correlation functions under target-present and target-absent conditions are shown in Fig. 1*I* by solid lines, which are in good agreement with the experimental results as shown in figure 2f in ref 7 and figure 3e in the supplementary text of ref. 7. Both our numerical results and experimental observation possess the following features: (*i*) there is no significant difference in the structures of spatial correlation between target-present and target-absent cases and (*ii*) there is a kink in the correlation structure, occurring at ∼500 μm, separating high correlations with a cusp-like shape at the center and a slow decay over long distances. In our model, we demonstrate that these features in the spatial correlation depend on the combined effect of local, isotropic AMPA connections and horizontal long-range connections with NMDA conductance. To show this, we reduce the long-range spatial scale (originally ∼1.5 mm) to that of local AMPA, which is ∼250 μm in the spatial kernels (*Methods*) and find that the slow-decay part in the correlation disappears and only the cusp-like part remains as shown in Fig. 1*I* (dashed lines). Therefore, the slow-decay part of the spatial correlation in our MCN indeed arises from the long-range NMDA synaptic coupling. This result suggests that NMDA horizontal connections may serve as a possible underlying neurophysiological mechanism for the correlation structures as observed in V1 (7, 15).

As seen in Fig. 1 *B*, *F*, *H*, and *I*, our MCN can well capture the neuronal population responses as observed in experiments (7, 9) whose distributions are Gaussian with long-range spatial correlations. In the work of Chen et al. (7), monkeys’ performance of detecting various Gabor stimuli with different contrasts was examined. It addressed how well the spatially correlated population responses in V1 could be used for the detection of visual signals and what is the optimal decoding rule that can take advantage of these correlated structures to maximize detection accuracy (7). Notice that the response in each site is a Gaussian distributed random variable and the covariance matrix, **Σ,** of responses between sites is stimulus-independent, as shown in Fig.1 *F*, *H*, and *I*. Therefore, by using the maximum likelihood method as demonstrated in ref. 7, we can also obtain the optimal pooling rule based on the statistical structures of responses in our MCN as follows. First, we can construct the prior probability density function (pdf) as a 2D Gaussian: *i* = 1 or 2, where **x** corresponds to the observed response, and represent the target-absent and target-present class, respectively, and are the corresponding relative mean responses across trials (the mean response in target-absent trials was subtracted from the responses in target-present trials in our analysis). Second, the decision rule in this two-alternative forced choice can be determined by the relative value between and , where corresponds to the posterior pdf. Finally, using Bayes’ rule, one can define the standard discriminant function . Then, the optimal decision rule is equivalent to the determination of the sign of , where **W** is the optimal set of weights (18), given by . Here, the inverse of covariance can be obtained by with the whitening filter **K**, which can be computed from the power spectrum of the responses (18). We use the averaged responses of neurons over the spatial scale around 0.25 mm in our MCN to represent the measured responses for each site in experiments. The spatial size for data analysis in our MCN is also chosen to be the same as that in experiments (7). Fig. 1 *J* and *K* display the 1D whitening filter (red dashed line in Fig. 1*J*) and the optimal set of weights **W** as computed from the responses in our MCN, respectively. It can be clearly seen that there is a spatially antagonistic center-surround structure in , leading to a corresponding “Mexican-hat”organization in **W**. These features are in agreement with the experimental results as shown in figure 3 c and d in ref. 7. Unlike commonly observed center-surround structures in receptive fields, we would ask why the optimal decision rule for decoding neuronal population responses in the cortex would possess such antagonistic structures. If these structures were truly embedded in the dynamics of the cortex, what would be their underlying functional significance? However, it turns out that these antagonistic structures in both **K** and **W** are data-processing artifacts. This can be shown as follows.

Following exactly the same data-processing procedure as in the experiments (7), the whitening kernel is computed bywhere is the whitening filter in frequency domain and is obtained through the inverse of square root of spectrum . By the Wiener–Khinchin theorem, is the Fourier transform of the correlation function . As shown in figure 2f in ref. 7 and figure 3e in the supplementary text of ref. 7, the measured 1D spatial correlation was smoothed to obtain the best fit to the form where is a fitting parameter. Then, the 2D function is constructed using this fitted exponential under the assumption of isotropy. Following the above procedure, we can calculate analytically and obtain that , if and , if . Hence, , if and , if that is, the analytically exact 1D whitening filter possesses no center-surround structures and the solid black line in Fig. 1*J* demonstrates this point. As a matter of fact, the appearance of the center-surround feature in ref. 7 arises from the finite discretization in the numerical computation of the integration in Eq. **1**. Fig. 1*J* shows the whitening filter computed using different grid sizes. As can be clearly seen, the center-surround structure shrinks, as the grid is refined, toward the correct limit, which has no center-surround structure. Therefore, this center-surround structure is not an intrinsic physiological property in either our MCN or the V1 experiment. It is simply caused by finite discretization in numerics. Using the analytically exact , we can obtain the corresponding analytically exact optimal set of weights **W** as shown in Fig. 1*L*. It can be clearly seen that the Mexican-hat structure disappears in the form of the exact **W**. Therefore, these antagonistic structures seen in numerically computed whitening filter and optimal set of weights are consequences of finite discretization used in data processing and are not any features of intrinsic decoding mechanism of the V1 circuitry. Note that in ref. 7 the results of and **W** are obtained from the analytical form of the exponential function of , which was a fit to experimental data. This fitted was the starting point in their analysis to obtain and **W** and the original experimental data were no longer used from this point on. Therefore, our above analysis is precisely repeating what was carried out in ref. 7 and our mathematical conclusion also holds for their results: The exponential fitting procedure along with finite discretization to handle data would always give rise to antagonistic structures in and **W**. Finally, we point out that, as the pixel is refined, the SD is no longer independent of contrast (Fig. 1*G*); as a result, the linear pooling strategy is no longer optimal for signals obtained by averaging over smaller pixels.

### Temporal Dynamics.

To further understand the architectural and synaptic consequences of long-range horizontal connections, we also investigate temporal properties of target-evoked responses in our MCN. Here, the temporal dynamics are spatially averaged as in refs. 9 and 10. We measure the mean and the SD of the population response in time across trials because they determine the reliability of neural population coding. In our model dynamics, just as observed in experiments (7, 9), the mean response depends on the target contrast, whereas the SD is relatively constant in time and largely stimulus-independent, as seen in Fig. 1*F*.

Fig. 2*A* shows the time course of the normalized responses under different contrast in the central annulus (a circle with a 0.5-mm radius in the center of 2D Gaussian fit as in refs. 9 and 10) in our MCN. Based on the methods in ref. 10, we divide the time course of responses into two parts. The first part (i.e., the rising stage) is the response in the first 200 ms after the stimulus onset; the second part is the falling stage, which is the remaining time course of responses. The response of each stage is fitted with a logistic function as , where (r and f stand for the rising and falling stages, respectively), corresponds to the time when the response reaches half of the peak value, describes the slope of the response, and *P* is the peak value of response in each contrast. Fig. 2*B* shows logistic fits for the dynamics of these two parts with square symbols indicating . Clearly, increases as stimulus contrast is reduced, whereas, is independent of stimulus contrast. Our numerical results qualitatively agree well with the experimental results (10). As in ref. 10, the latency of the rising stage is defined as the time for the fitted response to reach 10% of its peak amplitude. Similarly, the latency of the falling stage is the time at which the response falls from the peak by 10% after the stimulus offset. Fig. 2 *C* and *D* display the latencies and slopes, respectively, for both rising and falling stages of responses in our MCN. It is clear that and are relatively independent of contrast, whereas decreases and increases as the stimulus contrast is enhanced. Again, our results are consistent with experiments as shown in figures 2 c and d in ref. 10.

As will be discussed below, the dynamics in the falling stage in our MCN are dominated by the long-range NMDA lateral interactions because it is the only long-lived (∼80 ms) conductance type after stimulus offset in our model. Therefore, we obtain the same latency and slope in the falling stage of the dynamics for different contrasts. In addition, in experiments (9), there is a large and long-lasting temporal correlations (with a time scale of ∼50 to ∼120 ms) in the VSD responses, and this correlation seems to be similar in both target-absent and target-present trials. We note that this is consistent with the time scale of NMDA conductance. However, the dynamics in the rising stage are dominated by both the local AMPA excitation and the long-range NMDA lateral connections after the onset of stimulus. As the contrast increases, the reaction time for the V1 network dynamics will be accelerated (*Discussion*) similar to what was observed in single neuron studies (19, 20). Therefore, the latency decreases and the slope increases accordingly.

Fig. 2*E* summarizes the spatiotemporal structures of the responses in our MCN, in which we follow exactly the same procedure as in ref. 10. The spatiotemporal response for each contrast is shown as a separate subplot as a function of time; our results are consistent with experimental observation as shown in Fig. 2*F*. These intricate spatiotemporal patterns will be further discussed below.

## Discussion

Now, we discuss the underlying mechanisms for the above spatiotemporal dynamics in our MCN, which may explain the phenomena observed in experiments (7, 9, 10). We first address the temporal properties of dynamics in the rising and falling stages. The governing equation (Eq. **2**) for exponential integrate-and-fire (EIF) neurons can be written as , where is the total conductance and is an effective reversal potential. The initial driving force that makes the neuron’s voltage increase is from the excitatory input current as , where consists of both AMPA and NMDA receptors, that is, . Once the neuron’s voltage has crossed the spike initiation threshold , it will proceed to grow rapidly to produce an action potential due to the internal spike-generating current . Therefore, at the onset, the initial rising time scale of voltage responses is determined by the conductance . This may underlie the experimental observation: After the stimulus onset, the response latency in V1 neurons is strongly dependent on stimulus contrast (21, 22), the response latency being ∼50 ms when the stimulus contrast is relatively low (below 15%) and decreasing to ∼30 ms when the stimulus contrast is high (above 25%) (21, 23), as shown in Fig. 2*F*. Clearly, the spontaneous response before the onset of evoked cortical activity should be independent of contrast, as indicated by the dark blue color in the initial stage of each subplot in Fig. 2*E*. In addition, the AMPA conductance will be larger when the drive from the lateral geniculate nucleus (LGN) is stronger as the contrast increases. However, the larger AMPA conductance will make more neurons fire, which will also induce a strong shunting inhibition, as observed in experiments (24). Therefore, becomes larger under the higher contrast. Although the AMPA and GABA conductances are spatially short-ranged, the long-range NMDA will also be induced through the orientation-specific synaptic couplings from the center region. This is consistent with the long-range spatial correlation observed in both our MCN and experiments (7). All these combined effects give rise to faster response for all locations (due to long-range connections) in the rising stage under higher contrast because the response time scale is inversely proportional to the total conductance (25), hence the observed shortened latency and larger slope for higher contrast in the rising stage in Fig. 2 *C* and *D*. After the stimulus offset, the responses initially stay at the peak value and will eventually decay to spontaneous level. Because the time scales of AMPA and GABA conductances are very short (∼10 ms), they decay rapidly and only the long-lived and long-range NMDA conductance is left to maintain the response of dynamics in the falling stage. Therefore, the dynamics of the falling stages are similar for all locations and for all contrasts, with the decay time scale controlled by the time scale of the NMDA conductance. This decay feature is observed in both our MCN and the experiments as shown in Fig. 2 *E* and *F*.

Next, we address the mechanism underlying the spatial properties for the dynamics of rising and falling stages. The small Gabor stimulus maps to only ∼0.5 mm on the cortex. Because the spatial kernels of all conductance are Gaussian, the amplitude of the conductance decreases as the distance from the center increases, the amplitude of short-range AMPA conductance is higher in the center region of cortex than that in regions far away from the center (>1 mm). Meanwhile, the NMDA conductance persists in both nearby and far away from center regions through the long-range interactions. Taking into account that the response time scale is inversely proportional to the total conductance and spatial decaying of NMDA conductance over long distances, the response rises at a slower rate for larger distances away from the center, as shown in Fig. 2 *E* and *F* (i.e., the rising front in each space-time activity subplot is not perpendicular to the time axis). In contrast, after the stimulus offset, the long-lived and long-range NMDA conductance is the only component remaining to control the decay of the normalized response, as discussed above; therefore, under different contrast, there is a simultaneous onset of decay for all target-evoked locations with the same decay rate as shown in Fig. 2 *E* and *F* (i.e., the falling edge is perpendicular to the time axis). Therefore, it is clear that the long-range and long-lasting NMDA synaptic coupling in our MCN plays a central role in controlling the spatiotemporal responses.

Finally, we point out that (*i*) in the real V1, it is possible that the effect of the NMDA conductance could coexist with some other components that also give rise to similar long-range and long-lasting dynamics; (*ii*) the structured feedback from higher cortical regions or a correlated LGN background firing are not considered in our model; and (*iii*) the background inputs to our network have no spatial or temporal structures and are independent Poisson spike inputs to each neuron with the same constant rate. In addition, our results are robust to the parameters used in our model. For instance, (*i*) the spatial scale of local excitation or inhibition can be chosen from 100 μm to 300 μm, (*ii*) the percentage of NMDA receptor Λ in long-range interactions can be chosen from 40% to 100%, (*iii*) the background firing rate can be from 2 to 15 spikes per second, and (*iv*) the orientation spread projected by the long-range couplings can be chosen from to .

## Methods

Following ref. 12, we model a patch (∼8 × 8 mm) of V1 that reaches the same size as the evoked population responses in cortical areas as observed in VSD experiments (7⇓⇓–10). The patch is a 2D lattice that contains approximately 1 million coupled, conductance-based, EIF neurons. The model cortex consists of 75% excitatory and 25% inhibitory neurons that are uniformly distributed over the lattice. The *i*th neuron with spatial coordinate has a preferred orientation angle and the hypercolumn/pinwheel structure is established by LGN projections (26, 27), as shown in Fig. 1*A* (see ref. 12 for details). The dynamics of the *i*th EIF neuron is governed by

where , characterizing the spike-generating current (28). For physiological values, we use normalized units as in ref. 11: the leaky conductance and voltage . The excitatory and inhibitory reversal potentials are and , respectively. A spike is recorded when reaches the threshold , then is reset to with a refractory period ms. and are the excitatory and inhibitory conductances, respectively. consists of type with a sparse connection and a spatial scale of ∼250 μm (27, 29). consists of both AMPA type and NMDA type (30). The AMPA , where comes from excitatory neurons with isotropic, short-range (≤ 250 μm) sparse connections and from excitatory neurons with orientation-specific long-range (∼1.5 mm) connections. Here Λ denotes the percentage of long-range NMDA receptors (31). models the drive inputs from (*i*) LGN, which is driven by a small Gabor target with size (0.1° to ∼0.3°) and spatial frequency (1.5 to ∼2.5 cycles per degree), as in experiments (7⇓⇓–10), and (*ii*) the background noise in the V1 network. The NMDA conductance is from excitatory neurons with both orientation-specific long-range and isotropic short-range connections. The long-range interactions from both AMPA and NMDA receptors can be written as with , , *H* is a Heaviside function with and (i.e., range from −10° to +10°). The short-range interactions from both AMPA and GABA receptors can be written as with . The spatial kernels are modeled by Gaussian functions with a short-range scale for and long-range scale for (29, 32). The temporal kernels are modeled by *α* functions with rising time constant , , and and decay time constant , , and (27, 33). The LGN input to the *i*th neuron is modeled as a Poisson process with the input rate , where is the visual image and is a spatiotemporal convolution kernel as in ref. 34. is the constant rate of independent background Poisson spike trains to each neuron to maintain a spontaneous firing rate (approximately three spikes per second) per neuron in our MCN. We note that our network is an effective or “lumped” model of V1 because we do not include the detailed laminar structure of V1 in our modeling. This is sufficient for our modeling because the signal of the real-time VSD optical imaging measures the sum of all potential changes in an imaged area, which may reflect information of subthreshold synaptic potentials and dendritic action potentials in neuronal arborizations originating from neurons in all cortical layers whose dendrites reach the superficial layers (35). To obtain quantitative comparison for spatiotemporal dynamics with experiments (7⇓⇓–10), we model the response at the *i*th pixel using the synaptic input, which is proportional to averaged over neurons within a bin of size ∼0.25 mm as used in refs. 7 and 9 and smaller bin size (∼0.04 mm) to study the onset–offset dynamics as used in ref. 10. The sampling frequency (∼100 Hz) and the cortical area (8 × 8 mm) for data analysis in our MCN are also chosen to be the same as those in the experiments.

## Acknowledgments

D.Z. is supported by Shanghai Pujiang Program (Grant 10PJ1406300) and National Science Foundation in China (Grants 11101275 and 91230202); D.C. is supported by National Science Foundation Grant DMS-1009575. D.Z., A.V.R., and D.C. are supported by New York University Abu Dhabi Research Grant G1301.

## Footnotes

- ↵
^{1}To whom correspondence may be addressed. E-mail: cai{at}cims.nyu.edu or david.mclaughlin{at}nyu.edu.

Author contributions: D.Z., A.V.R., D.W.M., and D.C. designed research, performed research, contributed new reagents/analytic tools, analyzed data, and wrote the paper.

The authors declare no conflict of interest.

## References

- ↵
- ↵
- Felleman DJ,
- Van Essen DC

- ↵
- ↵
- Frostig RD

- ↵
- ↵
- ↵
- ↵
- Palmer C,
- Cheng SY,
- Seidemann E

- ↵
- Chen Y,
- Geisler WS,
- Seidemann E

- ↵
- ↵
- McLaughlin D,
- Shapley R,
- Shelley M,
- Wielaard DJ

- ↵
- Rangan AV,
- Cai D,
- McLaughlin DW

- ↵
- ↵
- ↵
- Arieli A,
- Sterkin A,
- Grinvald A,
- Aertsen A

- ↵
- ↵
- ↵
- Duda RO,
- Hart PE,
- Stork DG

- ↵
- Carandini M,
- Heeger DJ

- ↵
- ↵
- Gawne TJ,
- Kjaer TW,
- Richmond BJ

- ↵
- Albrecht DG,
- Geisler WS,
- Frazor RA,
- Crane AM

- ↵
- Carandini M,
- Heeger DJ,
- Movshon JA

- ↵
- Reich DS,
- Mechler F,
- Victor JD

- ↵
- ↵
- ↵
- ↵
- Fourcaud-Trocmé N,
- Hansel D,
- van Vreeswijk C,
- Brunel N

- ↵
- ↵
- ↵
- Myme CIO,
- Sugino K,
- Turrigiano GG,
- Nelson SB

- ↵
- Angelucci A,
- et al.

- ↵
- Koch C

- ↵
- ↵

## Citation Manager Formats

## Article Classifications

- Biological Sciences
- Neuroscience