## Abstract

Neighboring neurons in cat primary visual cortex (V1) have similar preferred orientation, direction, and spatial frequency. How diverse is their degree of tuning for these properties? To address this, we used single-tetrode recordings to simultaneously isolate multiple cells at single recording sites and record their responses to flashed and drifting gratings of multiple orientations, spatial frequencies, and, for drifting gratings, directions. Orientation tuning width, spatial frequency tuning width, and direction selectivity index (DSI) all showed significant clustering: pairs of neurons recorded at a single site were significantly more similar in each of these properties than pairs of neurons from different recording sites. The strength of the clustering was generally modest. The percent decrease in the median difference between pairs from the same site, relative to pairs from different sites, was as follows: for different measures of orientation tuning width, 29–35% (drifting gratings) or 15–25% (flashed gratings); for DSI, 24%; and for spatial frequency tuning width measured in octaves, 8% (drifting gratings). The clusterings of all of these measures were much weaker than for preferred orientation (68% decrease) but comparable to that seen for preferred spatial frequency in response to drifting gratings (26%). For the above properties, little difference in clustering was seen between simple and complex cells. In studies of spatial frequency tuning to flashed gratings, strong clustering was seen among simple-cell pairs for tuning width (70% decrease) and preferred frequency (71% decrease), whereas no clustering was seen for simple-complex or complex-complex cell pairs.

- clustering
- degree of tuning
- local circuitry
- tuning width
- visual cortex

neurons in the cerebral cortex of many mammals show columnar organization: certain response properties are similar among neurons across the depth of cortex at any given cortical position and change gradually with tangential movement across cortex (Hubel and Wiesel 1962; Mountcastle 1957; reviewed in Horton and Adams 2005; Van Hooser 2007). In cat primary visual cortex (V1), properties that show such organization include preferred stimulus orientation (e.g., Albus 1975; Berman et al. 1987; DeAngelis et al. 1999; Hetherington and Swindale 1999; Hubel and Wiesel 1962, 1963; Lee et al. 1977; Maldonado et al. 1997; Maldonado and Gray 1996; Ohki et al. 2005, 2006), preferred stimulus direction of movement (e.g., Berman et al. 1987; DeAngelis et al. 1999; Kim et al. 1999; Ohki et al. 2005; Payne et al. 1981; Tolhurst et al. 1981), and preferred spatial frequency of a drifting grating stimulus (e.g., DeAngelis et al. 1999; Issa et al. 2000; Kim et al. 1999; Maffei and Fiorentini 1977; Mallik et al. 2008; Tolhurst and Thompson 1982). Columnar organization implies a more easily measured property, spatial clustering, meaning that nearby neurons are on average more similar in these properties than pairs of neurons chosen from arbitrary locations.

Here we examine whether neurons spatially cluster by their degree of tuning, i.e., by the width of their orientation or spatial frequency tuning curves or the degree of direction selectivity, using single-tetrode recording in cat V1 to simultaneously isolate multiple cells at single recording sites. Better characterization of the structure of spatial clustering—which features are clustered, and how strongly?—can provide strong constraints for models of both the development and the function of columnar cortex. Models of circuit development typically explain clustering as arising from processes of “Hebbian” or correlation-based synaptic plasticity that lead nearby neurons to develop correlated responses and more distant neurons to develop uncorrelated or anticorrelated responses (e.g., Kaschube et al. 2010; Miller 1996), or else as inheritance from upstream structures (Paik and Ringach 2011; but see Hore et al. 2012; Schottdorf et al. 2014). Different models will have different predictions as to the clustering of response properties that will result. The function of columnar organization and more generally of spatial clustering of response properties is unclear (Horton and Adams 2005), although it has been argued that they may enable processing of stimulus features with minimal wiring length (Chklovskii and Koulakov 2004). Rodents and lagomorphs appear to lack spatial clustering (Van Hooser 2007), suggesting that there may be two different basic cortical designs, one involving columnar organization of derived response features (i.e., those beyond the topography established in the sensory periphery) and one lacking such organization. Both the differences in developmental mechanisms that lead to these different structures and the differences, if any, in their function remain unknown. Understanding of which properties are clustered, and how strongly, in animals with columnar organization will help inform studies of these differences.

There have been few previous studies of clustering of tuning strength. A handful of studies have found evidence for clustering in orientation tuning width (Hetherington and Swindale 1999; Nauhaus et al. 2008), but until recently (Martin and Schroder 2013) the strength of clustering had not been measured quantitatively, and numbers in that study were small (36 cell pairs). No significant clustering has been found in direction selectivity index (DSI) (DeAngelis et al. 1999; Martin and Schroder 2013) or spatial frequency tuning width (Martin and Schroder 2013), but numbers in those studies were small (20–50 cell pairs) so it is unclear whether clustering was absent or simply too weak to detect. We reexamine the clustering of all of these measures of degree of tuning with much larger numbers (many hundreds of cell pairs) for responses to both drifting and (except DSI) flashed gratings. In addition to assaying clustering of tuning strength, we also, for comparison, characterize the clustering of preferred stimuli by the same methods. Abstracts of this work have appeared previously (Emondi et al. 2000, 2006).

## METHODS

### Animal Preparation

#### Surgery and anesthesia.

All experimental recordings were conducted under a protocol approved by the University of California, San Francisco, Committee on Animal Research. We recorded from area 17 of 19 anesthetized male and female cats aged from 3 to 6 mo. Several hours before surgery, dexamethasone (0.5–5 mg/kg iv) was given to reduce anticipated cerebral edema. Initially the cat was anesthetized with 1.5–5% isoflurane. After intubation was performed, the isoflurane was discontinued and the animal was given pentobarbital (Nembutal, 25 mg/kg iv initial dose) for the remainder of the experiment. The animal was switched to external ventilation with an O_{2}-nitrous oxide (up to 1:2) mixture, with the overall flow chosen to keep lung pressure within physiological limits (8–12 cmH_{2}O). The animal was mounted in a stereotaxic apparatus, and a small craniotomy was made above the visual cortex. The dura was removed, and care was taken thereafter to keep the cortical surface in good physiological condition. The surface was kept moist and/or was protected by applying agar and covering the agar surface with silicone oil to maintain moisture. The animal was then paralyzed with a continuous intravenous infusion of gallamine (10 mg·kg^{−1}·h^{−1}) diluted with lactated Ringer solution and 2.5% dextrose (5–10 ml·kg^{−1}·h^{−1}). Every 6–12 h, atropine sulfate (0.04 mg/kg sc) was given to reduce tracheal secretions. Every 12 h, an antibiotic (ceftizoxime, 10 mg/kg sc) was administered to prevent infection. Contact lenses were used to protect the corneas and to focus the eyes on a tangent screen at a viewing distance of 30–40 cm.

#### Vital parameter monitoring.

A set of vital parameters—heart rate, EEG, respiratory rate, lung pressure, O_{2} saturation, expiratory CO_{2}, and body temperature—was continuously monitored throughout the entire experiment. The entire set of parameters was displayed numerically and represented graphically by vital monitor software, which also sampled those parameters once a minute and stored them in a database. Instantaneous values and long-term trends of heart rate, CO_{2} expiration level, and EEG spectrum were used to control proper anesthesia level, along with monitoring of heart rate responses to noxious stimuli (paw pinch).

#### Recording locations.

We recorded with one tetrode at a time. Most of our recordings were from penetrations down the medial bank of V1. This means that tetrodes typically ran roughly parallel to the cortical layers. As a result, we have no reason to think that we have uniformly sampled the cortical layers. In total, we recorded from 317 sites in 60 different electrode penetrations in 18 animals.

### Apparatus

#### Tetrodes.

We used different types of tetrodes to record from cat striate cortex. Some were made by Thomas Recording Scientific Resources (Giessen, Germany). Most were fabricated in our lab, made of either 12.5-μm NiCr or 7.5-, 12.5-, or 25-μm tungsten with recording tips gold plated, by a technique similar to Wilson et al. (1992) and Gray et al. (1995). In this procedure, a heat gun is used to meld the plastic coatings around the four wires. We used the heat gun extensively with the goal of preventing splaying of the tips, subsequently checking interelectrode impedances to ensure no shorts. Relatively frequent microscopic inspection of tetrode tips showed tips intact without splaying, although we cannot rule out that splaying may have occurred in some cases. In all tetrodes, single-electrode impedance varied from 0.7 to 1.4 MΩ measured at 1 kHz.

#### Data acquisition.

The tetrode was connected to a custom-made head stage amplifier (based on the INA110 chip by Burr-Brown) providing a gain of 10, DC coupled. The signal was further amplified and band-pass filtered by a CyberAmp 380 (Axon Instruments) with the following settings: gain of 1,000, AC coupling at 300 Hz, a Bessel-type fourth-order high-cut filter at 3,000 Hz, and a notch filter at 60 Hz. The voltage traces were sampled at 20 kHz with 12-bit resolution and streamed to disk. Recording to a file was initiated ∼1 s before the stimulus onset and was terminated ∼1 s after the stimulus presentation was completed.

#### Visual stimulator.

Gamma-corrected visual stimuli were shown on a FlexScan FX-E8 21-in. color display monitor (model MA-21A2, Nanao) with 120 Hz frame update rate. The monitor was calibrated and gamma-corrected in software. Synchronization between the data acquisition system and the visual stimulator running on a different computer was achieved by generating pulses in software at the start of each frame and recording them on a separate channel along with the four tetrode voltage traces at a common sampling rate.

### Visual Stimuli

Each site was driven by one or more of the following stimulus sets: drifting gratings sets, *S*_{orient} and *S*_{spat}, and flashed gratings, *S*_{flashed}.

#### Drifting gratings.

In the drifting gratings set of experiments, we studied responses of V1 neurons to drifting full-field (on a monitor that covered 43° × 55° of the animal's visual field), sinusoidally modulated luminance gratings. Two sets of stimuli, *S*_{orient} and *S*_{spat}, were used for these experiments. *S*_{orient} was used to estimate direction selectivity and orientation tuning at a site. It included drifting gratings (100% contrast) of 72 different directions of movement covering 0–360° in 5° steps, with spatial frequency 0.5 cyc/°. Using 5° steps gives essentially the same statistical power as using 10° steps with twice the number of repetitions but helps ensure that even the occasional very narrowly tuned cell will be well stimulated. *S*_{spat} was used to assess spatial frequency tuning. It consisted of drifting gratings (100% contrast) of 10 spatial frequencies, spaced approximately evenly on a logarithmic scale from 0.1 to 4.0 cyc/°, and of multiple orientations chosen according to the preferred orientations measured at the site using *S*_{orient} (see next paragraph). For both sets, temporal frequency was set to either 2 or 3 Hz. A block consisted of multiple repetitions of *S*_{orient} or *S*_{spat}, with all of the gratings in a block shown in pseudorandom order. Each grating was shown for 4 s, with successive gratings separated by a blank period (gray screen at same mean luminance as the gratings, 45 cd/m^{2}) of ∼1 s.

At every cortical site, we first tested whether the site responded reasonably to visual stimulation. For this we used a set of drifting gratings of different orientations (0–360°, 10° steps) with spatial frequency 0.5 cyc/° and temporal frequency 2 or 3 Hz. If the site was stable and showed reasonable visual responses, then orientation tuning was measured by showing a block of two to four repetitions of *S*_{orient}, for which each grating was presented for a few seconds (between 4 and 12 cycles). The entire block lasted between 6 and 15 min. To correctly measure spatial frequency tuning, it is important that a cell be studied at its preferred orientation (e.g., Issa et al. 2000). We did a preliminary online sorting and analysis of the data from the *S*_{orient} stimulus presentations to estimate the preferred directions of the various cells isolated at the site; this took <5 min. If preferred directions for any pair of extracted cells were within 5° of one another, then the cells were considered as sharing the same preferred direction. We then showed a block consisting of four repetitions of the *S*_{spat} stimulus set at each preferred direction (typically between 1 and 6) found at that site. The entire *S*_{spat} block contained 4 repetitions × 10 spatial frequencies × (number of different directions) gratings, each lasting 4 s followed by a 1-s blank period, and lasted from 7 to 20 min depending on the number of different grating directions.

#### Flashed gratings.

In the flashed gratings set of experiments, we studied the responses of V1 neurons to stationary sinusoidal luminance gratings (square 8° × 8° frames, centered on the site's receptive field), which we call *S*_{flashed}. The set of stimuli consisted of gratings of 36 orientations (in 5° steps), 10 spatial frequencies (evenly spaced on a logarithmic scale from 0.05 to 4.0 cyc/°), and either four or eight spatial phases (spaced 90° or 45° apart, respectively). This kind of stimulus ensemble is sometimes called a Hartley subspace stimulus set (Ringach et al. 1997). The monitor was updated with a new grating at a frame rate of 10 Hz, 24 Hz, or 60 Hz (monitor refresh rate was 120 Hz). The set *S*_{flashed} consisted of multiple presentations of this set of frames, each time in a pseudorandom order. To minimize the effect of electrode drift recording times were limited to 10–20 min, and so the full set of gratings was presented 4 times (for 10 Hz frame rate), 10 times (for 24 Hz frame rate), or 16 times (for 60 Hz frame rate).

### Data Processing

At the next stage, data were analyzed (off-line) in order to detect spikes from the continuous voltage traces and then sort the spikes from the unsorted multiunit pool into clusters corresponding to isolated single cells. The sorting of spikes was conducted separately for the responses to the *S*_{orient}, *S*_{spat}, and *S*_{flashed} blocks of stimuli. We did not attempt to draw correspondences between the cells studied with the different stimulus sets, that is, to determine which cluster in the responses to one type corresponds to the same cell as a given cluster in the responses to another type (Emondi et al. 2004). We analyzed the responses to the *S*_{orient} stimulus for orientation and direction tuning, the responses to the *S*_{spat} stimulus for spatial frequency tuning at the preferred direction, and the responses to the *S*_{flashed} stimuli for both orientation and spatial frequency tuning. The only correspondence that was assumed was that all of the preferred directions of the cells isolated in response to the *S*_{spat} stimulus block were actually shown in that block, that is, that new cells with different preferred directions had not appeared in the time between the showing of the *S*_{orient} block and the showing of the *S*_{spat} block.

#### Spike detection and sorting.

The detection and sorting of spikes were performed with custom-made MATLAB software.

##### SPIKE DETECTION.

The four-dimensional voltage trace *v*(*t*) (from the 4 channels of the tetrode) was high-pass filtered at 300 Hz with a 1st-order Butterworth filter. We calculated the 4 × 4 cross-channel covariance matrix, *C*, of this four-dimensional voltage data, using randomly selected time segments. We then marked a spike event whenever the Mahalanobis distance, , of the trajectory exceeded a threshold θ. We computed *C* from regions of the voltage trace without spikes, as follows. We first computed an initial *C* using the full voltage trace and then marked potential spikes using a conservative threshold of θ = 5 and omitted all segments of *v*(*t*) from 1 ms before to 3 ms after each potential spike. We then recomputed *C* from the trace with potential spikes omitted and finally marked spikes using this *C* with θ = 8.

For each spike, we recorded its time of occurrence as the time at which the Mahalanobis distance reached a maximum during that particular spike, as well as the surrounding waveform from 0.9 ms before the negative peak on each channel until 1.2 ms after it (43 samples from each channel at 20 kHz). The waveforms from each channel were upsampled by a factor of 10 using Fourier interpolation with the surrounding 80 samples, aligned by the negative peak amplitude on each channel, and then downsampled again.

##### SPIKE SORTING I: FEATURE EXTRACTION, CLUSTERING.

For each spike, we concatenate the four waveforms from each channel, creating a 43 × 4 = 172-dimensional vector. Since the voltage signals in a tetrode recording are highly correlated across channels, we performed “cross-channel whitening” to transform these vectors to a basis in which the redundancy across the four channels was eliminated (Emondi et al. 2004). This means, in essence, that differences between voltages in any direction in the four-channel space are always measured in units of the intrinsic variability in that direction. We then used the graph-Laplacian feature (GLF) algorithm (Ghanbari et al. 2011) [a modified version of principal components analysis (PCA), designed for clustering applications such as spike sorting] to reduce the dimensionality of the spike vectors from 172 dimensions down to 8. We used this algorithm with k (the parameter that determines number of nearest neighbors calculated for each spike) set to 15. These eight-dimensional spike vectors were sorted into clusters automatically with the KlustaKwik program (klustakwik.sourceforge.net), which fits a Gaussian mixture model to a distribution of data points (spikes). We ran the program with most of the default parameters, except that we set minClusters = 10 and nStarts = 5. This results in a larger number of random initializations (105 instead of the default of 11), which increases the probability of finding the cluster arrangement with a globally maximum likelihood.

##### SPIKE SORTING II: CLUSTER “PRUNING.”

To “clean up” the clusters and remove contaminating spikes from other cells, we reduced the size of the clusters by eliminating spikes that violate the cell's refractory period (i.e., they occur <1–2 ms from another spike in the cluster). These pairs of spikes that violate the refractory period are indicators of the presence of spikes from multiple cells, and so reducing the size of the cluster so that one of the two spikes in each pair is removed may reduce the contamination from other cells. For purposes of pruning, we represented each spike in the four-dimensional channel-whitened voltage space described above. The pruning was done by cutting this space of spikes with a hyperplane chosen to eliminate one refractory-violating spike while removing as few spikes as possible from the cluster and repeating this process until refractory violations were eliminated. If this procedure eliminated more than a third of the spikes in the cluster, the cluster was discarded. This procedure focuses on removing spikes that are as far as possible from the main densities of spikes in the cluster, since the main density would be most likely to come from a single cell. In tests of this algorithm against a data set of tetrode recordings in which “ground truth” was known for one intracellularly recorded cell (using the data set from Harris et al. 2000), we have found that it performs well, eliminating a far higher percentage of the spikes that did not belong to the cell than of spikes that did come from the cell. We describe the algorithm and these tests more fully in a separate paper.

##### SPIKE SORTING III: CLUSTER MERGING.

The KlustaKwik spike sorting program fits a Gaussian mixture model to the set of spike vectors and typically tends to “overcluster” the spikes (i.e., we usually get many more clusters than there are cells). This happens because the clusters of spikes are typically not exactly Gaussian shaped, but rather have longer tails (Harris et al. 2000). To complete the spike sorting, it is necessary to inspect the output clusters, identify which of the output clusters came from the same cell, and update their cell labels accordingly. This inspection and merging of the clusters was performed manually with custom-built MATLAB software, described in unrefereed supplemental materials, Supplemental Section S1 (see endnote).

We computed a measure of clustering quality, the isolation distance (ID) (Schmitzer-Torbert et al. 2005), for each isolated cell. For a cluster with *N* spikes, the ID is the square of the Mahalanobis distance of the *N*th closest noncluster spike. We only used cells that were reasonably well isolated as judged by an ID of at least 10. A total of 87% of our isolated cells satisfied this criterion. Of our isolated cells, the median ID was 26.0, with a 5th–95th percentile range of 12–197.

##### INTERCLUSTER GAUSSIAN OVERLAP.

To test whether neurons that showed more similar response properties showed more similar spike shapes (a possible sign of contamination in spike sorting), we computed an overlap measure between the spike clusters corresponding to different cells. We first fit a multidimensional Gaussian to each cluster (by calculating the sample mean and covariance over all its points in the 8-dimensional space in which they were originally clustered). For two clusters, we take their sample means (*M*_{1} and *M*_{2}) and their covariances (*C*_{1} and *C*_{2}) and calculate the integral *I*_{1,2} = *I*(*M*_{1},*C*_{1},*M*_{2},*C*_{2}) of the product of the two Gaussians:
where

The further apart the two clusters are, the smaller this integral will be. To factor out the effect of size and elongation of each cluster, we normalize this integral by dividing by the integral we would have obtained if the clusters had had the same mean:

The Gaussian overlap, *D*_{1,2}, between two clusters is defined as the negative log of the normalized integral:

### Measures of Response

A cell's response to a given stimulus was taken to be the average firing rate to the stimulus, minus the average spontaneous firing rate, with negative values set to zero (since some of our measures of tuning are only defined for nonnegative values). While some previous works have used spontaneous-subtracted response measures, as we do (e.g., Gizzi et al. 1990; Gur et al. 2005; Swindale et al. 2003), others, particularly those who have studied a “global” measure of orientation tuning (circular variance), have included spontaneous in the response measure (e.g., Alitto and Usrey 2004; Ringach et al. 2003). Inclusion of spontaneous shifts global orientation tuning widths and DSIs toward broader tuning but has very minor effects on measures of clustering. This is described further in unrefereed Supplemental Section S2.1 (see endnote). Hereafter, we will use “firing rate” to mean this spontaneous-subtracted, rectified firing rate unless otherwise indicated.

#### Drifting gratings.

For drifting gratings, the response to a particular stimulus orientation and spatial frequency was calculated as the mean firing rate during the period the grating was being presented. In studying the response of cells to drifting gratings, we were interested in steady-state responses (as opposed to flashed grating studies, in which we were studying transient responses). The first cycle often contained a response to the sudden increase in contrast (from a uniform screen at mean luminance to a drifting grating at full contrast) consisting of an increase in firing rate in the 30–60 ms after the beginning of the first cycle, so we omitted the first cycle of each drifting grating presentation from our analysis.

#### Flashed gratings.

For the flashed grating experiments, calculating the response was not as straightforward. Since stimuli of different orientations and spatial frequencies were interleaved, it was important to determine which flashed grating frame(s) was most likely to be responsible for evoking a particular spike, particularly in the experiments where the frame rate was 60 Hz (and each stimulus was displayed for only 17 ms). Since the strongest response to a preferred grating typically occurs in a window that can lie anywhere in a window from 30 to 80 ms after the stimulus is first displayed on the screen (and this epoch is different for each cell), we first had to determine the extent of this window so we could correctly attribute spikes to the corresponding stimuli that caused them. For simplicity, we used a single window for all stimuli for a given cell, even though in many cases the responses to different stimuli had slightly different latencies (for example, cells tend to have shorter latencies to gratings of lower spatial frequencies; Bredfeldt and Ringach 2002; Mazer et al. 2002).

To determine the extent of the optimal poststimulus window, we used the following algorithm, modeled after that of Mazer et al. (2002). The spikes following the onset of each stimulus were binned into 5-ms bins. To overcome the relatively high amount of sampling noise with short (5 ms) bins, we considered responses (spike count) in a set of five bins centered on a given poststimulus bin. Recall that each grating stimulus is presented multiple (between 4 and 16) times during the experiment. For each set of five bins, we constructed *r*_{odd}, the set of responses to all stimuli averaged over the odd-numbered presentations, and *r*_{even}, averaging over the even-numbered presentations. We called the central time bin “reproducible” if *r*_{odd} and *r*_{even} were significantly correlated (*P* < 0.001 with a 1-sided Pearson correlation test). The window for a particular cell was then simply the concatenation of all reproducible time bins.

We wanted our algorithm to select windows that consisted of a contiguous set of time bins. Usually, all the reproducible bins formed a single continuous sequence, so this was not an issue. However, on occasion we came across cells in which there were two separate sequences of reproducible time bins, sometimes relatively far apart (>15 ms). In such cases, we used the reproducibility of the entire window (as opposed to that of individual time bins) to arbitrate between three possible windows: the window consisting of *1*) the first set of bins, *2*) the second set of bins, or *3*) both sets of bins, including the nonreproducible time bins in between them. We chose whichever of these three windows was the most reproducible (i.e., produced the highest correlation coefficient between *r*_{odd} and *r*_{even}). Across all cells, the start times of the window (relative to the stimulus onset) were 33 ± 12 ms (mean ± SD) with 5/95 range (5th and 95th percentiles) of 21–55 ms. The width of a cell's response window was strongly correlated with the duration that each grating was present on the screen (i.e., 1/frame rate). For cells responding to gratings flashed at 60 Hz (frame duration 16.7 ms), window widths were 46 ± 20 ms with 5/95 range of 21–104 ms. For 24-Hz gratings (frame duration 41.7 ms), window widths were 60 ± 26 ms with 5/95 range of 25–113 ms. For 10-Hz gratings (frame duration 100 ms), window widths were 106 ± 35 ms with 5/95 range of 35–158 ms.

Since our window selection method can select a window that is longer than the duration of each stimulus, a concern is that spikes collected in response to a particular stimulus might have been driven by stimuli presented before or after. In practice, this is not an issue, both because the selection for reproducibility will tend to exclude stimulus bins with significant contamination from other stimuli and because the flashed grating stimulus space covers a large range of orientations, spatial frequencies, and phases while each cell typically responds to only a small subset of these stimuli. Thus the stimuli before and after any preferred stimulus are usually one of the cell's nonpreferred stimuli, leaving the response to preferred stimuli mostly uncontaminated by responses to other stimuli. However, to address this concern directly, we reran our analysis and constrained the window to be no longer than the length of presentation of each stimulus (as before, selecting its position to maximize reproducibility of the firing rates for each stimulus between odd/even trials). All response properties were essentially unchanged: the median changes in orientation and spatial frequency tuning widths (defined below) were all <1%. Thus allowing the analysis window to be longer than the frame presentation time does not seem to have affected our results.

### Studies of Orientation and Direction Tuning

#### Computing orientation and direction tuning curves.

We first studied the orientation tuning properties of cells responding to drifting gratings (using responses to the *S*_{orient} stimulus block) and to flashed gratings (using responses to the *S*_{flashed} block). For each cell, we computed the orientation (flashed gratings) or direction (drifting gratings) tuning curve, *r*_{k}, where *r*_{k} represents the firing rate elicited by gratings of the *k*th orientation or direction, θ_{k}, averaged over cycles (drifting gratings) or spatial phases (flashed gratings). For drifting gratings, direction tuning was measured at only one spatial frequency. However, for flashed gratings we have an orientation tuning curve for each of the 10 spatial frequencies used. Since preferred orientation is largely independent of stimulus spatial frequency (Issa et al. 2000; Webster and De Valois 1985), we averaged over the responses to all spatial frequencies when estimating the preferred orientation. Orientation tuning width, on the other hand, does depend on stimulus spatial frequency. We thus used the responses of the cell at its preferred spatial frequency when measuring orientation tuning width.

#### Cell selection criteria and the measure of preferred orientation.

Before using the computed orientation tuning curves for the estimation of orientation/direction tuning parameters, we used a set of three criteria to determine which cells had reliable responses to different orientations. Our criteria were that *1*) the cells showed orientation selectivity, *2*) they gave reproducible responses to the presented stimuli, and *3*) they were well fit by a Gaussian-shaped orientation tuning curve. This selection process is necessary because the tetrode samples all separable cells at a given site regardless of their responsiveness, but tuning could only be characterized in cells that showed reproducible and selective responses. In single-cell recordings, similar exclusion is likely to be done in the process of isolating a cell for study.

##### ORIENTATION SELECTIVITY.

To check whether the response of a cell was orientation selective, we first calculated the cell's preferred orientation. We took *r*_{k}, the vector of average firing rates in response to each of the *N* stimulus directions θ_{k}, and used it to calculate the orientation resultant vector, *R*, defined as
(1)
i.e., we associate each response *r*_{k} to stimulus direction θ_{k} with a vector of length *r*_{k} and orientation 2θ_{k} and calculate the vector sum of these vectors. Defining *R*_{x} and *R*_{y} as the components of *R*, the preferred orientation was then taken to be the half the angle of the resultant vector:
(2)
The angles are initially doubled so that opposite directions correspond to the same angle and are added together, while directions that are 90° apart correspond to opposite angles and are subtracted from each other. The division by 2 in *Eq. 2* is to undo the initial doubling.

We then determined whether the preference for this particular orientation was statistically significant, as follows. We calculated the magnitude of the orientation resultant vector, *s*_{0} = |*R*|. If a cell is tuned for orientation, the clustering of relatively high firing rates around the preferred orientation will result in *R* having a large amplitude in the direction 2θ_{k}. Conversely, if the cell is not tuned for orientation, the responses of the cell will be randomly distributed across the vector *r*_{k} and will mostly cancel out in their contribution to *R*, resulting in *R* having a small amplitude in a random direction. We used a Monte Carlo randomization test to determine whether *s*_{0} was significantly larger than what would be expected in the absence of tuning: the vector *r*_{k} was randomly shuffled, a new orientation resultant vector *R*_{rand} was calculated for this randomized *r*_{k} with *Eq. 1*, and we calculated its magnitude, *s*_{rand} = |*R*_{rand}|. This process was repeated *N*_{rep} = 10,000 times, yielding *N*_{rep} samples of *s*_{rand}. If there is some degree of orientation tuning, then *s*_{0} should be larger than what would be expected from a random rearrangement of *r*_{k}. The probability of obtaining the observed value *s*_{0} under the null hypothesis of no orientation tuning is given by
(3)
where *L* is the number of values of *s*_{rand} that are ≥*s*_{0}. [This formula for the probability is known as Laplace's rule, see Jaynes (2003), Chapter 18; for a simple derivation see Erwin and Miller (1999), Appendix A]. We accept the response of the cell as orientation selective if *P* < 0.01, chosen as a conservative criterion that should reject marginal cases.

##### ORIENTATION REPRODUCIBILITY.

To assess the reproducibility of the orientation tuning curves, we computed *r*_{odd} and *r*_{even}, the orientation tuning curves averaged over the odd- and even-numbered trials. A “trial” is the firing rate averaged over one cycle (drifting gratings) or one presentation of all phases (flashed gratings). We smoothed *r*_{odd} and *r*_{even} by convolving with a Gaussian of standard deviation 5°. Under the null hypothesis that the response to each orientation is not reproducible, we would expect that *r*_{odd} and *r*_{even} would be uncorrelated. We thus rejected the null hypothesis (and defined the cell as having reproducible responses) if the correlation coefficient between *r*_{odd} and *r*_{even} was positive and the one-sided *P* value was <0.01.

##### GOOD FIT.

The third criterion for acceptance was that the cell's orientation tuning curve was reasonably well fit (*R*^{2} > 0.4) by a Gaussian (see *Eq. 5* below). One of our measures of orientation tuning width was the width of the Gaussian that fit the orientation tuning curve; use of this measure of width implicitly assumes that the Gaussian tuning curve is a good description of the orientation tuning curve. In practice, most of the cells that had reliable and single-peaked tuning (multipeaked tuning might indicate poor clustering) were well fit by a Gaussian.

From a total of 968 (470) cells responding to drifting (flashed) gratings, 729 (317) cells passed the orientation selectivity criterion, 856 (339) cells passed the reproducibility criterion, and 720 (292) passed both of these criteria. The orientation tuning curves of 620 (232) of these orientation-selective and reproducible cells also had good fits to a Gaussian tuning curve. Thus a total of 64.0% (49.4%) of our cells passed all three selection criteria and were accepted for further analysis. These 620 (232) cells were recorded from 246 (101) sites, across 55 (24) tetrode penetrations, from 17 (8) animals.

These numbers are comparable to, although slightly lower than, those found in other tetrode recording studies in V1, which have found between 69% and 77% of cells to have reliable responses and good orientation tuning using drifting gratings (Gray et al. 1995; Hetherington and Swindale 1999; Maldonado et al. 1997).

#### Measures of orientation selectivity.

For orientation tuning width, we considered both a “global” and a “local” measure. A “global” measure, such as the circular variance (e.g., Alitto and Usrey 2004; Ringach et al. 2003), is one that is sensitive to all responses in the tuning curve. A “local” measure, conversely, is one that is sensitive only to the width of the tuning curve peak.

As a global measure of orientation tuning width, we use the ordinary (noncircular) standard deviation of the distribution of responses *r*_{k}:
(4)
where *d*(θ_{k} − θ_{pref}) is the shortest distance around a circle of 180°. In *Eq. 4*, the tuning curve *r*_{k} is treated as a distribution over orientations, with a mean of θ_{pref}. This measure is essentially equivalent to the circular variance but, as an ordinary standard deviation expressed in degrees, has a more intuitive meaning. The circular variance CV can be translated into a circular standard deviation in degrees σ as 180/2π (Batschelet 1981; Fisher 1996) (the denominator is divided by 2 to account for orientation being circular over 180° rather than 360°). We have found empirically for our data that this σ is almost perfectly correlated with *w*_{ORI}^{Global}: *σ* = 0.76*w*_{ORI}^{Global} + 2.59, *R*^{2} = 0.99.

As a local measure of the width of the peak of the tuning curve, we used the standard deviation σ of the Gaussian in a fit of the orientation tuning curve to a Gaussian plus a constant term:
(5)
Here θ is restricted to the range θ_{pref} ± 90°. The parameters *A*, *B*, and σ were simultaneously fit to minimize the least-square error between *f*(θ) and the orientation tuning curve. When fitting, *B* is constrained to be at least as large as the smallest value of *r*_{k} in the range θ_{pref} ± 90°. The constant term, *B*, in *Eq. 5* absorbs any response far from the preferred, and so the Gaussian ends up fitting the peak without much effect of far-from-peak responses. We refer to the local measure of orientation tuning width, σ, as *w*_{ORI}^{Local}. An example of the set of orientation tuning curves at a recording site (with an illustration of the orientation tuning widths) is shown in Fig. 1.

Because the difference between global and local measures involves responses far from the peak, we also tried to directly examine clustering of responses far from the peak (at 90° from the preferred). However, spontaneous activity is clustered and this contaminates these measurements, so we were not able to draw firm conclusions about clustering of such responses (described further in unrefereed Supplemental Section S2.2; see endnote).

#### Preferred direction and direction selectivity index.

For responses to drifting gratings, we defined the preferred direction as one of two values, θ_{pref} or θ_{pref} + 180°, depending on which of the two directions gave a greater “integrated response.” We define the integrated response *N*_{θ} associated with a specific direction θ to be the sum of the responses for all directions in the range (θ − 90°, θ + 90°), where the bounding directions θ − 90° and θ + 90° are excluded from the sum. The direction selectivity index (DSI) was then defined according to the formula
(6)

where *N*_{pref} is the integrated response associated with the cell's preferred direction and *N*_{opp} is the integrated response associated with the opposite direction. A DSI close to 0 indicates roughly equal responses to both the preferred direction and its opposite, while a DSI close to 1 indicates that the cell responds almost exclusively to the preferred direction.

#### Multiunits.

In addition to studying the individual cells at a site, we also considered the multiunit activity consisting of all spikes detected at a site that were not assigned to a cell (but excluding spikes that were “pruned” from clusters; see *Spike detection and sorting*). This was typically dominated by the many spikes too small to be sorted into clusters corresponding to distinct cells. We wanted to measure differences in preferred orientation between cells and the site, so we used the preferred orientation of the multiunit activity as an indicator of the preferred orientation of the site. We used the multiunits if they passed the orientation selectivity and reproducibility criteria, but we did not impose the third (goodness of fit to a Gaussian) criterion, since *1*) the multiunit response is made up of the responses of multiple cells, which may not fit a single Gaussian, and *2*) we were not trying to measure the tuning width of the multiunit response.

### Studies of Spatial Frequency Tuning

#### Response measures and cell selection criteria.

For each isolated cell studied with the *S*_{spat} (drifting gratings) and *S*_{flashed} (flashed gratings) stimulus sets, a spatial frequency tuning curve was computed. According to the experimental design, every cell had multiple spatial frequency tuning curves (at each of multiple orientations or directions). For each cell, we selected the spatial frequency tuning curve at the cell's preferred orientation/direction. For the drifting grating stimulus set, spatial frequency tuning was measured for each of the handful (1–6) of directions that were determined to give good responses, and so preferred direction was defined as the direction that had the highest firing rate averaged over spatial frequencies. For the flashed gratings, we had a spatial frequency curve for all 36 orientations, and so the preferred orientation was calculated with *Eq. 2*.

From the pool of selected spatial frequency tuning curves, we accepted for further analysis only those that met two criteria: *1*) they were reproducible across repetitions and *2*) they were well fit by a skewed lognormal function (*Eq. 7* below).

Reproducibility was assessed in a similar way as for the orientation tuning data: we computed *r*_{odd} and *r*_{even}, the spatial frequency tuning curves averaged over the odd- and even-numbered trials (presentations of the stimulus, averaged across phase). We then calculated the correlation coefficient between these two vectors and considered the spatial frequency tuning reproducible if the correlation coefficient was positive with a one-sided *P* value < 0.01.

Goodness of fit to a spatial frequency tuning curve was defined as having an *R*^{2} > 0.4 for the fit of the tuning curve to the skewed lognormal function (*Eq. 7* below) and having the peak of the fitted spatial frequency curve lie within the range of stimulus spatial frequencies (see *Measures of preferred stimuli and selectivity*).

Of the 557 (470) cells whose spatial frequency tuning was studied in response to drifting (flashed) gratings, 475 (204) had reproducible tuning curves and 444 (185) of these were well fit by a skewed lognormal function. Thus a total of 79.7% (39.4%) of our recorded cells passed our spatial frequency tuning criteria and were used for analysis. These 444 (185) cells were recorded from 144 (87) sites, across 45 (18) tetrode penetrations, from 16 (8) animals.

#### Measures of preferred stimuli and selectivity.

For cells that passed the selection criteria above, the preferred spatial frequency and the spatial frequency tuning width were obtained from the fit of the spatial frequency tuning curve to what we call a skewed lognormal (SLN) function. The SLN has the following form, which has sufficient flexibility to capture the essential features of the real tuning curves (e.g., either right- or left-sided skewness):
(7)
This functional form resembles a skewed version of a standard lognormal distribution with the exception that it might take negative values; examples of the fit are shown in Fig. 2. Here *r*_{max} is the maximum response, *r*_{bkg} characterizes a DC level (untuned component of response), *f* is spatial frequency, *f*_{pref} is the preferred spatial frequency, *w* is a width parameter, and *s* determines the degree of skewness for the curve. We tried a number of methods of characterizing the spatial frequency tuning curves and found all to give roughly consistent results, but fitting of the SLN curve seemed most robust.

Note that the SLN function is unchanged if the sign of both *w* and *s* is switched, so, without loss of generality, we constrained *w* to be positive; then, positive and negative values of *s* correspond to right-and left-sided skewness, respectively. To ensure that the curve fitting provided a plausible spatial frequency tuning curve, we constrained *r*_{max} to be at most 1.5 times the maximum response in the tuning curve. This constraint prevented the curve from taking on unrealistically large values to fit the peak of the tuning curve. Most of the time (more than ∼90% of the time), this constraint does not affect the fitted parameters. The example in Fig. 2*E* illustrates when this constraint is necessary: for tuning curves where there are very few (1 or 2) samples around the peak and/or the fitted value of *f*_{pref} falls close to midway between two sampled spatial frequencies, an unconstrained fit can sometimes take on unrealistically large values. In cases such as these, enforcing the constraint on *r*_{max} yields a more plausible fit to the tuning curve.

For four cells responding to *S*_{spat} that otherwise had reproducible tuning curves and were well fit by the SLN curve (by the *R*^{2} measure), the fitted value of *f*_{pref} (the peak of the fitted SLN curve) fell outside of the range of the spatial frequencies tested. Since we could not get a good sampling of the responses of these cells close to and around their preferred spatial frequencies, we excluded them from our analysis of preferred spatial frequency and spatial frequency tuning width.

The spatial frequency tuning width in octaves was given by *w*_{SF} = log_{2}(*f*_{hi}/*f*_{lo}), where *f*_{hi} and *f*_{lo} are the spatial frequencies higher and lower than *f*_{pref} that gave half-maximal responses in the fitted SLN curve. This is highly correlated with the width parameter *w* of the SLN: *w* = 0.39*w*_{SF} + 0.01, *R*^{2} = 0.98.

### Studies of F1-to-DC Ratio and Simple/Complex Cell Classification

To classify cells as simple or complex, we studied the F1-to-DC ratio (F1/DC) either of the curve of cycle-averaged firing rates vs. time during the stimulus cycle (drifting gratings) or the phase tuning curve (flashed gratings). Here, the F1 is the first harmonic, the amplitude of the Fourier transform of the curve at a frequency of one cycle per curve length, and the DC (or F0) is the zeroth harmonic, which is proportional to the mean response over the curve. Cells were classified as simple or complex if the F1/DC was >1 or <1, respectively (Skottun et al. 1991; but see Kagan et al. 2002; Martinez et al. 2005; Mechler and Ringach 2002; Priebe et al. 2004 for arguments against the adequacy of this classification).

Ideally, cells should be classified as simple or complex on the basis of their responses to their preferred orientation and spatial frequency. We do this for cells responding to flashed gratings (*S*_{flashed}) and for studies of spatial frequency tuning in cells responding to drifting gratings (*S*_{spat}). However, for studies of orientation tuning and direction selectivity to drifting gratings (*S*_{orient}), cells were studied at only a single spatial frequency, 0.5 cyc/°. This is an intermediate spatial frequency that drives most cells reasonably well. Using the responses of cells to flashed gratings, we estimated the “misclassification rate” using the F1/DC of the response at 0.5 cyc/° instead of the preferred spatial frequency for each cell and obtained a modest 9% error rate. Thus we proceeded to use the F1/DC of cells responding to *S*_{orient} (at their preferred orientation but at 0.5 cyc/°) to classify them as simple or complex.

### Jackknife Procedure for Computing Standard Errors

We used a jackknife procedure (Miller 1974) to determine the standard errors of our estimates of preferred stimuli and tuning widths. For a particular stimulus ensemble with *M* trials, we calculated the quantity (preferred stimulus or tuning width) *M* times, each time leaving out one of the trials, obtaining *M* estimates θ_{m} of the given quantity, whose overall estimated value was θ. The standard error of the estimate is then given by
(8)
In the case of preferred orientation, (θ_{m} − θ) is replaced by *d*(θ_{m} − θ), the shortest distance around a circle of 180° (for orientation) or 360° (for direction) between θ_{m} and θ. For spatial frequencies, we use the difference between the two spatial frequencies in octaves: log_{2}θ_{m} − log_{2}θ.

To ensure that studied response properties were well characterized, we excluded a measured response property of a cell from study if its standard error δθ was at or above a conservative threshold (results were extremely similar in all respects whether or not this exclusion was done). For preferred orientation and orientation tuning width, we used a threshold of 5°. For *w*_{ORI}^{Global}, this excluded 6%/9% (for drifting/flashed gratings) of cells (starting from the set of cells that passed cell selection criteria); for *w*_{ORI}^{Local}, 15%/13%; and for preferred orientation, 4%/4%. For DSI, we used a threshold of 0.1, which removed 13% of cells. For preferred spatial frequency and spatial frequency tuning width, we used a threshold of 0.5 octaves, which removed 6%/2% (preferred) or 22%/16% (width) of cells. For studies of F1/DC and for classifying cells as either simple or complex, we used a standard error threshold of 0.25, which excluded 19%/27% of cells for studies of orientation tuning and 8%/23% of cells for studies of spatial frequency tuning (before cells were additionally excluded for excessive jackknife error in the measured tuning properties). We included or excluded cell responses individually for each measure, so that a given cell might be included for some measures and excluded for others.

### Measures of Clustering

We measured the degree to which a given response property (for example, orientation tuning width) tends to spatially cluster. For a particular response property, *x*, we first calculated the pairwise differences in *x* between all pairs of cells. As above (*Eq. 8*), “difference” is the shortest distance around a circle of 180° for pairs of orientations or around a circle of 360° for pairs of directions; the difference in octaves for pairs of spatial frequencies; and the absolute value of the difference for other quantities. We can then determine *P*_{W}(*x*), the distribution of differences in *x* among pairs of cells recorded at the same site (the “within-site” distribution), and *P*_{B}(*x*), the distribution of differences among pairs of cells recorded from different sites (the “between-site” distribution). As in DeAngelis et al. (1999), our measure of clustering was the median ratio, the median of *P*_{B}(*x*) divided by the median of *P*_{W}(*x*).

To determine the precision of our estimates of median ratios, we also calculated a bootstrap estimate of the standard error of each median ratio, expressed as a 68% confidence interval around the estimated value (Efron and Tibshirani 1986). For a given measure, the within-site and between-site distributions consist of *N*_{W} and *N*_{B} pairwise differences, respectively. A single bootstrap sample of the median ratio is obtained by taking *N*_{B} samples (with replacement) from the between-site distribution and *N*_{W} samples (with replacement) from the within-site distribution and taking the ratio of the medians of these two sets of samples. We calculated 1,000 such sample median ratios and determined the standard error margin that bracketed the central 68.3% of these samples. To avoid problems with the different samples of pairwise differences not being independent (because a single cell might contribute to multiple pairs), we did not use these standard errors to test significance of clustering but instead used a randomization test, described next.

### Tests for Significance of Clustering Measures

A Monte Carlo randomization procedure was used to assay the statistical significance of clustering, similar to the cluster index introduced by DeAngelis et al. (1999). Suppose we have calculated *P*_{W}(*x*), the within-site distribution of the pairwise differences for a particular property *x*, and its median *m*_{W}. To test the null hypothesis that there is no spatial clustering, we randomly permuted the cells among the recording sites (leaving the number of cells at each site unchanged), recalculated the distribution of the within-site difference measure using the new rearrangement of cells, and computed the median again for this new randomized distribution. We repeated this randomization process *N*_{rep} = 10,000 times, obtaining a set of *N*_{rep} randomized “within-site” distributions, {*P*_{W-rand}(*x*)}, and a set of medians from each one, {*m*_{rand}}. We then computed the probability *P*, under the null hypothesis of no spatial clustering, of finding a median of *P*_{W}(*x*) as small as (or smaller than) that observed: *P* = (*L* + 1)/(*N*_{rep} +2), where *L* is the number of random medians in the set {*m*_{rand}} that are ≤*m*_{W} (Laplace's rule, *Eq. 3*). The median ratios that we report are the ratio of the true between-site median to the true within-site median (without randomization), while the *P* values are the probabilities computed using this randomization procedure. The use of *N*_{rep} = 10,000 gives our significance tests a “resolution limit” of *P* ≈ 1/*N*_{rep} = 10^{−4}.

Our randomization procedure (randomly shuffling cells across sites and then computing the new “within-site” distribution) is conceptually very similar to the procedure of DeAngelis et al. (1999), who built the randomized within-site distribution simply by repeated draws of *N*_{W} samples from the between-site distribution *P*_{B}(*x*), where *N*_{W} is the number of pairs in the within-site distribution. Our method differs in controlling for the fact that cells at sites with more cells contribute to more pairs in the within-site distribution than do cells at sites with fewer cells, by shuffling the cells across sites without changing the number per site and then recomputing the within-site distribution. In practice, this was not an issue in the study of DeAngelis et al. (1999), because they recorded with a single electrode and almost always had the same number of cells (2) at a site (at only 2 sites did they have 3 cells). In contrast, with tetrodes we recorded a mean of ∼2.5 cells per site, with 4 or more cells in ∼20% of sites and as many as 7 or 8 in ∼5% of sites. Thus it was important in our case to use the random-cell-shuffling control.

Given that differences between within-site and between-site distributions are visible in cumulative distributions (see Figs. 4, 5, 7, 8, 10, and 11), one might think that a natural test of significance would be the Kolmogorov-Smirnov (KS) test, which tests whether two distributions significantly differ based on differences in their cumulative distributions. While a KS test does show significant clustering when our randomization test does, its application to our data is not valid. The KS test requires that the different data points be independent samples from the respective distributions, a condition that does not apply because different pairwise differences involving the same cell cannot be said to be independent. The randomization test avoids this problem, as it randomizes while preserving the statistics of the number of pairs contributed to by single cells.

For studies of simple and complex cells, letting “S” stand for simple and “C” for complex, we computed the median ratios for each studied response property separately for SS pairs, CC pairs, and SC pairs. To calculate the significance of the median ratios for each of these pair types, we performed *N*_{rep} = 10,000 randomization trials, where in each randomization the simple cells were randomized among themselves and the complex cells were randomized among themselves across sites. For each pair type, we set *L* to be the number of randomized median ratios of that pair type that were as large or larger than the one observed and then obtained the *P* value using Laplace's rule as above. We also determined the significance of the differences in ratios for a given response property between any two of these three pair types (SS, SC, or CC), under the null hypothesis that all three pair types showed equal clustering of the given response property, by using the following randomization test. We shuffled the labels “simple” and “complex” randomly across cells *N*_{rep} = 10,000 times, leaving all other cell response properties untouched. This preserves the degree of clustering in the within-site distribution while randomizing which cells are labeled simple or complex, and thus which pairs are SS, SC, and CC. For each randomization, we computed the difference in median ratios between the two pair types. Setting *L* to the number of randomization trials in which the difference in ratios had an absolute value that was as large as or larger than the absolute value of the observed ratio difference, we then computed the significance level, using Laplace's rule as above.

### Spatial Extent of Clustering

To study the spatial extent of clustering, we tested the significance of clustering of a response property as a function of distance between recording sites within a tetrode penetration. We collected the between-site, within-penetration (BS-WP) distribution, consisting of the difference in the response property between all pairs of cells that were recorded at different sites but within the same tetrode penetration. Each pairwise difference was tagged by the distance between the two associated recording sites, with 5th and 95th percentiles of distances given by 75 and 1,430 μm, respectively. We ordered the differences from smallest to largest recording distance and calculated the median difference in sliding windows of 1,000 (500) pairs of cells for drifting (flashed) gratings, with a window step size of 200 (100) pairs. We define a window's distance to be the median distance between sites for all the pairs it contains. To assay at what distance the BS-WP distribution was not significantly different from the full between-site (BS) distribution, we needed to control for the fact that cells would appear in many more pairs in penetrations with many recorded cells vs. penetrations with fewer cells. For this, we used a randomization test, modified slightly from the version described in *Tests for Significance of Clustering Measures*, as follows: to collect one sample of the randomized BS-WP distribution, we randomized cell response properties across all sites (across all penetrations and animals), while leaving the identity of the recording site (and its associated tetrode penetration) untouched. Thus, in each sample of the randomized BS-WP distribution, the number of pairs to which a cell contributes has the same distribution as in the observed BS-WP distribution. We generated *N*_{rep} = 2,000 such samples of the BS-WP distribution. We then partitioned each of the randomized BS-WP distributions using the same set of sliding distance windows as for the actual data. For each distance window, we calculated the median of the observed BS-WP distribution and the medians of the *N*_{rep} randomized BS-WP distributions and used Laplace's rule to obtain a two-sided *P* value, corresponding to the probability that a random median for that window would be as extreme as the window's observed median. We found the index *j* of the first sliding window for which this *P* value was >0.01. We defined the spatial extent of clustering to be the midpoint between the distance of the *j*th and *j*−1th windows, that is, midway between the distance of the first window to fail to show *P* < 0.01 and the distance of the previous window, which did show *P* < 0.01.

### Estimating Effects of Measurement Errors on Median Ratios

We consider the possible effects of errors in our measurements of response properties on our estimates of clustering. To make a ballpark estimate of the magnitude of the effect of measurement error, we assume that, for a given response property of a cell, the distribution of measurement error (that is, the probability distribution for the cell's true value of the property, centered about the measured value) can be taken to be a zero-mean Gaussian with a standard deviation given by its corresponding jackknife standard error. We further simplify by ignoring the diversity among cells in measurement error, taking all cells to have identical Gaussian distributions of measurement error with standard deviation equal to the median jackknife standard error across all cells, σ_{m} (see Table 1, “Median SE”). Suppose the true distribution of within-site signed differences in the response property is an independent Gaussian with variance σ_{ws}^{2}, so the observed within-site distribution has variance σ_{ws}^{2} + 2σ_{m}^{2}. Since the median absolute deviation of a Gaussian is about *c* = 0.675 times its standard deviation, the observed median within-site difference would be *c*. Thus the observed ratio, *X*, between the median within-site difference and the median jackknife standard error is given by *X* = *c*/σ_{m}. Rearranging, we obtain σ_{m} = *q*σ_{ws}, where *q* = 1/. Note that this approach fails if *X* < 0.95, because in that case the equation *X* = *c*/σ_{m} cannot be satisfied, since the right side is *c* > *c* = 0.95. Correspondingly, for *X* < 0.95 the equation for *q* gives an imaginary result.

Now, assuming *X* > 0.95, suppose that the true overall distribution (over all sites) of a response property is a Gaussian with variance σ_{r}^{2}, so the measured between-site distribution of differences is a Gaussian with variance 2(σ_{r}^{2} + σ_{m}^{2}). Then the observed median ratio—the ratio of the medians of the between-site and within-site distributions—should be equal to the ratio of the standard deviations:
while the true median ratio is

Putting these results together, we find that the observed median ratio is
For measurement errors to cause the median ratio to be underestimated by >10% (i.e., to have *MR*_{obs} < 0.9*MR*_{true}), the term inside the square root must be <0.81. Since this term is always >1/(1 + 2*q*^{2}), such an underestimation can only happen if 1/(1 + 2*q*^{2}) < 0.81, i.e., *X* < 2.19. For *X* < 2.19, underestimation by >10% only occurs when the true median ratio is sufficiently large:

### Simulations of Orientation Preference Differences

In discussion, as one assay of the likely “seeing distance” of our recordings, we present simulations of the expected median difference, as a function of the seeing distance of our tetrodes, between preferred orientations of cells recorded at a site and the dominant preferred orientation at the site (measured as the preferred orientation of the site's multiunits). The simulations were conducted as follows. We constructed orientation maps with an orientation period of 1.2 mm (a typical such period for cat V1; Kaschube et al. 2002) and realistic statistics according to the prescription of Kaschube et al. (2010): we let **x** represent two-dimensional cortical position and define the complex function
where *i* = , *l*_{j} is randomly chosen to be 1 or −1 for each *j*, ϕ*j* is a random phase chosen from a distribution uniform on [0, 2π), *k*_{j} is the two-dimensional vector *k*_{c}(cosθ_{j}, sinθ_{j}), with and (in radians) and *N* = 50. Expressing *z*(**x**) = *r*(**x**)*e*^{i2θ(x)} for real numbers *r* (the modulus of *z*) and θ (1/2 the argument of *z*) at each position **x**, the map of preferred orientations is given by θ(**x**). We identify θ(**x**) with the preferred orientation of the multiunits at site **x**.

The data for the distribution of orientation differences in our data set were taken from 265 recording sites across 17 of the animals (an average of 15.6 recording sites per animal). To approximately replicate this in our simulations, we generated 17 maps, each with 16 recording sites, as follows. Each orientation map was created as a 5.0 × 5.0-mm square, sampled at intervals of 10 μm. We placed 16 “recording sites” on a 4 × 4 grid, spaced 1 mm from each other and from the edge of the map. We then determined the simulated distribution of orientation differences as follows.

The preferred orientations of cells at **x** vary about θ(**x**) with a median absolute difference of 5° (Ohki et al. 2006). We modeled this by assuming that the preferred orientation of a cell at **x** was drawn from a normal distribution that was centered on θ(**x**) and had a standard deviation of σ_{θ} = 5°/0.6745 ≈ 7.4° [truncated at ±90° about θ(**x**)], which gives a median difference of 5°; we write this distribution as *N*(θ(**x**),σ_{θ}). We further assumed that the probability of recording from a cell at a distance Δ**x** from the recording site is given by a Gaussian seeing distance function parameterized by its standard deviation, σ: *G*_{σ}(Δ**x**).

For a given seeing distance σ, we built a distribution *P*_{i}^{σ} of orientation differences at the *i*th recording site, **x**_{i}, as follows. For a given cell location **x**, we calculated the distribution of differences in orientation (calculated as shortest distance around a 180° circle) between θ(**x**_{i}) and the set of orientations with distribution *N*(θ(**x**),σ_{θ}) and weighted this by *G*_{σ}(|**x**_{i} − **x**|). We summed these weighted distributions over all **x** within 500 μm of **x**_{i} (in 3 dimensions) and normalized the resulting sum of weighted distributions to have unit integral and thus again represent a probability distribution. Because of the assumption that the probability of sampling depends on three-dimensional distance via a Gaussian seeing distance function and the fact that preferred orientations depend only on horizontal (*x*,*y*) position and not on vertical (*z*) position, the same result is obtained whether summing over three dimensions or only over the two horizontal dimensions, so in practice we did the latter.

We pooled all the distributions of differences *P*_{i}^{σ} from all recordings sites *i* together to give a distribution *P*^{σ}. We calculated the median of this pooled distribution, med(*P*^{σ}), and also *p*_{out}(*P*^{σ}), the proportion of “outliers,” i.e., differences in orientation >45° (see *Local scatter of preferred orientation and “outlier” cells*). We repeated this process for a range of different values of σ, the standard deviation of the Gaussian seeing distance function, to generate the curves med(σ) and *p*_{out}(σ) for this map. We repeated this process for *n* = 17 maps, and took med(σ) and *p*_{out}(σ) to be the average of the curves obtained for all maps. We interpolated the curve med(σ) to find σ_{best}, the value of σ that gave the experimentally observed value med(σ) = 8.7° (the median difference between cells and their multiunits in response to drifting gratings that we observed in our data set; *Local scatter of preferred orientation and “outlier” cells*). We calculated *p*_{out-best} as *p*_{out}(σ_{best}). We then repeated this entire process a total of *N*_{rep} = 500 times, each time generating a new set of *n* = 17 maps, calculating med(σ) and obtaining a sample σ_{best}^{α}, and *p*_{out-best}^{α} for each iteration α. We used these repeats to get more robust estimates of σ_{av} = <σ_{best}>_{α} and *p*_{out-av} = <*p*_{out-best}>_{α}, where the brackets indicate average over the 500 iterations, as well as to calculate 95% confidence intervals for these quantities. From σ_{av} we calculated *R*_{50}, the radius in which we see 50% of the cells (i.e., the radius for which the integral of the Gaussian seeing distance function is 0.5), which for a three-dimensional Gaussian is given by *R*_{50} = 1.538 σ_{av}.

## RESULTS

We used single-tetrode recordings to record multiple neurons at single recording sites (from 2 to 8 cells/site; mean ∼2.5). After sorting spikes into clusters corresponding to cells (methods), we only accepted cells for further study if they were reasonably well isolated as assayed by an isolation distance (ID) (Schmitzer-Torbert et al. 2005) of at least 10. We also required that cells passed certain selection criteria for reproducible and tuned responses (see methods; briefly summarized below). From among these cells, we only included in the study of a given response property those cells (and corresponding pairs) for which the property was well estimated, as determined by a jackknife estimate of the standard error of the property being less than a threshold (see methods; summarized below).

We assayed clustering of each studied response property by determining whether it was, on average, significantly more similar in pairs of neurons recorded simultaneously from the same site than in pairs of neurons recorded at different sites (at different times). To do this, we formed two distributions: the within-site distribution *P*_{W} of pairwise differences in the response property among simultaneously recorded cell pairs and the between-site distribution *P*_{B} of all pairwise differences in the response property among pairs of cells recorded at different sites. To quantify the degree of clustering, we formed the median ratio: the ratio of the median of *P*_{B} to the median of *P*_{W} (DeAngelis et al. 1999). Ratios > 1 indicate that pairs recorded at the same site tended to have more similar values than pairs selected at random. We used a Monte Carlo randomization procedure to determine the significance of ratios > 1 (see methods).

We summarize the distributions of all measured cell properties in Table 1; the distributions of within-site and between-site pairwise differences for each measured property in Table 2; and the median ratios and their significance for each measured property in Table 3.

### Studies of Orientation and Direction Selectivity

For studies of orientation and direction selectivity, we used the responses of cells to *1*) *S*_{orient}, a stimulus set consisting of drifting grating stimuli at 72 directions, randomly interleaved, at 1 spatial frequency, and *2*) *S*_{flashed}, a stimulus set containing flashed gratings at 36 orientations, 10 spatial frequencies, and either 4 or 8 spatial phases, randomly interleaved. There was a significant overlap between the identities of the cells in these two sets of stimuli, as in many instances both *S*_{orient} and *S*_{flashed} were presented while recording from the same site. However, we did not try to determine correspondences between cells recorded with drifting gratings and those recorded with flashed gratings. We only analyzed cells that were significantly orientation selective and had a reproducible orientation tuning curve that was well fit by a Gaussian (see methods). For studies of orientation tuning width and preferred orientation, we required that the jackknife standard error in our estimate of the corresponding property be <5°. For DSI (which ranges from 0 to 1), we required jackknife standard error to be <0.1.

#### Local scatter of orientation tuning width.

For responses to drifting gratings, we only studied orientation tuning of cells at a single spatial frequency, 0.5 cyc/° (see methods). The orientation tuning width in cat V1 depends on the spatial frequency of the stimulus, becoming narrower with increasing spatial frequency (Hammond and Pomfrett 1990; Issa et al. 2000; Jones et al. 1987; Lampl et al. 2001; Vidyasagar and Sigüenza 1985). Thus, ideally, orientation tuning width should be assessed at the cell's preferred spatial frequency. Nonetheless we can examine the distribution of tuning widths that we observed at this single spatial frequency. For orientation tuning of cells in response to flashed gratings, we used a range of spatial frequencies, and were thus able to study orientation tuning at each cell's preferred spatial frequency.

We considered two measures of orientation tuning width. One is a “global” measure, meaning that it is influenced by the entire orientation tuning curve, not only by the shape of the region around the peak. For this, we used the standard deviation of the orientation tuning curve, treated as a distribution. We refer to this as *w*_{ORI}^{Global}. The second measure is a “local” measure, meaning that it reflects only the width of the peak of the orientation tuning curve. For this we used the standard deviation of the Gaussian in a fit of a Gaussian plus a constant to the orientation tuning curve. We refer to this as *w*_{ORI}^{Local}. As found in previous work (e.g., Alitto and Usrey 2004; Ringach et al. 2003), it is generally true that *w*_{ORI}^{Global} > *w*_{ORI}^{Local} (a scatterplot of *w*_{ORI}^{Global} vs. *w*_{ORI}^{Local} roughly fills out the triangle below the diagonal; the few points above the diagonal are close to it). The distributions of *w*_{ORI}^{Global} and *w*_{ORI}^{Local} across recorded cells are summarized in Table 1.

The pairwise distributions of the orientation tuning measures across simultaneously recorded cell pairs can be seen in Fig. 3 (Fig. 3, *A* and *C*: *w*_{ORI}^{Global}; Fig. 3, *B* and *D*: *w*_{ORI}^{Local}). The correlation between the values of *w*_{ORI}^{Global} recorded on nearby cells can be seen in the fact that the upper right and lower left corners in Fig. 3, *A* and *C*, are more dense while the upper left and lower right corners are more sparse, indicating that very widely and very narrowly tuned cells tended to be found at separate sites. In the plots of *w*_{ORI}^{Local} (Fig. 3, *B* and *D*), the dense clump of cells with smaller values shows some tilt (very weakly for flashed gratings), indicating some correlation.

To quantify the degree of clustering, we examined the within-site and between-site distributions of differences between tuning widths (either *w*_{ORI}^{Global} or *w*_{ORI}^{Local}) of pairs of cells. The cumulative within-site and between-site distributions are shown in Fig. 4, and the statistics of these distributions are given in Table 2. As can be seen, the within-site distribution is generally shifted to smaller values than the between-site distribution. To determine whether these shifts are significant and to quantify the degree of clustering, we calculated the median ratio (the ratio of the median of the between-site distribution to the median of the within-site distribution), with significance of ratios > 1 determined by a randomization procedure. Both global and local measures of orientation tuning width show significant clustering (Table 3). For *w*_{ORI}^{Global}, the median ratio is 1.55 (*P* < 10^{−4}) for drifting and 1.17 (*P* = 0.03) for flashed gratings, while for *w*_{ORI}^{Local} the ratios are 1.41 (*P* < 10^{−4}) for drifting and 1.34 (*P* = 0.002) for flashed gratings.

For *w*_{ORI}^{Local}, we found a small, but significant, positive correlation between the difference in tuning width and the difference in preferred orientation for two cells recorded at the same site [correlation coefficient (cc) = 0.10, *P* = 0.02 for drifting gratings; cc = 0.24, *P* = 10^{−4} for flashed gratings]. This means that cells with similar preferred orientations at a site were more likely to have similar tuning widths (using the local measure). For *w*_{ORI}^{Global}, this correlation was present for flashed gratings (cc = 0.24, *P* = 10^{−4}) but not for drifting gratings (cc = 0.01, *P* = 0.9).

Hetherington and Swindale (1999) found a correlation between orientation scatter (the SD of preferred orientation at sites with at least 3 cells) and orientation tuning width (the mean orientation tuning width at the site): cc = 0.7, *P* < 0.001. We observed a similar trend in our data set, which, although weaker, was still statistically significant. Using *w*_{ORI}^{Global} as the measure of orientation bandwidth, we found cc = 0.38, *P* = 10^{−4} for drifting gratings and cc = 0.45, *P* = 0.01 for flashed gratings. For *w*_{ORI}^{Local}, we found cc = 0.31, *P* = 0.008 for drifting gratings and cc = 0.34, *P* = 0.07 for flashed gratings.

#### Local scatter of preferred orientation and “outlier” cells.

The distribution of preferred orientations across simultaneously recorded cell pairs is shown in Fig. 5*A* (for drifting gratings) and Fig. 5*C* (for flashed gratings). As expected, there is obviously strong clustering, with nearby cells tending to prefer similar orientations. This is also apparent in the comparison of the cumulative within-site and between-site distributions for differences of preferred orientation, which are dramatically different (see Fig. 5, *B* and *D*, and Table 2). This is confirmed quantitatively by our clustering measures: preferred orientation had a median ratio of 3.13/3.35 (for drifting/flashed gratings, respectively), which are highly statistically significant (*P* < 10^{−4}; see Table 3).

Interestingly, in addition to the cluster of cell pairs with smaller differences, we observed a number of pairs with differences in orientation > 45° (Fig. 6, *A* and *C*; see also Table 2). We examined the distribution of the differences between an individual cell's preferred orientation and the preferred orientation of the multiunits recorded at the same site (Fig. 6, *B* and *D*). The multiunits consist of the many spikes too small to be sorted into clusters. Thus the preferred orientation of the multiunits provides a natural “reference” orientation at a site, so that differences from the preferred orientation of the site's multiunits can be interpreted as differences from the site's preferred orientation. Most cells had a preferred orientation that was <45° from the site's multiunits, but 23/566 or 4.1% of cells responding to drifting gratings (12/213 or 5.6% of cells responding to flashed gratings) had preferred orientations that differed by >45°. We refer to this small set of cells as “outliers” and the remaining cells as “typical cells.” The tail of the distribution of pairwise differences in orientation (Fig. 6, *A* and *C*) is mostly attributable to these outliers: of the pairs of cells that have differences in preferred orientation >60°, 44 of 49 (89.8%) for drifting gratings and 23 of 29 (79.3%) for flashed gratings contain at least one of these outlying cells. An example of such an outlier is *cell 2* in Fig. 1.

These outliers as assessed from responses to drifting gratings were distinguished in other ways. Most prominently, they had significantly wider orientation tuning by both the global and local measures. The median *w*_{ORI}^{Global} of the outlier cells for drifting gratings was 36% wider than that of typical cells (median for outlier cells: 29.6°, typical cells: 21.7°; *P* = 10^{−4}, *U*-test; for flashed gratings, the median was 52% wider: 38.1° vs. 25.1°, but this difference was not statistically significant: *P* = 0.21, *U*-test), while the median *w*_{ORI}^{Local} was 42% wider (23.8° vs. 16.7°, *P* = 0.03, *U*-test; for flashed gratings, the median was 59% wider: 27.0° vs. 17.0°, *P* = 0.01). At least for drifting gratings, these wider orientation tuning widths were not a reflection of greater uncertainty in estimating the width: the jackknife estimate of the standard error was not significantly different between typical and outlier cells (median SE for outlier cells: 2.0°, for typical cells: 1.5°; *P* = 0.15, Wilcoxon rank sum test). For flashed gratings, the standard error was slightly higher for outlier cells (median SE for outlier cells: 1.7°, for typical cells: 0.9°; *P* = 0.02, Wilcoxon rank sum test). The outlier cells in responses to drifting gratings also had slightly broader waveforms measured from the negative to the positive peak, with median peak-to-peak width of 0.42 ms vs. 0.31 ms for nonoutlier cells (*P* = 10^{−4}, *U*-test; 0.36 vs. 0.29 ms; *P* = 0.42 for flashed gratings). Scatterplots illustrating these various outlier properties can be found in unrefereed Supplemental Section S2.3 (see endnote).

Might these cells be poorly clustered, despite our care to consider only well-clustered units (see methods), so that their wider tuning stems from mixing the tuning of two or more cells from multiunit activity? We do not think this is the case, for two reasons. First, we calculated a measure of cluster quality, the ID (Schmitzer-Torbert et al. 2005), for all cells. This quantity measures how well separated the cell's spikes are from the rest of the recorded spikes. The outlier cells were not significantly less well isolated than the rest of the cells (median ID for outlier cells responding to drifting gratings: 29.9; for typical cells: 25.1, *P* = 0.18, *U*-test; for flashed gratings, outlier cells actually tended to be slightly better isolated, with a median ID of 37.2 vs. 26.3 for typical cells, *P* = 0.04). Second, if badly separated clusters of action potentials from multiple cells were mixed together, one might expect the predominant contribution to be from orientations near that preferred by the site's multiunits (which represent the many small action potentials that are too small to cluster). Instead, these cells had preferred orientations differing strongly from that of the multiunits.

#### Local scatter of direction selectivity.

The degree of direction selectivity of a cell was parameterized by computing the DSI (see methods) and could only be measured for responses to drifting gratings. The distribution of DSI across simultaneously recorded cell pairs is shown in Fig. 7*A*. The cumulative distributions of differences in DSI between cell pairs show that the differences tend to be slightly smaller for cells recorded at the same site than for randomly chosen pairs (Fig. 7*B*). The median ratio was 1.31 (*P* < 10^{−4}), indicating that clustering is highly significant.

We found a significant negative linear correlation between the DSI and global orientation tuning width, indicating that there was a weak tendency for more direction-tuned neurons also to be more orientation selective. The correlation coefficient between DSI and *w*_{ORI}^{Global} was −0.32 and was highly significant (*P* = 10^{−13}). This might reflect the fact that higher baseline firing will both increase *w*_{ORI}^{Global} and reduce DSI. For *w*_{ORI}^{Local}, however, which is not affected by baseline firing level, the effect was absent, (cc = −0.07, *P* = 0.14). Scatterplots of these relationships can be seen in unrefereed Supplemental Section S2.4 (see endnote).

#### Local scatter of preferred direction.

The distribution of preferred directions across simultaneously recorded cell pairs is shown in Fig. 8*A*. As for preferred orientation, there is clearly strong clustering. Two strong lines of points can be seen, representing cell pairs with similar preferred orientations with either similar preferred directions (central line) or roughly opposite preferred directions (2 peripheral lines). Between these lines can be seen a weak scattering of points, representing pairs of cells that differ in preferred orientation. The strong tendency to prefer similar or opposite directions is also apparent in the cumulative distributions (Fig. 8*B*) and in a histogram of the differences in preferred directions (shortest distance around a circle of radius 360°; Fig. 9*A*).

Because of the bimodal distribution of pairwise differences, the median ratio is not a good quantitative measure of the strength of clustering of preferred direction. We formulated an alternative quantitative measure as follows. As we did for preferred orientation, we took the preferred direction of the multiunit activity to represent the site's overall preferred direction and studied the differences between this and the preferred directions of the individual cells recorded at the site (Fig. 9*B*). We call a cell “aligned” with its site's preferred direction if the difference between its preferred direction and that of the site's multiunits is <45° (380, or 67.1% of cells), “anti-aligned” if this difference is between 135° and 180° (163, or 28.8% of cells), and “unaligned” otherwise (the 23, or 4.1%, of cells that were “outliers” of preferred orientation). We indicate the probability that a cell is aligned with the site's multiunits by *p*_{a}, the probability that it is anti-aligned by *p*_{z}, and the probability that it is unaligned by *p*_{u}; empirically, these probabilities are *p*_{a} = 0.67, *p*_{z} = 0.29, and *p*_{u} = 0.04. We computed confidence intervals for these probabilities, by calculating the Bayesian posterior probability that each *p* has a certain value, given the observed data and the assumption of uniform prior probabilities for the *p*s. We found that the 95% confidence interval for *p*_{a} is [0.63, 0.71]; for *p*_{z}, it is [0.25, 0.33]; and for *p*_{u} it is [0.02, 0.06]. The large size of *p*_{a} (0.67) compared with *p*_{z} (0.29) gives a quantitative measure of the strength of clustering of preferred direction, since these quantities would be equal if preferred direction was not clustered.

Given the proportions of aligned and anti-aligned cells across all sites, we found that they were distributed randomly across sites, without a tendency for aligned or anti-aligned cells to cluster together. To test this, we randomly shuffled the labels “aligned” and “anti-aligned” among all such cells across all sites, leaving the unaligned cells untouched, determined the proportions of mixed-type (*az*) pairs in the randomized distribution, and repeated this *N*_{rep} = 10,000 times. We found that the proportion of mixed pairs observed in the original data 40.8% (272/667) did not differ significantly from those found in the randomized distributions (mean 39.0%, *P* = 0.38 with the randomization test).

Along the lines of Hetherington and Swindale's analysis mentioned above (measuring the correlation between orientation scatter and orientation tuning width), we found a significant correlation between proportion of cells aligned at a site and mean DSI at that site: cc = 0.46, *P* = 10^{−5}, indicating that sites with greater diversity in preferred direction (i.e., a lower proportion of aligned cells) also tended to have broader direction selectivity (lower mean DSI).

### Studies of Spatial Frequency

For studies of spatial frequency tuning, we used the responses of cells to *1*) *S*_{spat}, a stimulus set consisting of drifting grating stimuli of 10 spatial frequencies (randomly interleaved), spaced approximately evenly on a logarithmic scale from 0.1 to 4.0 cyc/°, and of multiple orientations chosen according to the preferred orientations measured at the site, and *2*) the *S*_{flashed} stimulus set, which as described above contained flashed gratings at 36 orientations, 10 spatial frequencies, and either 4 or 8 spatial phases, randomly interleaved. Again, we did not try to determine correspondences between cells recorded with drifting gratings and those recorded with flashed gratings. We only analyzed cells that had a reproducible spatial frequency tuning curve that was well fit by a skewed lognormal function (*Eq. 7*; see methods). From among these cells, we again included, in studying a given property, only cells with jackknife standard error of the property less than a threshold, which was ½ octave both for spatial frequency tuning width and preferred spatial frequency.

#### Local scatter of spatial frequency tuning width.

Spatial frequency tuning width was measured in octaves as the log_{2} of the ratio of the high- and low-frequency values that gave half-maximal responses in a fitted tuning curve (see methods). Spatial frequency tuning was narrower for flashed than for drifting gratings: the median width was 1.5 octaves for drifting gratings but only 1.2 octaves for flashed gratings (Table 1).

The distribution of spatial frequency tuning widths of simultaneously recorded cell pairs shows little visually obvious clustering for the drifting gratings (Fig. 10*A*), although the tilt in the scatterplot is somewhat visible for flashed gratings (Fig. 10*C*). Differences can more clearly be seen, however, in the cumulative within-site and between-site distributions (Fig. 10, *B* and *D*). Examined quantitatively, clustering is weakly significant for responses to drifting gratings (median ratio 1.09, *P* = 0.05) and not significant for flashed gratings (median ratio 1.04, *P* = 0.3). However, when we used the F1/DC to classify cells as simple or complex (see methods; discussed further below) and restricted our analysis to pairs of simple cells, we found very strong clustering of spatial frequency tuning width in response to flashed gratings (median ratio 3.31, *P* = 0.004). This clustering was absent in pairs of complex cells (*P* = 0.8) or in simple-complex pairs (*P* = 0.4) (see Table 4). Thus it appears that, in responses to flashed gratings, simple-cell pairs show strong clustering of spatial frequency tuning width but the effect is washed out when all pairs of cells are considered.

As with orientation tuning, we found a correlation between the difference in spatial frequency tuning width and difference in preferred spatial frequency for two cells at the same site (cc = 0.14, *P* = 0.003 for drifting gratings; cc = 0.18, *P* = 0.03 for flashed gratings). This means that cells with similar preferred spatial frequencies were more likely to have similar spatial frequency tuning widths.

#### Local scatter of preferred spatial frequency.

Our recordings were within the central 10° of eccentricity. The median preferred spatial frequencies of the recorded cells were 0.50 cyc/° for drifting gratings and 0.48 cyc/° for flashed gratings, consistent with previous measures of preferred spatial frequency at these eccentricities (Movshon et al. 1978).

Preferred spatial frequency shows clustering that is much weaker than for preferred orientation or direction. We measured differences in spatial frequency in octaves, i.e., as |log_{2}(*f*_{1}/*f*_{2})|, where *f*_{1} and *f*_{2} are the two frequencies. In scatterplots (Fig. 11, *A* and *C*), preferred spatial frequencies are shown on a logarithmic scale, so that the distance from the 45° line of a data point corresponding to a simultaneously recorded cell pair represents this difference between their two preferred spatial frequencies in octaves. The clustering tests showed a median ratio of 1.35 for drifting gratings, which was highly significant (*P* < 10^{−4}). For flashed gratings, the median ratio was only 1.12 and was not statistically significant (*P* = 0.1) when all pairs of cells are considered. As with spatial frequency tuning width, however, when we restricted our analysis of flashed grating responses to pairs of simple cells we found very strong clustering, with a median ratio of 3.49 (*P* < 10^{−4}), whereas no clustering was seen in pairs of complex cells (*P* = 0.5) or simple-complex pairs (*P* = 0.9; see Table 4). Here again, it appears that despite the strong clustering of preferred spatial frequency in responses to flashed gratings for simple-cell pairs, the effect is washed out when all pairs are considered because of the lack of clustering in simple-complex and complex-cell pairs.

Analysis of clustering of preferred spatial frequency should take into account the eccentricity of cells at a particular recording site, since preferred spatial frequency depends on eccentricity (e.g., Movshon et al. 1978). However, we do not have precise eccentricity data for our recordings. Therefore, the above approach will necessarily overestimate the degree of local clustering, as we will be comparing the local variability of preferred spatial frequency at a single recording site at one eccentricity to the variability across all recording sites over a range of eccentricities.

As a partial corrective to this, we also calculated the between-site distribution using only pairs of cells across sites recorded from different depths of the same tetrode penetration. Our penetrations were down the medial bank of V1 and hence were moving largely horizontally across cortex, so eccentricities changed across a penetration, but we expect that on average eccentricities varied less within a penetration than between penetrations. On the other hand, since preferred spatial frequency does spatially cluster in cortex on local scales (e.g., Issa et al. 2000), cells in the same penetration will have less diversity in preferred spatial frequency than cells sampled from different penetrations at the same eccentricities (see *Spatial Extent of Clustering Measures*), in which case this approach will underestimate the median ratios. The median ratios obtained by comparing the within-site distribution with the between-site, within-penetration distribution are, indeed, lower than those we obtained by comparing with the full between-site distribution (1.22 for drifting gratings, 1.12 for flashed gratings) but are still highly significant for drifting gratings (*P* = 0.0007 for drifting gratings, *P* = 0.09 for flashed gratings, using all pairs of cells). (These significance values were obtained with a similar randomization scheme as before, except that we only shuffled cells among sites in the same penetration).

As we did for orientation and direction, we applied Hetherington and Swindale's analysis (measuring the correlation between scatter and tuning width) to spatial frequency and tested whether sites with greater diversity in preferred spatial frequency [SD of (log) preferred spatial frequency] also had broader tuning (mean spatial frequency tuning width). In contrast to orientation and direction, however, we found no correlation, either for drifting gratings (cc = −0.02, *P* = 0.9) or flashed gratings (cc = 0.24, *P* = 0.33).

### Correlation Coefficients Between Response Curves

As a separate method to assess clustering of tuning properties, we examined the cc between tuning curves of cells at the same site. This is a measurement that combines correlation in preferred stimulus and tuning width. The mean cc between pairs of orientation tuning curves was 0.32 and 0.40 for drifting and flashed gratings, respectively. This is in comparison to the between-site distribution, for which the mean cc was 0. The mean cc between pairs of spatial frequency tuning curves at a site was 0.50 for both drifting and flashed gratings, compared with the mean of the between-site distribution, which was 0.39 (we did not expect this to be zero, since the distribution of preferred spatial frequencies is not uniform). For flashed gratings, we measured responses for all combinations of orientations and spatial frequencies (36 orientations × 10 spatial frequencies = 360 stimuli, averaged over spatial phase). The mean cc between these 360-element vectors between pairs of cells at the same site was 0.33, compared with 0.09 for cells between sites (this was also nonzero because of the nonuniform distribution of preferred spatial frequencies).

### Simple vs. Complex Cells

We have already noted the distinction between simple and complex cells in our studies of spatial frequency tuning in response to flashed gratings. Here we more generally examine the role of this distinction in our results. To classify cells as simple or complex, we used the ratio of the first harmonic (F1) to the mean (DC) of the phase tuning curve of the cell (for a drifting grating, this is just the average temporal response over a stimulus cycle), measured at its preferred orientation and spatial frequency. Cells were then classified as simple (F1/DC > 1) or complex (F1/DC < 1) (Skottun et al. 1991; but see Kagan et al. 2002; Martinez et al. 2005; Mechler and Ringach 2002; Priebe et al. 2004 for arguments against the adequacy of this classification). We only accepted cells for study of the F1/DC and simple/complex classification if the jackknife estimate of standard error of the F1/DC was <0.25. As expected, we found a bimodal distribution of the F1/DC, with one mode located around 0.4 and the other around 1.6 in responses to both drifting gratings and flashed gratings (histograms of F1/DC and example phase tuning curves can be found in unrefereed Supplemental Section S2.5; see endnote). For studies of orientation and direction tuning in response to drifting gratings we studied F1/DC at a fixed spatial frequency (0.5 cyc/°) rather than the preferred spatial frequency, but we estimate that this should not greatly distort results (see methods). Flashed gratings are not typically used to classify cells as simple or complex. However, the F1/DC of the phase tuning curves in response to flashed gratings has been shown to be strongly correlated with the ratio obtained in response to drifting gratings, at least when both are characterized at the drifting grating preferred orientation and spatial frequency (Nishimoto et al. 2005). Based on this, we applied our analysis of simple/complex cells to responses to flashed gratings as well.

As expected from the different proportions of simple vs. complex cells in different layers of V1 (Gilbert 1977; Hubel and Wiesel 1962), there is a tendency for simple or complex cells to cluster together at specific recording sites. For cells responding to drifting gratings in studies of spatial frequency tuning (for which F1/DC was characterized at the preferred spatial frequency), *x* = 50% of cells were simple and *y* = 50% were complex. Thus we would expect 2*xy* = 50% of pairs to be of mixed simple-complex (SC) type if simple and complex cells did not tend to co-occur at the same site, whereas only 29% were. To test the statistical significance of this, we randomly shuffled cells across the recording sites *N*_{rep} = 10,000 times, without regard to their identity as simple or complex, and never observed 29% or fewer SC pairs. Thus the clustering of simple and complex types is highly significant (*P* < 10^{−4}). Similarly, for cells responding to flashed gratings, 39% were simple and 61% were complex, so we would expect 48% of pairs to be SC without clustering, but only 28% were. Again, under randomization we never observed 28% or fewer to be SC (*P* < 10^{−4}). This tendency for simple/complex cells to cluster also manifests in clustering of the F1/DC: the distributions of differences in F1/DC (for the cells studied for spatial frequency) had median ratios of 1.48 and 1.84 for drifting and flashed gratings, respectively, which were both highly significant (*P* < 10^{−4} and *P* = 0.003).

We examined the relationship between simple/complex cell classification and clustering of response properties (preferred orientation/spatial frequency, and orientation/spatial frequency tuning widths and DSI). For each response property, we separated the original within-site and between-site distributions into subgroups that contained SS (simple-simple), CC (complex-complex), or SC (simple-complex) pairs and then computed the median ratios separately for SS pairs, CC pairs, and SC pairs. We then used a randomization test to calculate the significance of the differences in ratio between any two of these pair types (methods), under the null hypothesis that all three pair types clustered to the same degree. In most cases, ratio differences between pair types were not significant at the *P* = 0.05 level. The notable exception was in studies of spatial frequency in response to flashed gratings, for which, in accord with the presence of significant clustering for SS but not SC or CC as described above, we found significant differences between SS and SC and between SS and CC but not between SC and CC (see Table 4).

There was one other set of significant differences between clustering of pair types, for *w*_{ORI}^{Global} for drifting gratings, but its meaning was equivocal: SC pairs had significantly lower median ratios than SS or CC pairs, but this might have been an artifact of simple cells having narrower median *w*_{ORI}^{Global} than complex cells; this median difference might have dominated both between-site and within-site pairwise differences, leading the two distributions to be more similar for SC than for SS or CC. There were no other significant differences between the three pair types. We also confirmed earlier reports (e.g., Henry et al. 1974; Rose and Blakemore 1974; Wörgötter et al. 1991; but see Gizzi et al. 1990) that simple cells are more sharply tuned than complex cells, although the effect was present only for drifting grating responses for orientation tuning width (and DSI) and only for flashed grating responses for spatial frequency tuning width. Because none of these results gave clear or compelling new conclusions, we leave details of these findings on simple/complex differences to unrefereed Supplemental Section S2.6 (see endnote).

### Spatial Extent of Clustering Measures

We asked over what distance the measured clustering extends. To address this, we considered sites recorded within the same tetrode penetration. We binned the set of between-site, within-penetration cell pairs according to the distance between the two sites at which the cells were recorded. For each response property, we defined the spatial extent of clustering to be the intersite distance at which the distribution of pairwise differences in the response properties became statistically indistinguishable (*P* > 0.01) from the full between-site distribution of such differences (which includes between-penetration as well as within-penetration pairs) (see methods). The spatial extents of clustering of all measured cell response properties are listed in Table 5, and the dependencies on distance are shown in Fig. 12 and Fig 13.

Figure 12*E*, for example, shows the median difference in preferred orientation for drifting gratings as a function of distance between recording sites, using a sliding window of 1,000 pairs (with pairs ordered by their between-recording-sites distance). For small distances between sites (∼100 μm) the median difference is similar to that of the within-site distribution and significantly (*P* < 0.01) different from the randomized control. For distances > 465 μm, however, the median difference ceases to be significantly different from the control, and thus the spatial extent is defined as 465 μm. Interestingly, however, the median differences again become significant after ∼800 μm. This is consistent with our tetrode penetrations being mostly horizontal, parallel to the cortical surface (as expected, since most of our penetrations were down the medial bank of V1): since preferred orientation in cats is periodic across the cortical surface with an average period of ∼1.2 mm (Kaschube et al. 2002), two cells that are ∼1 mm apart horizontally will be more likely to have a similar preferred orientation than two cells picked at random. Thus our measure of the spatial extent of clustering largely reflects the dependence on horizontal distance; our methods did not allow an assay of dependence on vertical distance or on layer.

For most response properties, the spatial extent of clustering was typically ∼200 μm. It was considerably greater only for drifting grating responses for preferred orientation, orientation tuning width, and F1/DC (the large spatial extent for F1/DC likely reflects the well-known tendency of simple and complex cells to predominate in different layers, e.g., Martinez et al., 2005), while no intersite clustering could be detected at any distance for spatial frequency tuning width or, for drifting gratings, for preferred spatial frequency. Since spatial frequency measures for flashed gratings were significantly clustered only among simple cells, we would have liked to analyze the distance dependence of SS pairs but lacked sufficient data. It is possible that the lesser amount of data for flashed gratings, which led us to use sliding distance windows containing half as many points for flashed as for drifting gratings, might have played some role in our failing to detect significance at further distances for flashed gratings for some of the orientation measures or F1/DC.

The observed spatial extents would not be greatly changed by taking into account the fact that tetrodes sample cells over a finite distance. In discussion we estimate that, if we imagine that the probability of sampling cells as a function of distance from the tetrode is described by a three-dimensional Gaussian distribution, then it would have a standard deviation of ∼σ_{seeing} = 85 μm. If we also imagine that the actual strength of clustering decreases as a Gaussian function of distance with standard deviation σ_{cluster}, then in terms of the observed spatial extent *s*, and assuming *s* > σ_{seeing}, we would have σ_{cluster} = *s*. A measured spatial extent of 200 μm would be corrected to 181 μm under this transformation.

### Potential Confounds

For some of the tuning properties we measured, clustering was very weak. A worry is that this weak clustering might be artifactual. We considered three potential sources of contamination of clustering measures: bad cell isolation, noise correlations, and variations across animals. We also considered a confound that might have made clustering appear artifactually weak: errors in estimating response properties.

#### Cell isolation.

Similarity in response properties of cells at a single recording site might be artifactually induced by poor spike sorting, in one of at least two ways: *1*) tuning curves of two cells might appear more similar than they really are, if spikes from one cell are incorporated in the tuning curve of the other cell, or *2*) a cell might be compared with itself, if the spike distribution from a cell is bimodal in the clustering space and is mistakenly identified as two separate cells.

We only included in our study cells that were reasonably well isolated, as judged by the isolation distance (ID) (Schmitzer-Torbert et al. 2005), a measure of cluster isolation for each cell. Clusters with higher IDs tend to have fewer misclassified spikes, but there is no “correct” criterion for the minimum acceptable ID: the appropriate threshold will depend on the level of contamination acceptable for each specific study. We only used cells that had an ID of at least 10.

Nonetheless, if similarity in a pair's response properties correlated to low ID or to similarity in the two cells' spike waveforms, this could be an indication that poor spike sorting contributed to apparent clustering. We examined this possibility in three different ways, all with negative results. *1*) We calculated the Spearman (rank) correlation between the measured pairwise differences in a tuning property (e.g., differences in orientation tuning width) and the minimum ID of a cell pair, across all cell pairs. We did this for differences in all tuning properties and checked for correlation coefficients that were significantly greater than zero. No significant trends were found for any tuning properties (*P* > 0.3 for all rank correlation coefficients). *2*) We used a Wilcoxon rank sum test to test whether the median pairwise differences in tuning properties differed between the third of cell pairs with the smallest minimum IDs and the third of cell pairs with the largest minimum IDs. Again, no significant trends were found (*P* > 0.1 in all cases). *3*) Instead of using the minimum ID, which represents the cluster quality of the less well isolated in a pair of cells, we also more directly examined the separation between two cells in the clustering space, either using the simple Euclidean distance between the cluster centers or by examining a measure of the overlap of the two clusters when each is modeled as a multidimensional Gaussian (Gaussian overlap, defined in *Spike detection and sorting*). As with the minimum ID, we tested whether cells that were closer together in the clustering space (and were thus more likely to have spikes misclassified from one to the other) tended to have smaller differences in tuning properties. We used both Spearman's rank correlation coefficient and the Wilcoxon rank sum tests described above and found no significant trends (*P* > 0.1). In sum, we do not think that bad cell isolation contributed to the clustering reported here.

#### Noise correlations.

Another potential confound in measuring correlation is to mistakenly identify noise correlations (coordinated responses due to connections between cells or shared inputs) as signal correlations (correlations due to similar tuning properties). That is, because our response tuning measures are based on finite samples, it is possible that correlated noise fluctuations might contribute significantly to the measured tuning curves of a cell pair, particularly to the portions of the tuning curves corresponding to weak responses, causing the tuning curves to appear more similar than they really are. Since noise correlations are less likely to affect the peaks of the tuning curves, where responses are stronger, and since many of our measures depend only on the regions around the peaks, we did not think this scenario likely on its face.

Nonetheless, we tested this possibility as follows. We repeated our analysis but this time divided the experiment into odd- and even-numbered trials (as we did for the reproducibility tests and for the measures of variability) and only compared tuning curves that were obtained during different halves of the experiment. For example, for a particular measure of tuning width for cell *i*, we call the odd-trial-averaged and even-trial-averaged tuning widths *w*_{odd-i} and *w*_{even-i}. Then, to calculate the difference in tuning width between cells *i* and *j*, we calculate both |*w*_{odd-i} − *w*_{even-j}| and |*w*_{even-i} − *w*_{odd-j}| and place both of these values into the within-site (or between-site) distribution. (Our original distributions simply contained the value |*w*_{i} − *w*_{j}| for each pair of cells *i* and *j*). These response differences, calculated from different sets of trials for the two cells, cannot contain any contributions from noise correlations, which only can affect two cells on the same trial. We compared the median ratios we obtained with these odd/even-trial-averaged response properties to our original all-trial-averaged response properties. There were minor differences (some slight increases in the median ratio, some slight decreases) but no changes in significance levels, defined to mean whether a given clustering measure met zero, one, or both of the two (somewhat ad hoc) conditions *P* ≤ 0.05 and *P* ≤ 0.01. Thus noise correlations do not seem to have introduced spurious correlations that significantly affected our results.

#### Variations across animals.

A potential problem with using the full between-site distribution when calculating the median ratio is that it includes pairs of cells from different animals. Clustering of a response property means that there is less variability in the property at individual sites in an animal than across that animal's V1. To the extent to which there is between-animal variability in the distributions of these properties, comparing the within-site distribution of differences in the property to the full between-site distribution (including differences between pairs of cells recorded in different animals), as we have done, would tend to exaggerate the degree of clustering. We did not have sufficient data per animal to directly assay for between-animal differences in response property distributions, e.g., for studies of orientation tuning in response to flashed gratings we had, on average, 29 cells/animal, 12.6 sites/animal, and 3 penetrations/animal (see *Cell selection criteria and the measure of preferred orientation*), with cell properties generally correlated within sites and between nearby sites in a penetration as we have seen above.

Instead, as an additional control, we sought to redo our calculations using the between-site, within-animal (BS-WA) distribution as the control instead of the full between-site (BS) distribution. However, a sizable proportion of the BS-WA distribution consists of cell pairs from the same penetration, and a sizable proportion of those pairs comes from recording sites close enough together to show more similar cell response properties than expected by chance (see *Spatial Extent of Clustering Measures*). Including these pairs with smaller differences in the control distribution would tend to underestimate the median ratio. Thus we excluded within-penetration (BS-WP) pairs and considered only the between-penetration, within-animal (BS-BP-WA) pairs. We recalculated the median ratios using this (BS-BP-WA) distribution as a control and checked whether the median ratios were significantly different from those obtained with the full BS distribution. We found that, with one exception, median ratios were only slightly affected, with statistically significant ratios changing by <7%, including both increases and decreases, standard error ranges of the ratios overlapping with those found previously, and little or no change in significance level (unrefereed Supplemental Section S2.8 and Supplemental Table S6; see endnote).

The one exception was global orientation width in response to flashed gratings, for which clustering disappeared (median ratio 1.0, SE range 0.93–1.09, *P* = 0.5) when the BS-BP-WA distribution was used as a control (clustering was weak but present with the full BS distribution: median ratio of 1.17, SE range 1.10–1.30, *P* = 0.03). However, the gradual decrease in clustering with increasing inter-recording-site distance of global orientation width for flashed gratings (Fig. 12*B*) strongly suggests true clustering of this property. The BS-BP-WA distributions have many fewer pairs than the full BS distributions (for flashed grating global orientation width: BS-BP-WA distribution, 1,564 pairs, with a mean of 195 pairs/animal; BS distribution, 22,144 pairs). Thus it seems likely that the reduction in median ratio in this case may reflect a lack of adequate sampling of this response property within individual animals, which would artifactually reduce the median ratio, rather than being due to the removal of interanimal variations in this response property, which would indicate a truly smaller ratio.

#### Errors in estimating response properties.

Errors in our measurements of response properties could artifactually reduce the strength of clustering that we measure, by adding apparent within-site variability that was not present in reality. We estimated the precision with which we can measure response properties by the jackknife estimate of standard error in the measured property. For most response properties, the median within-site difference (Table 2) was considerably larger than the median jackknife standard error for measurements of that property (Table 1), indicating that the within-site variability we measured was well above any limit placed by measurement errors and hence that there was little underestimation of the median ratio. The ratio, *X*, of median within-site difference to median jackknife error was always >2.2, and it was >4 for all but a few properties: F1/DC for orientation stimuli (ratios 2.9/2.2 for drifting/flashed gratings), F1/DC for flashed spatial frequency stimuli (ratio 2.5), and spatial frequency tuning width (ratios 3.8/2.2). (However, the ratio was considerably smaller for spatial frequency tuning width and preferred spatial frequency for flashed gratings for simple-simple pairs, as further discussed below.)

To make a ballpark calculation of the quantitative effect of measurement errors, we made simplifying assumptions that the true within-site distribution, the true between-site distribution, and the distribution of measurement errors are all Gaussians, with the standard deviation of the measurement error distribution given by the median of the jackknife standard errors that we observed over cells (see *Estimating Effects of Measurement Errors on Median Ratios*). We then examined the effect of measurement error on our estimate of the median ratio, which is our measure of the strength of clustering. The result is that, for measurement error to have caused our observed median ratio to underestimate the true median ratio (call it *MR*_{true}) by >10%, it must be the case that *X* < 2.19 and *MR*_{true} > 1/. These conditions exclude all of our data from Tables 1–3. If *X* > 4, as is the case for most of our measurements, the median ratio will be underestimated by at most 3%. This calculation, though based on simplified assumptions, supports the idea that measurement error did not greatly reduce our estimates of clustering for any of the quantities in Table 3.

However, the effect of measurement error was likely to have been more significant for the strong median ratios we observed for spatial frequency tuning width and preferred spatial frequency of simple-simple pairs responding to flashed gratings (Table 4). For these pairs and response properties, the median jackknife standard error (tuning width: 0.21; preferred: 0.16) was comparable to the median within-site difference (tuning width: 0.11; preferred: 0.21). For preferred spatial frequency, with the observed ratio *X* = 1.3 and an observed median ratio of 3.49, the above approach yields an estimate that the true median ratio was 4.97, i.e., the observed median ratio underestimated the true value by 30%. For spatial frequency tuning width, for which *X* = 0.52, the above approach fails (see methods). However, the fact that the median standard error was greater than the median within-site difference for tuning width suggests that observed within-site variability was roughly as small as possible given measurement error, i.e., clustering was very strong. Thus preferred spatial frequency and spatial frequency tuning width are likely to be considerably more strongly clustered for simple-simple pairs responding to flashed gratings than we have estimated, even though our estimates already showed strong clustering.

## DISCUSSION

We have found that neurons in cat V1 show modest but statistically significant clustering by degree of tuning, as summarized graphically in Fig. 14 (all quantitative measures of clustering summarized in Table 3). To measure clustering, we used the median ratio, the ratio of the median difference between the value of the property for two cells from different sites to that between two cells recorded at the same site. In line with previous studies, we found a strong clustering of orientation tuning width: for the local measure *w*_{ORI}^{Local}, for example, we found a median ratio of 1.41 (*P* < 10^{−4}), corresponding to a 29% decrease in the median difference between pairs from the same site relative to pairs from different sites. For flashed gratings, the median ratio was smaller (1.34, corresponding to a 15% decrease in the median difference) but still statistically significant (*P* = 0.002). Median ratios for global orientation width, *w*_{ORI}^{Global}, were similar. In contrast to previous studies, we also found highly significant clustering of direction selectivity (DSI), with a median ratio of 1.31 (*P* < 10^{−4}). We found weak clustering of spatial frequency tuning width in response to drifting gratings, with a median ratio of 1.09 (*P* = 0.05). For flashed gratings, clustering of spatial frequency tuning width was strong when considering only simple-cell pairs (median ratio 3.31, *P* = 0.004; as discussed above, this is likely to be an underestimate of the true strength of clustering in this case, due to errors in estimating the width), but no clustering was seen in simple-complex or complex-complex pairs. These measures of clustering of tuning widths should be compared to the generally much stronger clustering of preferred orientation, which for flashed gratings had a median ratio of 3.35 (*P* < 10^{−4}), corresponding to a 68% reduction in the median difference at a site, and of preferred direction, for which we found that ∼61% of cells at a site are aligned with one preferred direction and 26% with its opposite (the median ratio is not a good measure for preferred direction because of the typical presence of anti-aligned as well as aligned cells). The clustering of preferred spatial frequency for drifting gratings was small compared with that of preferred orientation and comparable to the other measures of tuning widths, with a median ratio of 1.35 (*P* < 10^{−4}). For flashed gratings preferred spatial frequency, like spatial frequency tuning width, showed strong clustering for simple-simple pairs (median ratio 3.49, *P* < 10^{−4}; in this case also, this is likely to be an underestimate of the true strength of clustering) but no clustering for simple-complex or complex-complex pairs. Clustering had a horizontal spatial extent of ∼200 μm for the majority of properties, but there was a considerably larger spatial extent of clustering in responses to drifting gratings for preferred orientation, orientation tuning widths, and F1/DC, while there was no detectable spatial extent of clustering for spatial frequency tuning width or, for drifting gratings, preferred spatial frequency.

As noted in the introduction, previous studies failed to find significant clustering of direction selectivity (DeAngelis et al. 1999; Martin and Schroder 2013) or spatial frequency tuning width in responses to drifting gratings (Martin and Schroder 2013) but only had access to small numbers of pairs and, in the case of DeAngelis et al. (1999), used a different methodology for assaying DSI. It seems likely that these small numbers lacked the statistical power to detect the relatively modest clustering we found for these measures.

For many response properties we observed differences between responses to drifting vs. flashed gratings in tuning strengths or in the strength or spatial extent of clustering. Particularly notable were that *1*) for flashed but not drifting gratings, significant clusterings of spatial frequency tuning width and preferred spatial frequency were seen only for simple-simple pairs, and these simple-simple clusterings for flashed gratings were much stronger than the clusterings observed for these properties for drifting gratings (Table 4); and *2*) the spatial extents of clustering for preferred orientation, orientation tuning widths, and F1/DC were restricted to ∼200 μm for flashed gratings but were 1.5–4 times larger for drifting gratings (Table 5). Some of these differences might reflect differences between steady-state, adapted cellular and network responses, which are measured in response to drifting gratings, and transient, unadapted responses, which are measured in response to flashed gratings.

The local variability we found in preferred stimuli, which serves as a form of calibration of our measurements, is similar to that found in previous studies (Albus 1975; Berman et al. 1987; DeAngelis et al. 1999; Hetherington and Swindale 1999; Lee et al. 1977; Maldonado and Gray 1996; Martin and Schroder 2013; Tolhurst et al. 1981; Tolhurst and Thompson 1982), except that we found more local variability of preferred spatial frequency to drifting gratings than previous studies. A full discussion of this comparison to previous studies is available in unrefereed Supplemental Section S3 (see endnote).

### Origins of Local Response Variability

While we found significant clustering in orientation and spatial frequency tuning width and in direction selectivity, there is still a large degree of diversity at a site, as is apparent from scatterplots of each property (Figs. 3, 7, 10, and 11). The local diversity we observe in response properties could have at least two causes. First, it may reflect genuine disorder in the spatial arrangement of these properties. Second, they might be ordered on a much finer scale than the scale over which the tetrode samples. Our results show that measures of selectivity (orientation and spatial frequency tuning width and direction selectivity) must be much more weakly clustered and/or clustered on a much finer spatial scale than preferred orientation or direction.

#### Tetrode seeing distance.

Our findings of weak clustering for many response properties indicate diversity in these response properties over the tetrode sampling distance but do not address whether stronger clustering may exist on finer scales. Thus it is important to estimate the tetrode sampling distance.

Mechler et al. (2011) modeled the spiking neuron as a current dipole, as suggested theoretically (Mechler and Victor 2012). They concluded from measurements in cat and macaque V1 that the recorded spike corresponds to a dipole that extends some tens of micrometers along the primary dendrite that most closely points to the electrode tip. They calculated a seeing distance (*R*_{50}, defined as the radial distance from the tetrode containing half of isolated source dipoles) that they noted scales approximately linearly with the separation Δs between the tetrode wires. They found *R*_{50} = 100 μm for tetrodes with Δs = 45 μm and *R*_{50} = 80 μm for Δs = 35 μm. In most of our experiments we used smaller tetrodes, with Δs ranging from 8 to 25 μm [separations could be larger if tips were splayed (Chelaru and Jog 2005; Jog et al. 2002), but we believe our tips were generally intact; see methods]. Extrapolating from their data, their approach would estimate the *R*_{50} of most of our tetrodes to be in the range 25–60 μm. However, this is the distance over which the nearest primary dendrite of 50% of recorded cells would be found. Given that dendrites of cells in cat V1 typically extend horizontally at least 100 μm (e.g., Kelly and van Essen 1974), this suggests that the seeing distance in terms of cell bodies of recorded cells would be considerably larger.

We estimate the seeing distance in terms of cell bodies using our observed distribution of differences in orientation preferences. In calcium-imaging studies in upper layers of area 18 in 4- to 5-wk-old kittens, Ohki et al. (2006) found that the median difference of a cell's preferred orientation from the preferred orientation expected at the cell's site was 9° within 65 μm of pinwheel centers and 5° further from pinwheels. For orientation maps with a period of 1,200 μm, a typical period for cat V1 (Kaschube et al. 2002), and a density of π pinwheels per (1,200 μm)^{2} (Kaschube et al. 2010), regions within 65 μm of a pinwheel occupy <3% of the cortical surface, so we neglect the greater diversity near pinwheels. We equate the preferred orientation of the multiunits at a site with the preferred orientation at the site and assume that the local preferred orientations at the site have a Gaussian distribution about the site's preferred with a median absolute difference of 5°. We then asked (see methods) over how wide an area of statistically realistic orientation maps, averaging over maps and map locations, we would have to sample to reproduce our observation that the median difference of a cell's preferred orientation from that of the multiunits at a site is ∼8.7° (for drifting gratings). We assumed that cells are sampled as a Gaussian function of distance with a standard deviation σ. As shown in Fig. 15, this produced the estimate σ = 85 μm (95% confidence intervals 77–91 μm), which for three-dimensional sampling corresponds to an *R*_{50} of 130 μm (95% confidence intervals 118–140 μm). With this seeing distance, the expected proportion of “outliers” in the simulations is 4.8% (see Fig. 15*D*), close to our observed 4.1%, which further corroborates this estimate of seeing distance. This result is similar to that of Hetherington and Swindale (1999). They ignored local diversity of preferred orientations and found that their measured distribution of pairwise differences in preferred orientations could be reasonably accounted for by sampling an orientation map in two dimensions, with sampling probability given by a uniform distribution over a circle of radius 30 μm convolved with a Gaussian distribution with σ = 60 μm (similar to the σ = 85 μm of the Gaussian seeing distance function in our simulations).

This estimate for seeing distance appears consistent with the distribution of differences in preferred orientation from Albus (1975). He measured the distributions of differences in preferred orientation between pairs of cells recorded by electrodes at different tangential distances from one another. His Figure 9b, showing the distribution of pairwise differences between cells recorded at sites 10–100 μm apart, closely resembles our Fig. 6*A* (direct comparison of figures is available in unrefereed Supplemental Fig. S9; see endnote), suggesting that our cell pairs and his had a similar range of separations. The range of intercell separations in Albus's recordings would have exceeded the 100-μm maximal distance between recording electrodes, since the seeing distances of his electrodes must also be taken into account. This produces a rough estimate of seeing distance that seems in line with our estimate from our simulations of orientation maps.

In summary, based on our simulations and comparison with some previous studies, we estimate the seeing distance (*R*_{50}, the distance over which ½ of cell bodies of recorded cells would be found) of our tetrodes to be ∼130 μm. Simultaneous tetrode recording and calcium imaging could provide precise answers as to tetrode sampling, at least in upper layers where calcium imaging is readily accomplished. While calcium imaging can directly give more precise answers to clustering within upper layers, tetrodes calibrated by calcium imaging in upper layers can assay clustering across the cortical layers.

#### Preferred orientation outliers.

As mentioned above, in our studies of orientation tuning we found that a small fraction of cells (4.1% responding to drifting gratings) had preferred orientations that differed from that of the site's multiunits by >45°. These cells tended to have broader orientation tuning curves and broader waveforms (measured from the negative to the positive peak). These “outliers” were responsible for the majority of pairs of cells with larger (>60°) differences in preferred orientation, which have also been seen in previous studies (Lee et al. 1977; Maldonado and Gray 1996). One possibility is that these cells may represent a different class of cell, which tends to have broader orientation tuning and slightly wider waveforms. This is consistent with subsequent reports of complex inhibitory neurons with poor orientation tuning (Hirsch et al. 2003; Nowak et al. 2008; see also Sohya et al. 2007). Another possibility is that these cells are simply the occasional cells that are sampled from sites of very different preferred orientation than the recording site, either when recording near pinwheels or from the small number of cells recorded at the furthest distances from our tetrodes. Indeed, when we sampled from simulated orientation maps, using a Gaussian seeing distance function that matched our observed median orientation difference of recorded cells from multiunits, we observed 4.8% outlier cells (95% confidence interval: 4.1–5.6%), which is close to the 4.1% that we observed (Fig. 15*F*). Cells at larger distances would be more likely to be sampled if they have more widespread dendritic arbors, according to the theory of Mechler and Victor (2012), and wider dendritic arbors might provide a mechanism underlying the greater-than-average orientation tuning width seen in these outlying cells.

#### Variations across layers.

Response properties that differ between layers would appear clustered by our analysis, since the within-site distribution (which would mostly sample single-layer cell pairs) would show less variability than the between-site distribution (which would sample many pairs in which the 2 cells came from different layers). Since we did not do postrecording histological analysis, we do not know in which layers each set of cells was recorded and cannot directly assay laminar differences in response properties. However, if clustering were due to laminar differences, it would be expected to persist over long horizontal distances. Our results showed that for most properties the spatial extent of clustering, which we argued is very largely a horizontal (within layer) distance, was <400 μm (Table 5, Figs. 12 and 13). The exceptions were, for responses to drifting gratings, preferred orientation, which does not vary between layers in cats (e.g., Hubel and Wiesel 1962) but which showed spatial extent of clustering consistent with the known periodic arrangement of orientation maps (e.g., Kaschube et al. 2002), and F1/DC, which is known to cluster by layers [e.g., Martinez et al. (2005) showed that simple cells predominate in cat layer 4, and Ringach et al. (2002) showed a corresponding result in monkeys for F1/DC for drifting gratings]. Thus we can conclude that, with the exception of the F1/DC, laminar differences in response properties were not likely to have contributed substantially to the clustering we have seen.

#### Effect of stimulus ensemble.

It is possible that some of the diversity we observed in orientation tuning width might actually represent underlying diversity in spatial frequency tuning and some of the observed diversity in direction selectivity might represent diversity in temporal tuning. Orientation tuning width of neurons in cat V1, as assayed by local measures, has been shown to decrease with increasing stimulus spatial frequency (Hammond and Pomfrett 1990; Issa et al. 2000; Jones et al. 1987; Lampl et al. 2001; Vidyasagar and Sigüenza 1985). At least in ferret V1 it is independent of temporal frequency (Moore et al. 2005). DSI in both cat and ferret V1 has been shown to decrease with deviation from the preferred temporal frequency (Moore et al. 2005; Saul and Humphrey 1992). For drifting gratings, we studied both orientation and direction selectivity at a single spatial and temporal frequency, inducing diversity from the varying preferred spatial and temporal frequencies of recorded cells that might contribute to diversity of orientation tuning width and DSI. However, we found similar diversity for local orientation tuning width when we used flashed gratings, for which cells were studied at their preferred spatial frequency, suggesting that the use of a single spatial frequency for drifting gratings may not have greatly impacted our results, at least for orientation tuning width. With regard to spatial frequency tuning, most studies have found that it appears to be largely independent of stimulus temporal frequency in cat V1 (Foster et al. 1985; Holub and Morton-Gibson 1981; McLean and Palmer 1994; but see Hammond and Pomfrett 1990), although the same is not true in monkey V1 (Priebe et al. 2006). Thus the fact that we studied spatial frequency tuning using drifting gratings at a single temporal frequency seems unlikely to have affected our results.

### Conclusion: Implications for Cortical Circuitry

Several studies have demonstrated that nearby neurons carry surprisingly independent information, despite their similarities in preferred stimulus features (Gawne et al. 1996; Martin and Schroder 2013; Montani et al. 2007; Reich et al. 2001; Vinje and Gallant 2000; Yen et al. 2007). The diversity we have found among nearby neurons in their degree of stimulus selectivity may be one important mechanism underlying this finding. Variability in tuning width or selectivity might be explained at least in part by variability in the overall strength of inhibition received by cells relative to the excitation they receive (e.g., McLaughlin et al. 2000; Troyer et al. 1998). Variability in stimulus preferences, on the other hand, would seem to require that nearby cells receive quite different patterns of excitatory input, even while each cell receives a specific pattern of input. In particular, local variations in preferred spatial frequency and more generally in space-time receptive fields (DeAngelis et al. 1999) raise important problems. Consider, for example, two nearby simple cells with similar preferred orientations receiving input from the same portion of lateral geniculate nucleus (LGN) but showing an alternation of ON and OFF subregions (Reid and Alonso 1995) corresponding to two quite different preferred spatial frequencies (for simple cells, preferred spatial frequency is linearly related to, though systematically higher than, the spatial frequency corresponding to the subregion alternation; Jones and Palmer 1987, Fig. 7). If receptive field development is driven by maximizing correlations among the inputs a cell receives, as suggested by developmental modeling (reviewed in Miller 1996), then one cell would have detected correlations oscillating with one frequency across the LGN input cells, while the other cell detected a quite different correlation structure, even while both cells “locked onto” the same orientation. Our finding that, in response to flashed gratings, simple cells show a high degree of clustering both of preferred spatial frequency and of orientation tuning width suggests that, locally, simple cells may all “lock on” to a common correlation structure, but the finding that this is not true for responses to drifting gratings casts doubt on that conclusion.

Assuming that nearby simple cells may often have very different preferred spatial frequencies, how can two cells detect such different aspects of the input activity structure, rather than each detecting some overall average activity structure, even while each detects the same orientation? Our current understanding of models of development under activity-dependent rules of synaptic plasticity does not provide a good answer to this question. Existing proposals for forcing different cells to learn differing aspects of the input include the following. *1*) Recurrent connections between cortical cells that are modified by anti-Hebbian learning rules (Foldiak 1990): such connections have not been observed, and they would not explain the learning of a common preferred orientation or other clusterings or periodic organizations of response properties. *2*) Combinations of plasticity of excitatory and inhibitory neurons: as implemented thus far, these either do not produce sufficient variability of properties (Kayser and Miller 2002) or do not explain clustering or periodic organization (King et al. 2013). *3*) Initial learning of orientation selectivity in phase-nonselective cells, which is then transmitted in some manner to phase-selective cells (discussed in Miller 1994) that otherwise develop in an uncorrelated manner and thus develop diversity in their properties other than preferred orientation (Antolik and Bednar 2011): as implemented thus far, this approach would not explain diversity of response properties in layer 2/3, i.e., in the cells that initially learn orientation selectivity.

The present study, by characterizing the degree of clustering and variability in multiple response properties, underlines the importance of determining both the circuit mechanisms (e.g., the connectivity patterns) that achieve the observed combinations of specificity in responses of individual cells with both clustering and diversity in the response properties of nearby cells and the developmental rules that can lead such circuits to form.

## GRANTS

This work was supported by gifts from the Swartz Foundation and by Grants R01-EY-13595 and R01-EY-11001 from the National Eye Institute. A. J. Ziskind was supported by National Institutes of Health Training Grant T32-NS-064929.

## ENDNOTE

At the request of the author(s), readers are herein alerted to the fact that additional materials related to this manuscript may be found at the institutional website of one of the authors, which at the time of publication they indicate is: http://dx.doi.org/10.7916/D8H41Q78. These materials are not a part of this manuscript, and have not undergone peer review by the American Physiological Society (APS). APS and the journal editors take no responsibility for these materials, for the website address, or for any links to or from it.

## DISCLOSURES

No conflicts of interest, financial or otherwise, are declared by the author(s).

## AUTHOR CONTRIBUTIONS

Author contributions: A.J.Z., A.A.E., A.V.K., S.P.R., and K.D.M. conception and design of research; A.J.Z., A.A.E., A.V.K., and K.D.M. analyzed data; A.J.Z., A.A.E., A.V.K., and K.D.M. interpreted results of experiments; A.J.Z. and A.A.E. prepared figures; A.J.Z., A.A.E., A.V.K., and K.D.M. drafted manuscript; A.J.Z., A.A.E., A.V.K., and K.D.M. edited and revised manuscript; A.J.Z., A.A.E., A.V.K., S.P.R., and K.D.M. approved final version of manuscript; A.A.E., A.V.K., and S.P.R. performed experiments.

## ACKNOWLEDGMENTS

We thank Tatyana Sharpee and Hiroki Sugihara for helping to collect the data. We also thank Steve Lisberger for helpful comments on early versions of the manuscript.

- Copyright © 2015 the American Physiological Society