## Abstract

In the primary visual cortex (V1), Simple and Complex receptive fields (RFs) are usually characterized on the basis of the linearity of the cell spiking response to stimuli of opposite contrast. Whether or not this classification reflects a functional dichotomy in the synaptic inputs to Simple and Complex cells is still an open issue. Here we combined intracellular membrane potential recordings in cat V1 with 2D dense noise stimulation to decompose the Simple-like and Complex-like components of the subthreshold RF into a parallel set of functionally distinct subunits. Results show that both Simple and Complex RFs exhibit a remarkable diversity of excitatory and inhibitory Complex-like contributions, which differ in orientation and spatial frequency selectivity from the linear RF, even in layer 4 and layer 6 Simple cells. We further show that the diversity of Complex-like contributions recovered at the subthreshold level is expressed in the cell spiking output. These results demonstrate that the Simple or Complex nature of V1 RFs does not rely on the diversity of Complex-like components received by the cell from its synaptic afferents but on the imbalance between the weights of the Simple-like and Complex-like synaptic contributions.

## Introduction

Hubel and Wiesel initially distinguished two types of cells in cat and monkey primary visual cortex (V1) with respect to the properties of their spiking receptive fields (RFs): Simple cells, which are selective to the position and contrast polarity of oriented stimuli; and Complex cells, which are selective to the orientation of contours but not to their polarity and precise location (Hubel and Wiesel, 1962, 1968). In terms of RF architecture, Simple RFs have been mostly described as a single spatiotemporal filter followed by a half-rectifying nonlinearity (LN model) (Mancini et al., 1990; Heeger, 1992; DeAngelis et al., 1993). In contrast, Complex RFs have been classically modeled as two parallel LN branches, with linear spatiotemporal filters 90° out of phase and fully rectifying output nonlinearities (energy model) (Adelson and Bergen, 1985; Watson and Ahumada, 1985).

Many studies have pointed out that a significant proportion of V1 cells have RF properties intermediate between these Simple and Complex cell archetypes (Movshon et al., 1978; Emerson et al., 1987; Szulborski and Palmer, 1990; Baker, 2001). Recently, the spike-triggered covariance technique (STC), which relies on a principal component analysis of the spike-triggered ensemble, has been used in Macaque V1 to show that Simple and Complex discharge fields cover up multiple Complex-like subunits inconsistent with a strict dichotomy between Simple and Complex RF architectures (Rust et al., 2005; Chen et al., 2007). However, a complete description of both the 2D spatial structure and the temporal dynamics of the Complex-like subunits is still missing. It also remains largely unknown how the properties of these RF subunits correlate with the degree of *Simpleness* of the cell RF. Moreover, the issue of the functional diversity of RF subunits in Simple and Complex cells of cat V1 is still unclear (Touryan et al., 2002, 2005; Chen et al., 2007).

In the present study, we investigated the selectivity of the Complex-like components of cat V1 RFs with 2D ternary dense noise (DN) stimulation. As spike-triggered methods remain limited by the spike threshold operation to fully explore the functional space covered by the synaptic afferent (Rust et al., 2005), we measured the spatiotemporal selectivity of Simple-like and Complex-like RF subunits directly at the subthreshold level. We show that the Complex-like components of cat V1 cells rely on a large diversity of stimulus feature selectivity, regardless of the *Simpleness* of the cell RF. In both Simple and Complex cells, we found excitatory and inhibitory Complex-like subunits covering a large spectrum of orientation and spatial frequency tuning. Furthermore, we show that this diversity of feature selectivity measured at the subthreshold level is expressed in the cell spiking output. These results demonstrate quantitatively that the Simple or Complex nature of RFs in cat V1 does not rely on the diversity of Complex-like components received by the cell from its afferents but on the imbalance between the weights of the Simple-like and Complex-like synaptic contributions.

## Materials and Methods

##### Electrophysiology.

Data presented here were obtained from area 17 of anesthetized (alfatesin, 3 mg/kg/h) and paralyzed (pancuronium bromide, 0.2 mg/kg/h) cats of either sex, according to the American Physiological Society's Guiding Principles in the Care and Use of Animals. Animals used in these experiments have been bred in the Central Centre National de la Recherche Scientifique Animal Care at Gif-sur-Yvette (French Agriculture Ministry Authorization B91–272-105) under the required veterinary and National Ethical Committee supervision. Intracellular electrodes were made from thick wall borosilicate pipettes filled with 2 m potassium methyl sulfate and 4 mm KCl. Electrode resistances ranged between 60 and 90 mΩ. Recordings were performed using an Axoclamp 2A amplifier. Two cells of our population were labeled with biocytin, and their axonal and dendritic arborizations were reconstructed.

##### Visual stimuli.

Visual stimuli were generated using in-house software (Elphy, Gerard Sadoc, UNIC-Centre National de la Recherche Scientifique) and presented on a γ-corrected monitor with a refresh rate of 150 Hz and a background luminance of 12 cd · m^{−2}. The ON and OFF regions of the RF were first characterized with sparse noise stimuli consisting of nonoverlapping white (23 cd · m^{−2}) or black (1 cd · m^{−2}) squares presented one at a time. Cells were then stimulated by a 2D ternary DN, which consisted of random sequences of 10 × 10 (15 × 15 for one cell) nonoverlapping squares, which could be either white (23 cd · m^{−2}) or black (1 cd · m^{−2}), or gray (12 cd · m^{−2}) with equal probability. The stimulus refresh rate ranged from 30 to 75 Hz. To maximize the size of the stimulus dataset, the seed used for initializing the random process was changed for each sequence.

##### RF characterization.

The Simple or Complex nature of the cell RF was estimated from the spatial antagonism between ON and OFF subregions mapped with sparse noise stimuli, using the Pearson's correlation coefficient between ON and OFF responses at time of maximal response (ρ(ON,OFF)) (Priebe et al., 2004; Martinez et al., 2005; Mata and Ringach, 2005). This measure was further scaled to a 0–1 range so that the index tended to 0 for Complex RFs and to 1 for Simple RFs:
This index was computed on ON and OFF spiking responses for the majority of the cells (*n* = 27 of 38) or at the subthreshold level (*n* = 11 of 38) when the number of spikes was insufficient to compute a reliable RF estimate. Cells were considered as Simple or Complex when their RF Simpleness was >0.5 or <0.5, respectively. The results presented in this study were globally unchanged if we considered for all cells the correlation index computed from subthreshold ON and OFF responses: only one cell that was classified as Simple from its spiking response appeared as Complex at the *V*_{m} level. The RF Simpleness, defined as in Equation 1, was highly correlated with other classical indices, such as the discreteness index (Dean and Tolhurst, 1983) (*r* = 0.90; *p* < 10^{−11}; slope = 0.93; intercept = 0.02) or the segregation index (Van Hooser et al., 2013) (*r* = 0.89; *p* < 10^{−11}; slope = 0.85; intercept = 0.22), and our results were largely insensitive to the exact metric used to classify RFs.

To further characterize in detail the functional properties of the Simple-like and Complex-like components of the RF, we estimated the first- and second-order Volterra kernels from the subthreshold response to 2D DN, assuming a current-based model of the cell membrane. The estimated kernels were then represented as a bank of parallel LN branches, based on the eigenvectors of the second-order kernel (Marmarelis, 1989; Westwick and Kearney, 2003). In the following sections of Materials and Methods, we first present the theoretical background for obtaining the bank of filters from the estimated Volterra kernels and then describe how these mathematical tools were used in practice to recover the filter bank of the subthreshold RF.

##### From Volterra kernel expansion to filter bank representation: background and theory.

The selectivity of a visual neuron to the first- and second-order correlations of a stimulus can be estimated by fitting the cell response (*R*(*t*)) with the sum of the contributions of the first- and second-order Volterra kernels, which correspond to weighting functions in the first and second-order stimulus space, respectively:
which can be reformulated as a matrix product:
where **S**(*t*) is the stimulus input vector at time *t*, and *h*_{0}, **h _{1st}**, and

**h**correspond to the zero-, first-, and second-order Volterra kernels, respectively.

_{2nd}In the present study, we estimated the first- and second-order Volterra kernels from V1 cell subthreshold responses evoked by 2D DN stimuli. Before detailing practically how Equation 2 was solved (see Filter bank decomposition of V1 subthreshold responses), we first explain here how the estimated Volterra kernels were represented as a bank of parallel filters. Because of the high dimensionality of the DN stimulus, the estimated *h _{2nd}* kernel is not the most optimal way to visualize and quantify the selectivity of the cell to the second-order correlations of the stimulus. One convenient alternative is to reexpress the estimated second-order kernel as an orthonormal basis of eigenvectors (Westwick and Kearney, 2003; Marmarelis, 2004). By definition, this set of eigenvectors respects the following equation:
where

**B**is the matrix containing the eigenvectors in its columns and

**Δ**is a diagonal matrix of eigenvalues, which indicate the contributions of each eigenvector to

*h*.

_{2nd}Therefore, by replacing **h _{2nd}** in Equation 2, we obtain:
which can be developed as follows:
where

*h*are the eigenvectors of

_{2.k}*h*and λ

_{2nd}*are their respective eigenvalues.*

_{k}The eigenvector decomposition of *h _{2nd}* thus operates a change of coordinates, which allows reformulating the second-order Volterra kernel of the RF as a parallel set of first-order-dimensioned filters

*h*(eigenvectors) whose outputs are squared and weighted with specific coefficients (eigenvalues).

_{2.k}Although the eigenvalues computed in the eigenvector expansion define the relative weight of each eigenvector contribution to the cell response, they are not indicative of their significance. In the present study, we used the same approach as described by Rust et al. (2005) and used bootstrapping to identify the eigenvectors whose contribution to the *h _{2nd}* variance was above chance level (

*p*< 0.05). More precisely, the estimated

*h*and

_{0}*h*contributions (the first terms of Eq. 5) were subtracted from the cell response. The residual response was then randomly time shifted relative to the stimulus sequence,

_{1st}*h*was estimated, and its eigenvalues were determined. By repeating this procedure 300 times, we obtained the distribution of eigenvalues expected when the stimulus and the response are not correlated. If the largest (or smallest) eigenvalue of

_{2nd}*h*was not in the confidence interval of this distribution, we concluded that the corresponding eigenvector contributed significantly to the second-order kernel output. This eigenvector contribution was then subtracted from the response, and we assessed whether the next largest (or smallest) eigenvalue contributed significantly to the rest of the response.

_{2nd}The selected significant eigenvectors *h _{2.k}* correspond to a minimal set of orthogonal vectors, which describe entirely the functional space covered by the significant components of the second-order kernel. Still, it is not unique, and any linear combination of these eigenvectors would be equivalent provided it covers the same functional space. It is noteworthy that, although the eigenvectors are orthogonal, this constraint can be achieved in many different ways along the spatial and temporal dimensions; in particular, there is no constraint imposing orthogonality of the eigenvectors in the orientation or the spatial phase domain: these special cases of spatial quadrature simply reflect that the second-order kernel contains distinct components selective to different ranges of orientation or spatial phase, and which can be reexpressed conveniently by a linear sum of two orthogonal functional axis spanning the appropriate domain of selectivity. Similarly, because of the orthogonality constraint, the excitatory or inhibitory nature of the recovered eigenvectors cannot be directly interpreted as purely excitatory or inhibitory inputs. However, the sign of a given eigenvalue is indicative of the net excitatory/inhibitory balance along the stimulus feature for which the eigenvector is selective.

This method, originally described by Marmarelis (1989) (see also (Westwick and Kearney, 2003; Marmarelis, 2004) as the Principal Dynamic Mode method is virtually equivalent to the principal component analysis of the covariance matrix performed in spike-triggered covariance methods (Touryan et al., 2002; Rust et al., 2005; Schwartz et al., 2006), except that the second-order Volterra kernel is estimated by a least-squares method rather than by reverse correlation. Moreover, for the purpose of determining subthreshold RF subunits, the first- and second-order kernels were interpolated at ∼1 ms resolution, to allow the estimation of the time constant parameter in the current-based model used to estimate the subthreshold filter bank (see Filter bank decomposition of V1 subthreshold responses). In practice, this interpolation is similar to classical peristimulus triggered averaging methods and relatively straightforward: by considering that every single frame of the stimulus sequence consists of one single 1 ms impulse in an interstimulus interval and is zero elsewhere, solving Equation 2 is indeed equivalent to solve independently the following set of equations:
where T correspond to the number of 1 ms time bins in the interstimulus interval, **S ^{1}**(

*t*) corresponds to the stimulus impulses, and

*R*(

^{i}*t*) corresponds to the response down-sampled at the stimulus frequency by taking the value at the

*i*1 ms time bin. The estimated

^{th}**h**and

_{1st}^{i}**h**thus give the value of the first- and second-order kernels at successive time bins of the interstimulus interval. The eigenvectors of the second-order kernel were obtained by concatenation of the eigenvectors of the successive

_{2nd}^{i}**h**kernels, multiplied by their respective eigenvalue. The resulting orthogonal vectors were then normalized to form an orthonormal eigenvector basis of the estimated

_{2nd}^{i}*h*kernel.

_{2nd}##### Filter bank decomposition of V1 subthreshold responses.

As the present study addresses the issue of the selectivity of the second-order RF components at the subthreshold level, we searched for an optimal way to apply the filter bank decomposition method described above to intracellular responses. One issue is that the filter bank recovered directly from the membrane potential response (by equating *R*(*t*) with *V*_{m}(*t*) in Eq. 2) is unlikely to match the functional properties of the synaptic inputs because the latter are filtered by the membrane time constant (Pillow and Simoncelli, 2003). To take into account the membrane filtering of the filter bank components, we reformulated the Volterra kernel expansion considering the following current-based model of the cell response:
in which *V*_{syn} corresponds to the synaptic component of the cell response (*V*_{syn} = *I*_{syn} × *R*_{mb}), *V*_{rest} corresponds to the resting potential, and τ_{mb} corresponds to the membrane time constant. After replacing *V*_{syn} by the second-order Volterra expansion formulated in Equation 2, Equation 7 becomes:
with *h _{0}* corresponding to

*V*

_{rest}and the outputs of

*h*

_{1st}and

*h*

_{2nd}corresponding to the evoked synaptic component

*V*

_{syn}.

Practically, the derivative of the membrane potential (*dV _{m}/dt*) was estimated from the Fourier transform of the membrane potential multiplied by

*e*

^{−2iπf}.

*V*and

_{m}*dV*/

_{m}*dt*were low-pass filtered at 75 Hz and resampled with 1 ms resolution. When necessary, we restricted spatially and temporally the stimulus positions over which we estimated the second-order kernel parameters, to keep a ratio of at least 1:5 between the number of kernel parameters and the number of data points. Equation 8 was then solved for different presumed values of time constant τ

_{mb}, using the Cholesky factorization implemented by the fast orthogonal algorithm of Korenberg (1988). We thus obtained a collection of Volterra kernel estimates (

*h*

_{1st}

^{τmb}and

*h*

_{2nd}

^{τmb}), each corresponding to one particular assumption of the membrane time constant (τ

_{mb}). The eigenvectors of the second-order kernel estimates were calculated, the significant eigenvalues were identified (

*p*< 0.05, see above), and we reconstructed the Volterra kernel output

*V*

_{syn}

^{τmb}according to Equation 5, considering only the eigenvectors that were significant. The reconstructed membrane voltage response

*V*

_{m}

^{τmb}corresponding to each assumption of time constant was then calculated by convolution with the membrane time constant: We then measured the goodness-of-fit as the percentage of variance explained by each reconstructed membrane voltage

*V̂*

_{m}

^{τmb}and identified the time constant value that gave the best fit (τ*

_{mb}). The

*h*and

_{1st}*h*kernels estimated with this assumption of time constant were selected as the best current-based estimate of the RF. Eventually, the RF filter bank was defined as the set of subunits corresponding to the first-order kernel

_{2nd}*h*and the collection of significant eigenvectors

_{1st}*h*of the second-order kernel.

_{2.k}##### Validation of filter bank decomposition on simulated RFs.

To validate our approach, we simulated current-based second-order RF models consisting of multiple LN branches (Fig. 1*A*), with one filter contributing linearly to the response (E_{0}) and several filters with positive or negative square output nonlinearities (E_{1}, E_{2}, I_{1},I_{2}). The filter bank output v_{syn} was transformed into membrane voltage fluctuations according to a current-based process (Eq. 7). The simulated response to DN stimulation was fitted with the subthreshold filter bank estimation method described above. The estimated filter bank was compared with the filter bank expected from the RF model parameters: the estimated Simple-like subunit *h _{1st}* was compared with the filter of the model whose output contributed linearly to the response (

*h*

_{1st}

^{model}=

*E*

_{0}); the estimated Complex-like subunits

*h*were compared with the eigenvectors of the expected second-order kernel (

_{2.k}*h*

_{2nd}

^{model}), computed from the filters of the model whose linear outputs were fully rectifying (Westwick and Kearney, 2003): The estimated filter bank perfectly recovered the expected RF subunits and their relative weights (eigenvalues) (current-based estimate, Fig. 1

*B*) and the estimated optimal time constant value (Fig. 1

*C*) always recovered exactly the time constant used in the simulation (Fig. 1

*D*). The success of this current-based RF estimation method relies on the reduction of dimensionality performed before assessing the goodness-of-fit, when selecting the significant eigenvectors to approximate the output of the second-order kernel: if this step is omitted and the full second-order kernel is used to reconstruct the response, no optimal solution appears over the set of assumed time constant values (Fig. 1

*C*, gray line). Filter banks were also estimated from these same simulated RFs with a direct

*V*-based method, by neglecting the membrane filtering in the Volterra kernel expansion (i.e., by equating

_{m}*R*(

*t*) with

*V*(

_{m}*t*) in Eq. 2;

*V*-based estimate, Figure 1

_{m}*B*). In this case, the recovered filter bank did not match the RF subunits expected from the model: their temporal profiles were altered, their relative weights were misestimated, and additional spurious components appeared

*h*

_{exc}

^{2.3}and

*h*

_{exc}

^{2.4}; Fig. 1

*B*.

##### Filter bank decomposition of V1 spiking responses.

Spike times were recorded by detecting the crossing of the membrane voltage with a −30 mV threshold. We then selected cells for which we had at least 25 times more spikes than stimulus dimensions (*n* = 12 of 38) for estimation of their spiking filter bank. The recorded spike train was down-sampled to the stimulus frequency, and the zero-, first-, and second-order Volterra kernels (*h*_{0}, *h*_{1st}, and *h*_{2nd}, respectively) were estimated by fitting Equation 2 to the spiking response. The significant eigenvectors of the spiking second-order kernel (*h*_{2nd}) were then identified by bootstrapping (*p* < 0.05; see above), and the filter bank was defined as the collection of filters consisting of the normalized first-order kernel *h _{1st}*, and the collection of significant eigenvectors

*h*of the second-order kernel. This method is equivalent to the classical spike-triggered covariance technique used in earlier studies (Touryan et al., 2002, 2005) because the covariance matrix estimated in STC methods is closely related to the second-order kernel estimated by cross-correlating the stimulus with the spike train (Marmarelis, 2004). For the sake of comparison, we computed the spike-triggered ensemble and estimated the spike-triggered average (STA) and the significant eigenvectors of the spike-triggered covariance matrix, as described previously (Rust et al., 2005; Schwartz et al., 2006) (data not shown). Compared with the Volterra kernel approach used in the present study, the STA/STC technique led to very similar results: across these 12 cells, the RF subunits recovered with both approaches were highly correlated (average correlation value: 0.91), and we found almost the same number of second-order components with the STA/STC approach (

_{2.k}*n*= 28) than with the Volterra approach (

*n*= 31).

##### Measures of RF subunit weights.

To avoid misestimating the RF subunits weights due to overfitting of the response by the estimated filters, the weighting coefficients of the *h _{1st}* and

*h*subunits were fitted to a set of data, which had not been used for kernel estimation. Because we had only single-trial data, this fitting dataset was constructed as follows: we estimated the filter bank on 95% of the total recording length, reconstructed the output of the estimated filters for the remaining 5%, and repeated this procedure until we completed the reconstruction of the full single-trial response. The weighting coefficients of the filter bank subunits were then estimated by fitting the following equation to the recorded response

_{2.k}*V*(

_{m}*t*): where

*V̄*

_{m}

^{1st}(

*t*) and

*V̄*

_{m}

^{2.k}are the reconstructed outputs of the filter bank subunits, filtered by the estimated membrane time constant τ

_{mb}*: where h̅1̅s̅t̅ denotes the first-order kernel normalized to unity Euclidian norm.

The relative weight of each eigenvector *h _{2.k}* (Complex-like subunit relative weight) was defined by normalization with the weight of the first-order kernel

*h*

_{1st}: where |α

_{2.}

*| denotes the absolute value of the coefficient.*

_{k}The balance between the inhibitory and excitatory components of the filter bank was defined as follows:

##### Predictive correlation.

The ability of each filter bank to predict the subthreshold response to DN (*Predictive correlation*) was quantified by the Pearson's correlation coefficient between the reconstructed filter bank output *V _{m}*), considering a set of single-trial data which had not been used for kernel estimation (same dataset as the one used to estimate the subunit weights, see above).

The ability of the filters *h*_{2}.* _{k}* to predict the second-order components of the response (

*Second-order predictive correlation*) was measured as the Pearson's correlation coefficient between the

*V*fluctuations reconstructed from the

_{m}*h*

_{2}.

*filters (*

_{k}*V*

_{m}

^{2nd}(

*t*) =

*V*(

_{m}*t*) − (

*h*

_{0}+ α

_{1st}×

*V̄*

_{m}

^{1st}(

*t*))).

##### STA of subthreshold filter bank components.

The correlation (STA) between the recorded spike train (Spike(t)) and the estimated subthreshold filter bank components was estimated by multidimensional linear regression, fitting the following equation for different time lags τ:
where *V̄*_{m}^{1st} and *V̄*_{m}^{2.k} correspond to the reconstructed outputs of the subthreshold RF subunits, filtered by the estimated membrane time constant (i.e., the unweighted *V _{m}* contribution of each RF subunit; see Eqs. 12 and 13, respectively) and

*b*is a constant parameter. To avoid overfitting, the STAs were estimated from a dataset of single-trial responses, which had not been used for the kernel estimation (same fitting dataset as described above). Only cells that showed at least 50 times more spikes than the number of subthreshold subunits were included in this analysis (

*n*= 24).

The estimated STA_{1st} and STA_{2.}* _{k}* correspond to the cross-correlation between the recorded spike train and the stimulus projected onto the set of RF subunits identified at the subthreshold level, filtered by the estimated membrane time constant. They are equivalent to spike-triggered averages and indicate the selectivity of the spike response to the components of the stimulus matching the feature selectivity of the subthreshold RF subunits. The amplitudes of the STAs at zero time lag (i.e., STA

*(0)) were measured as indicators of the weight of each subunit in the cell spiking response. The relative weights of the Complex-like subunits (STA*

_{x}_{2.k}) were then defined by normalization with the weight of the first-order subunit (STA

_{1st}): The statistical significance of the contribution of each RF subunit to the spiking response was measured by comparing the cross-correlation value estimated at zero time lag STA

*(0) to the distribution of cross-correlation values recovered at long correlation delays (*

_{x}*p*< 0.01).

##### Tuning properties of RF subunits.

The tuning properties of the estimated RF subunits were quantified from their 2D spatial power spectrum at time of maximal response: the preferred orientation and spatial frequency were deduced from the position of the maximum in the spatial power spectrum; the orientation tuning curve was measured at best spatial frequency and fitted by the sum of two von Mises distributions, forced to be at 90° one of the other; the half-width at half-height (HWHH) was measured from the fitted curve. Filters were considered to have similar orientation and spatial frequency preferences when they differed by <30° and 0.5 octave, respectively. Their preferred spatial phases were then measured from the phase of the 2D Fourier transform at optimal orientation and spatial frequency. The selectivity of the excitatory and inhibitory RF subunits was compared using the HWHH and the preferred spatial frequency measured in the excitatory and inhibitory pooled power spectra. The excitatory and inhibitory pooled power spectra were calculated for each filter bank as the sum of the power spectra of the excitatory and inhibitory subunits, respectively, weighted by the coefficients of the subunits (α_{2.k}).

##### Spatiotemporal envelope of RF subunits.

The significance of the filter bank subunits in each spatiotemporal position (*z*-score) was estimated by bootstrapping, using the SD of the filter values measured over the collection of *h _{1st}* and

*h*estimated when the response was time shifted relative to the stimulus sequence. Areas defined by contours delimiting 99% significant regions (

_{2.k}*z*-score ≥ 2.33) were plotted as a function of time for each filter of the bank. The maximal spatial extent of these contours (measured by their apparent visual diameter expressed in degrees of solid angle) as well as the timing of this maximum (peak latency) were measured for each filter of the bank. The onset latency was derived from the time derivative of the filter spatial extent: we defined a size threshold corresponding to half of the maximal spatial extent and went backwards in time until the derivative of the filter size fell to <10% of the derivative calculated at half-amplitude, for five continuous time steps (5 ms).

## Results

A total of 38 cells were recorded intracellularly in the primary visual cortex (V1) of the anesthetized cat. The Simple or Complex nature of their RF (*Simpleness*) was estimated from the spatial antagonism between ON and OFF subregions mapped with sparse noise stimuli (Priebe et al., 2004) (see Materials and Methods). We then investigated the selectivity of the Simple-like and Complex-like components of both Simple and Complex cell RFs from synaptic responses evoked by 2D ternary DN. We first focused our analysis on the RF components estimated at the subthreshold level and further assessed their contribution to the cell spiking response.

### Recovering the Simple-like and Complex-like RF subunits at the subthreshold level

The decomposition of the sensory response into kernel series provides a systematic method to identify the linear and nonlinear components of RFs. In the present study, we estimated the first- and second-order Volterra kernels of V1 subthreshold RFs from DN responses (Fig. 2*A*) with a least-squares method. In this Volterra expansion, the first-order kernel *h*_{1st} quantifies the linear components of the response, which are evoked along individual stimulus positions (Fig. 2*B*). It reveals the stimulus feature for which the cell receives net excitatory and inhibitory contributions in a *push-pull* arrangement and can be considered as the Simple-like RF component in the strict sense. In contrast, the second-order kernel *h*_{2nd} measures the nonlinear components of the response selective to pairwise correlations between stimulus positions (Fig. 2*B*). It can be considered as the Complex-like component of the RF: the diagonal elements of the second-order kernel (*h*_{2Diag}, Fig. 2*B*) correspond to contributions that are insensitive to the polarity of the contrast along individual positions of the stimulus, whereas the off-diagonal terms provide information on the stimulus features to which these polarity-invariant components are selective.

Because of the high dimensionality of the second-order space associated with 2D DN stimuli, the spatiotemporal selectivity of the Complex-like RF components is not easily readable from the structure of the second-order kernel. We thus searched for a more efficient and parsimonious representation of the second-order kernel and opted for a decomposition method similar to the principal component analysis of the spike-triggered covariance matrix described in other studies (Touryan et al., 2002; Rust et al., 2005). Practically, after the Volterra kernels were estimated, we computed the eigenvectors of the second-order kernel and identified those with eigenvalues above chance level (Fig. 2*C*, left; see Materials and Methods). This set of significant eigenvectors corresponds to an orthonormal basis of linear filters (*h _{2.k}*), dimensioned as the first-order kernel, which, by definition, account for the significant component of the second-order kernel provided a quadratic transformation of their linear output (Fig. 2

*C*, right; see Materials and Methods). The second-order kernel can thus be reformulated as a set of eigenvectors whose linear outputs are fully rectified and weighted with specific coefficients (eigenvalues). The excitatory or inhibitory nature of the eigenvector contribution is defined by the sign of the weighting coefficient (which determine the sign of the square output function).

Eventually, the RF is thus represented as a parallel bank of multiple LN input branches composed of one Simple-like subunit *h _{1st}* corresponding to the first-order kernel and a set of Complex-like subunits

*h*

_{2.}

*corresponding to the significant eigenvectors of the second-order kernel (Fig. 2*

_{k}*D*). This filter bank decomposition of subthreshold RFs is virtually equivalent to the spike-triggered covariance technique (Touryan et al., 2002; Rust et al., 2005), except that the second-order nonlinearities of the RF are initially estimated by a least-squares method rather than by reverse correlation. Moreover, to account for the membrane filtering of the synaptic inputs, we used a more realistic current-based model of the cell membrane under the assumption that the filter bank output is converted into membrane voltage fluctuations through a passive leak process (see Materials and Methods; Eq. 7). The membrane time constant parameter was estimated in the least-squares sense, concomitantly with the kernels coefficients (see Materials and Methods). For all cells, we found an optimal time constant value that best explained the response (Fig. 2

*E*). Over our cell population, the estimated time constant values ranged from 5 to 27 ms (average = 14 ms; Fig. 2

*F*). Although these values fit well with the range of membrane time constants measured at rest,

*in vivo*in V1 cells (Monier et al., 2008), they should not be taken as an accurate estimate of the membrane time constant of the cell under DN stimulation: they only correspond to an optimal approximation of the membrane time constant, assuming that the evoked membrane voltage fluctuations result from a current-based RF. Still, simulations of RF models showed that this time constant parameter reduces the spurious effect of the membrane filtering on the selectivity and the relative weights of the estimated RF subunits (Fig. 1; see Materials and Methods). Over our cell population, this current-based decomposition method always provided a better fit of the evoked response (in the least-squares sense) than a more straightforward

*V*

_{m}-based method neglecting the filtering of the synaptic inputs by the cell membrane (data not shown).

### Subthreshold subunits of Simple and Complex V1 RFs

Figure 3 shows the spatiotemporal profiles of the Simple-like and Complex-like RF subunits recovered from the subthreshold responses of two Complex cells (Fig. 3*A*; RF *Simpleness*: cell 1, 0.35; cell 2, 0.17). The Complex-like RF subunits are mostly shaped as typical Simple RFs, with several subfields of alternating polarities, indicative of their orientation selectivity. As reported previously in Macaque V1, we first notice that the Complex-like subunits of these Complex RFs cover a much larger functional space than expected from the energy model: both cells in Figure 3 reveal more than two Complex-like subunits, which can be excitatory or inhibitory and whose spatial profiles cover a variety of orientation preference and spatial frequency selectivity (Fig. 3*B*,*C*).

Figure 4 illustrates the Simple-like and Complex-like subunits of two typical Simple cells (Fig. 4*A*; RF *Simpleness*: cell 3, 0.81; cell 4, 0.82), a layer 4 spiny stellate neuron (Fig. 4*B*, cell 3) and a layer 5/6 pyramidal neuron (Fig. 4*B*, cell 4). If these cells were responding to DN as predicted by the LN model, our analysis should reveal only one complex-like excitatory subunit, with the same selectivity as the linear Simple-like component. In contrast, the filter bank of these Simple cells revealed several Complex-like subunits (Fig. 4*C*,*D*), both excitatory and inhibitory and whose spatial profiles are markedly different from that expressed by the Simple-like subunit.

Across the whole population of V1 cells that we recorded, ∼90% (34 of 38) revealed at least one Complex-like subunit (Fig. 5*A*); ∼50% of the Complex cells (13 of 27) showed three or more excitatory Complex-like subunits, and ∼65% of the Simple cells (7 of 11) had at least two excitatory Complex-like subunits. Also, ∼55% of the Complex cells (15 of 27) and ∼65% of the Simple cells (7 of 11) exhibited at least one inhibitory Complex-like subunit. Overall, we found that the total number of Complex-like RF subunits recovered at the subthreshold level is not correlated with the *Simpleness* of the RF (*p* = 0.78).

Although the recovered Complex-like subunits cannot be taken literally as the structural correlate of presynaptic cells, they still correspond to a minimal representation of the functional space covered by the Complex-like synaptic inputs of the cell; as such, they are indicative of an at least equivalent functional diversity among the effective presynaptic afferents. Our data are thus consistent with previous spike-triggered covariance analysis performed in Macaque V1 (Rust et al., 2005) and further show at the intracellular level that the characterization of cat V1 RFs as either Simple or Complex does not reflect the diversity of Complex-like contributions received from the synaptic afferents.

We further investigated how the relative weights of these Complex-like RF subunits scale with the *Simpleness* of the RF. To avoid overfitting of the cell response, the weighting coefficients of the filter bank were estimated on a stimulus dataset, which had not been used for kernel estimation (see Materials and Methods). The weights of the Complex-like subunits were normalized by the weight of the Simple-like subunit (Complex-like subunit relative weight *h*_{2.}* _{k}*/

*h*

_{1st}). We found that the relative weights of the Complex-like subunits were significantly correlated with the RF

*Simpleness*: both excitatory and inhibitory Complex-like subunits contributed less to the cell response relative to the Simple-like component as the cell RF was more Simple (Fig. 5

*B*; excitatory:

*r*= −0.56;

*p*< 10

^{−8}; inhibitory:

*r*= −0.6;

*p*< 10

^{−4}). Our data thus demonstrate that the Simple or Complex nature of cat V1 RFs does not depend on the number of functionally distinct Complex-like influences the cell receives but rather on the imbalance between the synaptic weights of the Complex-like and Simple-like contributions (Chance et al., 1999; Tao et al., 2004). For cells that revealed both excitatory and inhibitory Complex-like components (

*n*= 22), we measured the sum of the weights of the inhibitory subunits over the sum of the weights of all Complex-like RF subunits (excitatory + inhibitory). We found that this “inhibitory–excitatory” balance is significantly correlated with the degree of

*Simpleness*of the RFs (Fig. 5

*C*;

*r*= 0.54,

*p*< 0.01), suggesting that the Complex-like contributions received by Simple cells are composed of proportionally more inhibitory synaptic inputs than in Complex cells. These results were globally unchanged if we considered other classical metrics of RF classification, such as the discreteness index (Dean and Tolhurst, 1983) or the segregation index (Van Hooser et al., 2013) (data not shown).

We quantified the ability of the RF subunits to explain the membrane voltage fluctuations evoked by DN stimuli, by measuring the correlation between predicted and recorded responses on single-trial datasets, which had not been used for kernels estimation (predictive correlation; see Materials and Methods). The prediction of DN responses by the subthreshold filter banks showed an average predictive correlation of 0.31 (ranging from 0.08 to 0.80; Fig. 6*A*), similar to average scores reported previously from extracellular studies, which used stimuli of comparable dimensionality and included both Simple and Complex RF types (David et al., 2004; Chen et al., 2007). These scores are to be interpreted as lower bound estimates because they were established from predictions of single-trial responses (response variability was not averaged out). We found that the ability of the filter bank to predict DN responses largely depended on the *Simpleness* of the RF (*r* = 0.7, *p* < 10^{−5}). This correlation, however, mostly relied on the contribution yielded by the Simple-like RF component: when we considered only the Complex-like components of the DN responses, by subtracting the Simple-like contributions (see Materials and Methods), we found that response prediction was improved in a comparable fashion in both Simple and Complex RFs (second-order predictive correlation; Fig. 6*B*; *r* = 0.03, *p* = 0.89).

### Tuning properties of subthreshold RF subunits

To quantify in more detail the spatial tuning properties of V1 RF subunits, we measured their preferred orientation, orientation tuning width, and preferred spatial frequency from their spatial power spectrum at time of maximal response (Fig. 7*A*). The preferred orientation of the Complex-like subunits was measured relative to that of the Simple-like component (Relative Orientation). Over the cell population, 58% of the excitatory Complex-like subunits were tuned to orientations close (±30°) to that expressed by the Simple-like subunit, whereas 22% were tuned to oblique (±30–60°) and 20% to orthogonal features (±60–90°) (Fig. 7*B*, left). This distribution was in contrast with that measured from the inhibitory Complex-like subunits, which revealed a slight bias for orientation orthogonal to the preferred orientation of the Simple-like subunit (Fig. 7*B*, right; ±30° for 27% of the filters; ±30–60° for 27% of the filters; ±60–90° for 46% of the filters).

The orientation tuning widths (HWHH) of the Simple-like and Complex-like RF subunits were not significantly different (Fig. 7*C*; 34° and 35° on average, respectively, *p* = 0.82), and excitatory and inhibitory Complex-like subunits also had similar orientation tuning widths (35° and 34° on average, respectively, *p* = 0.86), consistent with previous reports based on direct measurements of excitatory and inhibitory conductances (Monier et al., 2003).

The comparison of preferred spatial frequencies (Fig. 7*D*) showed that the excitatory Complex-like subunits are tuned to higher spatial frequencies compared with the Simple-like subunits (paired *t* test, *p* < 10^{−6}) or with the inhibitory Complex-like subunits (paired *t* test, *p* < 0.01). This latter difference was mostly because 23% (*n* = 7) of the inhibitory Complex-like subunits were not selective to orientation (i.e., with a preference for null spatial frequencies in the Fourier domain), whereas only 2% (*n* = 2 filters) of the excitatory subunits were in the same case. These inhibitory Complex-like subunits with no orientation selectivity were mainly found in Simple cells (5/7, Fig. 7*D*) and generally revealed one single subregion, which was neither elongated nor shaped for any oriented feature detection but still covered several positions of the cell RF (Fig. 4, *C*,*D*).

To relate our findings to classical LN models, we measured the preferred spatial phase of RF subunits tuned to similar orientation (±30°) and spatial frequency (±0.5 octave) and identified those with closely matching spatial phase (±30°) or in spatial phase quadrature (i.e., 90° out of phase, ±30°). In the case of Simple cells, we found that ∼35% (4 of 11) had one dominant Complex-like subunit with a feature selectivity matching that of the Simple-like component, as expected from an LN model. Nevertheless, all of these Simple cells also exhibited additional Complex-like subunits (among which some were suppressive). In the case of Complex cells, ∼50% (14 of 27) showed two dominant excitatory Complex-like subunits in spatial phase quadrature, as expected from the energy model of Complex cells. Again, a majority of these Complex cells (11 of 14) also showed excitatory and/or inhibitory Complex-like subunits with additional feature selectivity. These complex RFs with multiple Complex-like components and a pair of dominant excitatory subunits consistent with the energy model are reminiscent of previous results obtained in awake monkeys with spike-triggered covariance analysis (Chen et al., 2007).

Over our cell population, the orientation or spatial frequency selectivity of the Complex-like RF subunits was not correlated with the degree of *Simpleness* of the RFs. Our analysis thus demonstrates quantitatively that the Complex-like synaptic inputs of cat V1 RFs cover a broad spectrum of orientation and spatial frequency selectivity, regardless of the Simple or Complex nature of the RF. Although we could not correlate across our cell population the diversity of Complex-like RF subunits to cell types or laminar positions, the two reconstructed cells of our dataset (cells 3 and 4) show that, in layer 4 or 6, typical Simple cells with spatially nonoverlapping ON and OFF subfields can also exhibit a large repertoire of feature selectivity in their Complex-like synaptic afferents.

### Spatiotemporal envelope of subthreshold RF subunits

We compared the spatiotemporal profiles of the Simple-like and Complex-like RF subunits by measuring their spatial extents and their relative timing. The maximal spatial extent covered by each RF subunit was measured as the largest area delineated by the 99% significant response contour. We found that, in most cells, the Complex-like RF subunits had significantly smaller size than the Simple-like component (Fig. 8*A*; paired *t* test, *p* < 10^{−12}): whereas the Simple-like components had their largest spatial spread ranging from 0.4 to 4.3 degrees of visual angle (average = 1.70 degrees), the size of the Complex-like subunits ranged from 0.4 to 2.9 degrees (average = 1.4 degrees). The size of the Simple-like and Complex-like subunits tended to increase and decrease, respectively, with the *Simpleness* of the RF (*h _{1st}*:

*r*= 0.40,

*p*= 0.02;

*h*:

_{2.k}*r*= −0.21,

*p*= 0.02), resulting in a significant shrinkage of the Complex-like subunits relative to the Simple-like subunits in more Simple RFs (

*r*= −0.48,

*p*< 10

^{−5}). No significant difference was found between excitatory and inhibitory Complex-like subunits regarding their spatial extent.

In the temporal domain, the Simple-like and Complex-like RF components revealed a systematic temporal interplay (Onset latency, Fig. 8*B*): the inhibitory Complex-like subunits were systematically delayed compared with the excitatory Complex-like subunits (average delay = 23 ms; paired *t* test, *p* < 10^{−4}), which themselves generally appeared with a slightly delayed onset relative to the Simple-like component (average delay = 3 ms; paired *t* test, *p* < 0.01). The onset latency of the Simple-like or Complex-like RF subunits was not significantly correlated with the RF *Simpleness*, but the delay between the Simple-like and the excitatory Complex-like subunits increased significantly with more Simple RFs (*r* = 0.35, *p* < 10^{−3}). A tight temporal arrangement was also observed when comparing the times to maximal spatial extent (Peak latency, Fig. 8*C*): the excitatory Complex-like subunits reached their maximal spatial spread before the Simple-like components (paired *t* test, *p* < 10^{−5}), whereas the inhibitory Complex-like subunits expressed their full spatial profile significantly later (paired *t* test, *p* < 10^{−5}). This temporal interplay between excitatory and inhibitory RF subunits suggests the existence of a Complex-like inhibition whose effective contribution would be systematically delayed and/or have longer time constants than the Complex-like excitatory components (Levy et al., 2013).

### Complex-like components of V1 cell spiking responses

For a subset of cells, we estimated the Simple-like and Complex-like subunits of spiking RFs, using the first- and second-order Volterra kernels measured from spiking responses (*n* = 12; 8 Complex cells and 4 Simple cells; conservative threshold: >25 times more spikes than spatiotemporal positions in the stimulus space; see Materials and Methods). Figure 9*A*, *B* shows the filter banks obtained from the spiking responses of cells 2 and 3, respectively. The filter bank estimated directly from their spiking response revealed a much lesser number of Complex-like subunits than observed at the subthreshold level (Figs. 3*C* and 4*C*, respectively).

Over the analyzed subset of cells, we recovered half as many Complex-like RF subunits at the spiking level than at the subthreshold level (31 vs 61; Fig. 9*C*). Moreover, the Complex-like subunits recovered from spiking responses tended to better match the LN and energy models of Simple and Complex cells than what we observed at the subthreshold level (Touryan et al., 2002). Nevertheless, it is likely that the discrepancy between subthreshold and spiking filter bank estimates is attributable, at least in part, to an insufficient number of recorded spikes, relative to the dimensionality of the stimulus (Rust et al., 2005). Therefore, on the sole basis of this comparison, we could not conclude whether or not the diversity of feature selectivity observed at the subthreshold level is represented in the cell spiking output.

To answer this question despite the relatively low number of collected spikes, we reduced the dimensionality of the stimulus space and computed STAs of the stimulus projected onto the RF subunits identified at the subthreshold level. Practically, for cells that showed a sufficient number of spikes (*n* = 24, see Materials and Methods), we estimated the cross-correlation (STA) between the recorded spike train and the reconstructed response of each subthreshold subunit, for different time lags (Fig. 10*A*). Figure 10*A*, *B* shows the STAs obtained for the same 4 cells as presented in Figures 3 and 4. In all cases, the excitatory Complex-like components showed a significant (*p* < 0.01) positive deflection at the time of spike emission, whereas the amplitude of the inhibitory Complex-like contributions decreased. The spiking response of these cells was thus significantly correlated with the inputs received along the stimulus feature expressed in the subthreshold RF subunits. The weight of each RF subunit contribution to the cell spiking response was measured from the amplitude of the STAs at time of spike emission (i.e., at *t* = 0 on the STA abscissa). The weights of the Complex-like subunits were expressed relative to that of the Simple-like component (Complex-like subunits relative weight, STA; see Materials and Methods). We found that the relative weights of the Complex-like subunits measured from the STAs were highly correlated with those we had previously estimated at the subthreshold level (Complex-like subunits relative weight, *V _{m}*) (excitatory subunits:

*r*= 0.82,

*p*< 10

^{−12}; inhibitory subunits:

*r*= 0.7,

*p*< 10

^{−3}; Fig. 10

*C*). This result demonstrates that the Complex-like RF subunits identifiable at the subthreshold level are represented in the spiking response in the same proportion as their relative subthreshold contribution. We conclude that the diversity of feature selectivity observed at the subthreshold level in cat V1 RFs is expressed in the cell spiking output and is thus propagated to downstream network targets.

## Discussion

The present study is based on a new analytical method never performed before at the intracellular level and provides new insights into the functional organization of cat V1 RFs. Our conclusions are based on a second-order analysis of V1 synaptic responses evoked by 2D ternary DN. The estimated subthreshold RFs were decomposed into a parallel set of multiple LN branches consisting of a Simple-like subunit, whose linear output accounts for the *push-pull* components of the RF, and a variable number of Complex-like subunits, which contribute in a full-rectified manner to the cell response.

Our data show that, at the synaptic level, the Complex-like RF subunits of cat V1 cells rely on a large diversity of feature selectivity, regardless of the Simple or Complex nature of the RF. In both Simple and Complex cells, we found multiple Complex-like RF subunits whose tuning properties can differ largely from that expressed by the Simple-like RF component. We show that, although excitatory Complex-like subunits are often tuned for orientations close to that expressed in the Simple-like subunit (±30°), a substantial proportion is selective for oblique and crossed orientations (42%). Complex-like inhibitory subunits also exhibit orientation preference, which differs notably from that of the Simple-like subunit. We also found a substantial proportion of Complex-like inhibitory subunit, which exhibited no orientation selectivity, especially in Simple cells. These inhibitory components constitute additional evidences favoring the existence of an orientation unselective Complex inhibition, which could participate, at least in some V1 cells (∼45%, 5 of 11; in our population), to the contrast-invariance of V1 cell selectivity (Hirsch et al., 2003; Lauritzen and Miller, 2003; Nowak et al., 2008; Liu et al., 2009; Shapley and Xing, 2013). Interestingly, our data further show that the balance between the relative weights of inhibitory and excitatory Complex-like subunits tends to increase with the RF *Simpleness*. It suggests that the Complex-like inputs to Simple cells consist of relatively more inhibition than those of Complex cells, consistent with earlier reports of a major role of intracortical inhibition in shaping Simple cell RFs (Sillito, 1975; Debanne et al., 1998).

Consistent with previous results obtained in Macaque V1 using STC techniques (Rust et al., 2005; Chen et al., 2007), we found that Simple and Complex RFs in cat V1 are not fully accounted for by the LN or the energy model respectively. Instead, our data show that, regardless of the *Simpleness* of their RFs, V1 cells receive a diversity of feature selectivity from their synaptic afferents; Simple and Complex RFs primarily differ with respect to the weight of their Complex-like subunits relative to the Simple-like component. The discrepancy between our data and previous STC analyses performed in cat V1 with DN or natural images (Touryan et al., 2002, 2005; Felsen et al., 2005) is most likely the result of an insufficient number of spikes relative to the dimensionality of the stimulus (as suggested earlier, see Rust et al., 2005). In our own data, we observed that the collection of RF subunits estimated from subthreshold responses was not readily retrievable from the Volterra kernels estimated directly from spiking responses. However, by mapping the spiking RF in the functional space defined by the RF subunits identified at the subthreshold level, we were able to reduce the dimensionality of the parameter space and assess the stimulus feature selectivity represented in the spike train. Our results showed that the Complex-like subunits of the subthreshold RF contribute to the spiking response in proportion to their relative weights. The functional diversity of Complex-like contributions observed in the synaptic afferents is thus expressed in the cell spiking output and transmitted to postsynaptic targets. This is consistent with recent results based on spiking responses to stimuli of lower dimensionality and showing, in cat V1, significant deviations of Simple cell RFs from the classic LN model (Levy et al., 2013). Cat and Macaque V1 are therefore comparable with respect to the computational richness of Complex-like nonlinearities encountered in Simple and Complex cell RFs (Rust et al., 2005; Chen et al., 2007).

It is noteworthy that the collection of RF subunits recovered in the present study probably still underestimates the functional space covered by the cell synaptic afferents: on the first hand, the signal-to-noise ratio of the response constrains the number of significant subunits that can be recovered in a given period of time; on the other hand, it is likely that the nonlinear contributions conveyed through recurrent or feedback connections of the cortical network are difficult to estimate with a catenary second-order RF architecture.

The diversity of feature selectivity measured in the subthreshold RF subunits is remarkably consistent with previous conductance measurements, also performed in cat V1, which demonstrated substantial excitatory and inhibitory conductance inputs for orientations far away from the cell preferred orientation (Monier et al., 2003, 2008). Our results further show that these synaptic contributions cover a functional space defined by multiple Complex-like components with distinct feature selectivity. In agreement with previous conductance measurements (Monier et al., 2003; Boudreau and Ferster, 2005), we also found a systematic temporal interplay between excitatory and inhibitory Complex-like subunits, inhibitory components arriving with an average delay of ∼20 ms relative to the excitatory components. Still, the exact functional role of this delayed Complex-like inhibition remains unclear. One possibility would be that it contributes to the dynamic sharpening of the orientation selectivity observed at the spike level in some V1 neurons (Ringach et al., 1997, 2002; Schummers et al., 2002; Shapley et al., 2003; Xing et al., 2005). It may also play a major role in adapting V1 neural code to the contrast of visual stimuli (Levy et al., 2013). As reported in V1 (Monier et al., 2008) and other sensory cortices (Wehr and Zador, 2003), this delayed inhibition contributes to sharpening the temporal precision of the spiking response, generating selective spiking opportunity windows for certain stimulus configurations (Pouille and Scanziani, 2001). Consistent with this view, in our data, the STAs of the subthreshold Complex-like components showed an almost systematic temporal sequence of excitatory and inhibitory contributions (Fig. 10*A*,*B*) in which the spike emission is preceded by an early increase of excitatory components concomitant to a decrease of inhibitory contributions and often followed by an inhibitory rebound. Although the estimated excitatory and inhibitory RF subunits cannot be strictly interpreted as pure conductance changes, this inferred sequence of an early excitation followed by a delayed inhibition at time of spike emission likely reflects equivalent changes in the sign of the balance between concomitantly activated excitatory and inhibitory circuits.

Although the Simple-like RF subunits may result from the *push-pull* arrangement of excitatory and inhibitory feedforward inputs selective for the same orientation, the diversity of feature selectivity expressed by the Complex-like RF subunits is not consistent with a strict iso-orientation preference rule for excitatory and inhibitory connections as generally posited (Ferster and Miller, 2000; Priebe and Ferster, 2012). Although the estimated Complex-like subunits do not necessarily correspond to the RFs of neurons presynaptic to the recorded cell, they bear a striking resemblance to the linear RF of V1 Simple cells, which suggests that they originate from within the cortex (Rust et al., 2005; Chen et al., 2007). Moreover, the distribution of preferred orientation that they reveal (Fig. 7*B*) is consistent with the distributions of excitatory and inhibitory axonal projections over adjacent orientation columns reported by anatomical studies (Kisvárday et al., 1997; Karube and Kisvárday, 2011). The diversity of orientation and spatial frequency preferences reported here thus supports the hypothesis that the Complex-like components of V1 RFs arise from lateral interactions between adjacent cortical columns. Under this assumption, our results are consistent with the proposal that the Simple or Complex nature of V1 RFs arises from the balance between feedforward and lateral connectivity (Chance et al., 1999; Tao et al., 2004).

Advancing our knowledge on the functional properties of second-order nonlinearities in V1 RFs may improve our understanding of how visual cortex processes visual cues. The diversity of Complex-like RF subunits recovered in the present study raises the question of the computational benefit of integrating over such a large repertoire of feature selectivity. Because features in natural images often combine multiple orientations and spatial frequencies in the same location, filtering the visual input through multiple Complex-like RF components, covering a large set of stimulus selectivity, may contribute to the detection of high-order correlations of the visual scene (Victor et al., 2013). The expression of adaptive gain control processes over this variety of Complex-like subunits could also account for contextual changes in the selectivity of V1 RFs (Fournier et al., 2011). Such distributed gain control would fit with adaptive shifts in orientation tuning curves (Müller et al., 1999; Dragoi et al., 2001; Felsen et al., 2002) or with long-lasting changes observed through imposed Hebbian supervision (Debanne et al., 1998). We propose that the wide functional spectrum of Complex-like synaptic contributions to Simple and Complex RFs constitutes a multicompetent substrate for adaptating V1 cells to the statistics of visual inputs.

## Footnotes

This work was supported by Centre National de la Recherche Scientifique, the Agence Nationale de la Recherche (Natstats and V1-Complex), and EC Grants BrainScales FP7-2010-IST-FETPI 269921 and Brain-i-nets FP7-2009-ICT-FET 243914. K.S. and Z.F.K. were supported by National Neuroscience Program (NAP) to Z.F.K. We thank Gérard Sadoc for invaluable technical assistance in developing stimulation software and kernel analysis library tools and Andrew Davison for comments and suggestions on the manuscript.

The authors declare no competing financial interests.

- Correspondence should be addressed to either Dr. Julien Fournier or Dr. Yves Frégnac, Unité de Neuroscience, Information et Complexité, Centre National de la Recherche Scientifique, Unité Propre de Recherche CNRS 3293, 1 Avenue de la Terrasse, Gif-sur-Yvette, F-91198, France, julien.fournier{at}unic.cnrs-gif.fr or yves.fregnac{at}unic.cnrs-gif.fr