We studied the relative accuracy of drifting gratings and noise stimuli for functionally characterizing neural populations using two-photon calcium imaging. Calcium imaging has the potential to distort measurements due to nonlinearity in the conversion from spikes to observed fluorescence. We demonstrate a dramatic impact of fluorescence saturation on functional measurements in ferret V1 by showing that responses to drifting gratings strongly violate contrast invariance of orientation tuning, a fundamental property of the spike rates. The observed relationship is consistent with saturation that clips the high-contrast tuning curve peaks by ∼40%. The nonlinearity was also apparent in mouse V1 responses to drifting gratings, but not as strong as in the ferret. Contrast invariance holds, however, for tuning curves measured with a randomized grating stimulus. This finding is consistent with prior work showing that the linear portion of a linear-nonlinear system can be recovered with reverse correlation. Furthermore, we demonstrate that a noise stimulus is more effective at keeping spike rates in the linear operating regime of a saturating nonlinearity, which both maximizes signal-to-noise ratios and simplifies the recovery of fast spike dynamics from slow calcium transients. Finally, we uncover spatiotemporal receptive fields by removing the nonlinearity and slow calcium transient from a model of fluorescence generation, which allowed us to observe dynamic sharpening of orientation tuning. We conclude that for two-photon recordings it is imperative that one considers the nonlinear distortion when designing stimuli and interpreting results, especially in sensory areas, species, or cell types with high firing rates.
- primary visual cortex
two-photon calcium imaging is a powerful technique that allows for in vivo recording from many cells in a given brain volume (Garaschuk et al. 2006; Göbel and Helmchen 2007; Kerr and Denk 2008; Svoboda and Yasuda 2006). When using calcium imaging to study neuronal output, it is important to account for the considerable transformation between the spiking activity of neurons and the observed calcium indicator fluorescence. Various factors associated with the calcium indicator and the given firing pattern influence whether or not the spike-to-fluorescence transformation can be considered linear (Hendel et al. 2008; Neher and Augustine 1992; Pologruto et al. 2004; Tay et al. 2007; Wallace et al. 2008; Yaksi and Friedrich 2006; Yasuda et al. 2004). Furthermore, the degree to which nonlinearity might alter the functional properties in a large population has not yet been studied, especially at firing rates indicative of V1 in larger mammals. If two-photon imaging is to correctly address how functional properties of different neuronal populations are shaped by the circuitry, then the measurements must accurately represent the neuronal responses. The first goal of this study was to demonstrate the impact of the nonlinearity encountered with the most commonly used calcium indicator, Oregon Green BAPTA 1-AM (OGB), on measurements of neuronal tuning curves using visual stimuli that elicit continuous firing rates. Specifically, we quantified the degree to which fluorescence saturation of the calcium indicator distorts tuning curve measurements by making use of the well-established contrast invariance of the orientation tuning of neurons in visual cortex. Contrast invariance implies that the shape of a neuron's orientation tuning curve should not change with contrast of the visual stimulus (Alitto and Usrey 2004; Nowak and Barone 2009; Sclar and Freeman 1982). Fluorescence saturation at high firing rates should cause an apparent violation of contrast invariance, such that high-contrast gratings will yield broader tuning than low-contrast gratings. When stimulating with drifting gratings, we indeed find that high-contrast tuning curves measured with two-photon imaging are significantly broader than low-contrast tuning curves. The strength of the observed contrast-dependent tuning with drifting gratings points to the measurement error that can be encountered in general when using calcium imaging to estimate neuronal firing rates.
Operating in the saturating regime of the fluorescence nonlinearity impedes the use of standard linear methods for reconstructing spike rates (Holekamp et al. 2008; Yaksi and Friedrich 2006). Removal of the fluorescence saturation to reconstruct precise spike times requires relatively sophisticated processing (Vogelstein et al. 2009). Here we demonstrate that the spike-to-fluorescence transformation can easily be accounted for when measuring receptive fields with stochastic stimuli. Linear-nonlinear systems can be fully identified with Gaussian white noise (Bussgang 1952; Chichilnisky 2001; Korenberg and Hunter 1990) and subspace noise (Benucci et al. 2009; Ringach et al. 1997a) stimuli by cross-correlating the input with the output. We applied the same framework to a biologically inspired model of imaging neuronal calcium responses as a simple means of removing the fluorescence nonlinearity from our tuning curve measurements. Consistently, our data show that tuning bandwidth is invariant to contrast when tuning curves are computed via cross-correlation with noise. We then take the analysis a step further by using the same model of fluorescence generation to estimate the nonlinearity, calcium transient, and spatiotemporal receptive field of the underlying spike rates. This method allows us to look at properties of rapid tuning dynamics across local populations.
In addition to accounting for distortion caused by a static nonlinearity, noise stimuli may improve signal-to-noise ratios by keeping the majority of responses from saturation. Unlike steady-state stimuli, a rapid stochastic stimulus will produce transient responses. The transformation from a spike train to intracellular calcium concentration ([Ca2+]) is usually modeled as the convolution between the spike train and an exponential decay function. We convolve a typical exponential function with extracellular spike data to show that the relatively sparse firing patterns elicited by stochastic stimuli are expected to maintain a far more linear operating regime than drifting gratings with two-photon imaging.
MATERIALS AND METHODS
Animal Preparation and Surgery
All procedures were conducted in accordance with guidelines of the National Institutes of Health and were approved by the Institutional Animal Care and Use Committee (IACUC) at the Salk Institute. Experiments were performed on male ferrets aged postnatal days 40–60. Ferrets were pretreated with atropine (0.05 mg/kg), and anesthesia was induced with ketamine (40 mg/kg) and maintained with isoflurane (during surgery: 1.5–2.5%, during recording: 0.75–1.25%). Ferrets were paralyzed with pancuronium bromide (0.2 mg·kg−1·h−1). Respiration was maintained with a ventilator (Inspira, Harvard Apparatus), adjusting breathing rate and volume to maintain the EtCO2 between 3.2% and 3.9%. Small silver wires inserted through burr holes over frontal cortex were used to monitor the EEG. EtCO2, SpO2, EEG synchronization, and heart rate were monitored continuously to assess the depth of anesthesia. Body temperature was maintained at 38.5°C with a heating pad. Subcutaneous injections of 5% dextrose Ringer solution (5–10 ml·kg−1·h−1) were administered to prevent dehydration. A small metal headplate was attached to the skull with dental acrylic and was screwed onto a holder mounted in a standard stereotax (David Kopf Instruments). We also attached a small well over visual cortex, inside which we performed a craniotomy (diameter 4 mm) and durotomy to expose areas 17 and 18. After dye injection, the exposed brain was covered with agarose [1.5% in artificial cerebrospinal fluid (ACSF); type III-A, Sigma-Aldrich] and a coverslip (World Precision Instruments) to reduce brain motion. Atropine and Neosynephrine eye drops were used to dilate the pupil and to retract the nictitating membrane. Animals were fitted with contact lenses.
All procedures were conducted in accordance with guidelines of the National Institutes of Health and were approved by the IACUC at the Salk Institute. Experiments were performed with two male C57BL mice ∼2.5 mo old. Surgical plane anesthesia was induced and maintained by isoflurane (2.5% induction, 1.25% surgery). A custom-made metal frame was mounted to the skull, centered on the visual cortex of the left hemisphere with stereotaxic coordinates. The bone was thinned to allow for reliable retinotopic mapping of primary visual cortex with intrinsic signal imaging. After the surgery, chloprothixene (1 mg/kg) was administered intramuscularly and isoflurane was reduced to 0.4–0.7% for visual stimulation and recording experiments. After intrinsic imaging was performed to identify the location of V1 relative to the blood vessel pattern, the remaining bone over V1 was removed to allow for two-photon imaging.
The mouse V1 boundary was identified with intrinsic imaging prior to two-photon imaging. We used a Dalsa 1M60 camera (Dalsa, Waterloo, ON, Canada). For image acquisition, we used a Matrox Helios Camera-Link frame grabber (Matrox, Montreal, QC, Canada) and the Matrox Imaging Library to create a custom interface with MATLAB to run experiments. For illumination, we used an Illumination Technologies (East Syracuse, NY) 3900 Smart-Lite illuminator. Vasculature images were obtained with a green filter, whereas intrinsic images were obtained with red illumination and a 700-nm filter with 20-nm bandwidth (58460, Spectra Physics, Mountain View, CA). Retinotopic maps were obtained by taking the phase of the Fourier component at the frequency of a drifting bar (Kalatsky and Stryker 2003). V1 was identified during dye injection by overlaying the vasculature with the functional maps.
Multicell Bolus Dye Loading and In Vivo Two-Photon Microscopy
We injected a dye containing 1 mM OGB with 10% DMSO, 2% Pluronic F-127, and 25% sulforhodamine 101 (SR101) (all from Invitrogen) in ACSF (Stosiek et al. 2003). A small volume of dye was pressure injected under visual guidance ∼200 μm below the cortical surface with a Picospritzer (Parker). Dye was injected through glass pipettes with a tip diameter of ∼3 μm. Data acquisition began ∼60 min after dye injection.
Two-photon microscopy was performed with a moveable objective microscope (Sutter Instruments) coupled to a Ti:sapphire laser (Chameleon Ultra II, Coherent). Data acquisition was controlled by a modified version of ScanImage (Pologruto et al. 2003). We used a ×40, 0.8 NA lens (Olympus) for imaging. The beam size was adjusted so that the back aperture of the objective was overfilled. For drifting grating experiments, images were acquired at a frame rate between 8 and 16 Hz. For flashed grating experiments, images were acquired at a frame rate between 16 and 32 Hz. OGB and SR101 fluorescence was excited at 920 nm, and emission was collected with a green (535 nm center, 50 nm band) and a red (610 nm center, 75 nm band) filter (Chroma).
Signals were acquired with a System 3 workstation from Tucker-Davis Technologies in conjunction with one of two types of electrode arrays from Neuronexus Technologies. We used either a 4-shank by 4-electrode configuration (A4x4-3mm-50-125-177-A16) or a 1-shank by 16-electrode configuration (A1x16-3mm-50-125-177-A16). Single units were isolated in the principal component dimensions of spike waveforms and verified based on the refractory period. On occasion, the same cell was recorded on adjacent electrodes. This was determined by observing the cross-correlogram of the spike trains.
Visual stimuli were generated with the Psychophysics Toolbox extensions for MATLAB (Brainard 1997; Pelli 1997) on a 17-in. LCD monitor with a refresh rate of 60 Hz. The monitor was gamma corrected with a Photo Research-701 spectroradiometer. The monitor distance was 25 cm for ferrets and 15 cm for mice.
For the mice, the drifting grating experiments were run at full field and always set at a spatial frequency of 0.04 cyc/° and a temporal frequency of 1.5 cyc/s. The temporal frequency shown to ferrets was always 3 cyc/s. For the ferret experiments alone, we ran a set of preliminary stimuli to help optimize the spatial frequencies and retinotopic location of the stimulus for each region of interest (ROI) that was imaged. This required three short experiments. The first was a large (90° × 90°) drifting square-wave grating stimulus at eight orientations, which gave the preferred orientation(s) within the ROI. In the next experiment, we presented sine-wave drifting gratings at four spatial frequencies (0.04, 0.08, 0.16, and 0.32 cyc/°), each at one or two directions. Finally, we ran a retinotopy stimulus, which consisted of a drifting grating within a rectangular window, at the optimal orientation and spatial frequency for the population. The size of the grating window was chosen such that it covered the entire extent of the screen along one axis and about one-fifth of the screen along the other. We systematically changed the position of these gratings along their shorter axis, thereby mapping out responses to all possible stimulus locations on the screen. This allowed us to quickly approximate the screen location of the receptive fields. The stimulus in subsequent experiments was between 15° and 25° in diameter.
Drifting gratings were shown at 2 contrasts and 12 or 16 directions, at the optimal spatial frequency for ferrets and at 0.04 cyc/° for mice. For the ferrets, the low contrast was usually set to be 15%, but in some ROIs it was necessary to increase it up to 25% to get responses at a signal-to-noise ratio that would yield accurate tuning curves. For mice, the low contrast was set to 22%. Both orientation and contrast were randomly interleaved across trials, with a gray screen shown between trials. A gray screen was also presented instead of a grating on a few random trials. Each stimulus was shown for 3–4 s.
The noise stimulus was randomly pooled from a set of gratings of different orientation, spatial frequency, and spatial phase. The domain of orientations and spatial frequencies were organized into Cartesian coordinates. The stimulus was updated every 333 ms, and each grating was shown for the entire time within this window. A gray screen was randomly shown instead of a grating 5% of the time. The spatial frequency range was defined by the preliminary experiments using drifting gratings. Each trial lasted 60 s and contained only one of the two contrasts.
For each ROI, we created a cell mask defining the boundary of every cell in the ROI. Creating cell masks consisted of first implementing an automated procedure, followed by a manual (un)selection tool created in MATLAB (MathWorks) to correct for any obvious errors. The automated procedure entailed first performing a local Z-score filter on the entire image, followed by thresholding the image to make it binary, and finally morphological “opening” to clean the mask. To account for glia, cells were removed that were present in both the red and green channels. For each cell, we computed its responses as the mean of the pixel values within its boundary. To account for slow x/y drift of the brain over the course of an experiment, we aligned the images on a trial-by-trial basis. The location of the brain on each trial was aligned to the location on the first trial by finding the x/y location that maximized the correlation between the mean images. The SR101 channel was used for this purpose.
Responses to drifting gratings at each orientation were averaged over repeats and within the stimulus presentation period (Fig. 1). The response for each trial was normalized by the response to the gray screen immediately before the stimulus came on. For example, the response time course for each trial was computed as [F(t) − Fo]/Fo, where Fo is the mean response to the gray screen prior to stimulus onset for the given trial.
Prior to computing the kernels for the flashed gratings, we used a conservative band-pass filter to help remove noise—the signal for each neuron was convolved with a difference-of-Gaussian function (σlowpass = 50 ms; σhighpass = 5,000 ms). Responses were then Z-scored on each 60-s trial. The kernels were computed as the mean fluorescence (normalized Z units) response, given the stimulus parameters of the presented grating. Unlike the drifting gratings, the noise stimulus had extra dimensions in spatial phase and spatial frequency. The kernel for each neuron was thus a function of all these parameters. To compute the orientation tuning curve for each neuron, responses were first averaged over spatial phase. Next, the resulting kernel was smoothed in the spatial frequency dimension. Finally, the orientation tuning curve was computed at the time delay and spatial frequency of maximum response.
Our main objective in selecting cells for further analysis was to ensure that the tuning curves from low-contrast gratings yielded accurate measurements. Cell selection was based on two criteria. The first criterion was that the low- and high-contrast tuning curves for a given cell had to have the same peak to be included. The second criterion was based on a metric for the signal-to-noise ratio of the low-contrast tuning curve. For low-contrast gratings, we computed the response mean (μ) and standard error (SE) both to the optimal orientation and to a gray screen that was randomly interleaved in all experiments. We included cells that had (μopt − μgray)/(SEopt + SEgray) > 1.
To quantify the degree of saturation with contrast invariance (Figs. 2, 3, and 6), the function H = A(1 − e−αL) was fit to the scatter of points for each ROI, where H and L are the responses to high and low contrast, respectively. A and α are the free parameters, which were fit by minimizing the sum of the errors. Error for each data point was computed as the minimum distance to the function. The function derivative at the origin is Aα, which was used to quantify the amount of saturation.
Fluorescence Saturation Distorts Orientation Tuning Curves Obtained with Drifting Gratings
We used contrast invariance of orientation tuning to demonstrate the impact of fluorescence saturation on tuning curves acquired with calcium imaging. Neurons in ferret V1 were loaded with OGB (Stosiek et al. 2003), and their orientation tuning curves were measured by presenting drifting gratings of two different contrasts and monitoring the induced changes in OGB fluorescence with a two-photon microscope. The high-contrast condition always used a contrast of 90%. The low contrast was usually set to 15%, as this has been shown previously to elicit half the maximum spike rate for most complex cells in ferret V1 (Alitto and Usrey 2004). However, the low contrast was increased to 25% in a few experiments where cells did not respond reliably to 15%. In addition to the two-photon experiments, we performed extracellular recordings in ferret V1 using metal electrodes. Again, orientation tuning curves were measured for gratings of two different contrasts.
For both electrophysiology and two-photon data sets, a Gaussian function was fit to each tuning curve, and σ from the fit was used to quantify tuning bandwidth (see Fig. 1C for an example fit for the calcium imaging data). For each cell, we then quantified the dependence of tuning width on contrast by computing the log of the bandwidth ratio for the two contrasts, log10(σhi/σlow). For the spike data, the resulting distribution of values shows that tuning widths are independent of contrast (Fig. 1, F and G) (geometric mean of 1.04; P = 0.45, t-test). Thus our data are in agreement with a previous study demonstrating contrast invariance of orientation tuning in ferret V1 (Alitto and Usrey 2004). Performing the same analysis on the calcium imaging data yields very different results. A visual comparison of bandwidth between the two contrasts across a population of 115 neurons from 4 animals indicates that the tuning curves are broader for the higher-contrast drifting gratings (Fig. 1D). Indeed, the distribution of the log of bandwidth ratio (Fig. 1E) shows that the average high-contrast bandwidth was 22% larger than the low-contrast bandwidth (geometric mean of σhi/σlow = 1.22; P < 10−10, t-test). The observed contrast-dependent change in tuning width in the calcium imaging data is as expected from fluorescence saturation at high firing rates.
Next, we used the same ferret V1 two-photon imaging data to quantify the amount by which high-contrast tuning curves are getting clipped by the saturating nonlinearity. This estimation rests on the assumption that spike rate tuning curves for the average cell are contrast invariant in overall shape, meaning that tuning curves for different contrasts are simply scaled versions of one another. Since firing rates near the base of the tuning curves are low, we additionally assume that these portions of the tuning curves can accurately be determined by calcium imaging, without distortion by the fluorescence saturation. Taken together, these assumptions mean that the portions of the calcium imaging tuning curves near the base will reflect the underlying scaling relationship that exists between tuning curves for different contrasts. We then quantify the clipping as the deviation from this linear relationship at the tuning curve peak. This idea can be made more concrete by plotting the OGB responses to high- and low-contrast gratings, measured in response to each stimulus orientation, against each other in a scatterplot. Figure 2 shows data for an example cell. These points will sit on a line if the cell is contrast invariant. However, if responses are saturating, they will form a concave curve where the trend near the origin represents the contrast invariance of underlying spikes. Because it is difficult to make a robust estimate of the saturation for each cell independently, we pooled the data for neurons assayed simultaneously within a given imaging ROI (Fig. 2, C and D). Within each ROI, the tuning curves for each neuron and contrast were first fit by a Gaussian, as before. Tuning curves were then normalized by subtracting the Gaussian baseline and dividing by the amplitude. To quantify the degree of saturation, the normalized high- and low-contrast responses of each neuron at each stimulus orientation were plotted against each other for all neurons within an ROI (Fig. 2C). These data were then fit with a saturating exponential function. The slope of the fit at the origin allows us to quantify the saturation (materials and methods). Without saturation, the slope will be unity because of the normalization. A higher slope indicates a more concave curve, and therefore a more dramatic saturation. For the example ROI shown, the slope is 1.40, which means that the high-contrast response would reach 1.4 (at a low contrast of 1.0) if it were not for clipping. This estimates that the high-contrast tuning curve amplitudes are getting clipped by 29% [i.e., 29 = 100 × (1.4 − 1)/1.4] of their actual values. To help illustrate this measure, we first plotted the mean tuning curve (with standard error), after normalization, for the low-contrast condition (Fig. 2D). The corresponding tuning curve for the high-contrast condition is shown after dividing by the slope of 1.4. As expected, the trends for the two contrasts match at the base of the tuning curve but diverge toward the peak. The saturation is quantified based on the degree of this divergence. This analysis was performed on the same four ROIs (4 animals) as the two-photon data from Fig. 1. The number of neurons for each ROI was 29, 11, 21, and 54. The four respective slopes were 1.40, 1.66, 1.85, and 2.11, corresponding to peaks that are clipped by 29%, 40%, 46%, and 51% of the actual tuning curve. Because these estimates of saturation are based on the relative difference in tuning shape for the two contrasts, any saturation at lower contrasts puts an upper bound on the measured saturation. In other words, these saturation values may be an underestimate.
Previous work has shown that the degree of fluorescence saturation can be used to acquire estimates of intracellular calcium concentration ([Ca2+]) (Maravall et al. 2000; Yasuda et al. 2004). We used this framework to provide estimates of [Ca2+] based on our computed clipping percentages. In particular, Yasuda et al. (2004) derived [Ca2+] as a function of a “nonlinearity” variable (NL): where FL is the linear approximation of the fluorescence using the slope at [Ca2+] = 0, FNL is the actual fluorescence, Fmin is the fluorescence at [Ca2+] = 0, and KD is the dissociation constant. We have reexpressed this equation as something slightly different from NL yet more directly related to the clipping percentages (CP) defined above: To compute [Ca2+] from the above equation based on the clipping percentage from our data (Fig. 2), we first make the assumption that Fmin is similar to the fluorescence during the baseline response. This is already implied by the fact that we have fit a saturation function to the data that does not have a baseline (i.e., it crosses the origin in Fig. 2C). By making this assumption, the fluorescence ratio in the above equation is effectively the same as the clipping percentages. Next, we assume a value of 0.21 mM for OGB's dissociation constant, KD (Yasuda et al. 2004). This yields [Ca2+] equal to 86, 140, 179, and 219 nM based on the clipping percentages computed from the four ROIs that were stimulated with drifting gratings (given above). More specifically, these [Ca2+] values correspond to the average concentration during the presentation of the optimal high-contrast gratings. On the basis of previous estimates of change in [Ca2+] for a single action potential (with Ca2+ buffer), we can also find loose estimates for the firing rates. From the [Ca2+] estimates during a constant barrage of spikes, the spike rate can be computed as [Ca2+]/(Δ[Ca2+] × τ), where Δ[Ca2+] is the amplitude of the change in [Ca2+] after an action potential and τ is the time constant of the [Ca2+] transient with OGB. We used 50 nM and 150 ms for these constants (Yasuda et al. 2004). This gives 11.5, 18.7, 23.9, and 29.2 spikes/s in response to the optimal orientation at high contrast, which are sensible values for ferret V1 (see Fig. 8).
Contrast-dependent changes in fluorescence responses to drifting gratings were not unique to the ferret visual cortex but were also observed in mouse V1 (Fig. 3), where published data indicate that spiking responses are contrast invariant (Niell and Stryker 2008). The two contrasts shown to the mice were 22% and 90%. Similar to the analysis performed on the ferret data, we quantified the change in tuning bandwidth using the log of the bandwidth ratio between the two contrasts (Fig. 3, A and B). For the two mice analyzed (66 cells), the geometric mean of the bandwidth ratio was 1.16 and significantly greater than 1 (P ∼ 0.0006). Identical to the analysis performed on the ferret data, we quantified the saturation by fitting an exponential to scatterplots that compare responses between the two contrasts. These response values are shown for an example ROI (Fig. 3, C and D), which had a saturation estimate of 30%. The slopes were computed independently for three ROIs (2 mice), which each had 13, 25, and 28 neurons. The respective slopes were 1.07, 1.42, and 1.56, which correspond to peaks clipped by 7%, 30%, and 36%. On average, these saturation values are not as high as those computed for the ferret. This difference may be due to lower firing rates in the mice.
Using Noise Stimuli to Measure Spatial Properties of Visual Receptive Fields with Two-Photon Calcium Imaging
We have shown that ignoring the nonlinearity inherent in the spike-to-fluorescence transformation can lead to inaccurate tuning curve measurements with calcium imaging. In the next set of experiments, we therefore aimed to establish a method to account for the nonlinearity. For this purpose, we used a subspace noise stimulus to measure tuning curves (materials and methods) (Ringach et al. 1997a). The stimulus consisted of a sequence of flashed gratings randomized for orientation, spatial frequency, and spatial phase. Each grating was on the screen for 333 ms and was immediately followed by the subsequent grating. To aid our interpretation of tuning curve measurements with this stimulus, we modeled the fluorescence response as the output of a linear-nonlinear cascade (Fig. 4A). The first linear element in this cascade is the neuron's receptive field as a function of the parameters in the stimulus ensemble, such as orientation. Our approximation of a purely linear receptive field is based on the slow presentation rate of our stimulus in relation to the temporal integration time of the receptive field (see Fig. 9, A and B). The tuning curve thus incorporates the cell's spike output nonlinearity. A more detailed justification for this linear approximation is provided below. The second linear element represents the impulse response of the calcium influx caused by spikes. The nonlinearity instantaneously converts the calcium signal into the observed fluorescence. The entire cascade can be treated as one linear stage followed by the nonlinearity, which allows for simple means of identifying the linear and nonlinear sections with noise stimuli. For example, computing the expected fluorescence response to each grating in the stimulus ensemble yields a kernel (Fig. 4B) that is theoretically proportional to the combined linear front end, which we will call the “Ca2+ kernel” (materials and methods; see also Benucci et al. 2009; Ringach and Malone 2007 for a similar analysis with subspace noise stimuli). A requirement for this technique is that there exists sufficient temporal and/or spatial summation within the Ca2+ kernel. This is achieved by the slow response dynamics of the Ca2+ signal, relative to the presentation rate of the stimulus (see Fig. 9, C and D). The Ca2+ kernel is the temporal convolution between the Ca2+ transient and the cell's receptive field. Because the original receptive field has been smoothed only in time, its spatial properties, such as orientation tuning bandwidth, are mostly preserved. Computing spike rate tuning curves from drifting gratings, for example, would provide similar information about the cell. We show below how to recover the spatiotemporal dynamics of the receptive field.
Similar to the analysis with drifting gratings, we used contrast invariance to measure the impact of the fluorescence nonlinearity on tuning curve measurements using the noise stimulus. By cross-correlating the noise stimulus with the fluorescence response, we estimated the Ca2+ kernel for two stimulus contrasts in four animals (four ROIs, 85 cells). Two of the four ROIs are the same as those used for the analysis with drifting gratings (Fig. 1, D and E). Orientation tuning curves were computed by taking a slice of the kernel at the optimal spatial frequency and time delay (Fig. 5C). We tested for contrast invariance of tuning curve bandwidth by applying the same analysis used for the drifting gratings. The average tuning bandwidth from the noise stimulus was invariant to contrast (geometric mean of σhi/σlow ∼ 0.98; P ∼ 0.51) (Fig. 5, D and E). We also used tuning curves from these same experiments to estimate clipping, as was done for the drifting grating data. Figure 6 plots the data for an example ROI—the same ROI used in Fig. 2, C and D, for the analysis based on drifting gratings. The contrast independence is clear from the linearity observed in the scatterplot (Fig. 6A) and the similarity between the mean tuning curves (Fig. 6B). The slope of the fitted function at the origin is 1.05 (∼5% saturation). The slope for each of the ROIs analyzed was 1.05, 0.97, 1.06, and 1.10 (N = 17, 28, 18, and 15 neurons, respectively).
Single-unit recordings using the same noise stimulus also yielded contrast-invariant orientation tuning: We recorded single-unit activity in response to the same flashed grating noise stimulus that was used with two-photon imaging and compared tuning functions for high- and low-contrast gratings. To provide a more accurate comparison against the tuning curves computed from calcium imaging data, each point on the tuning curves was measured as the average response within a temporally broad window. We chose 40 ms to 200 ms as the averaging window, as this generously accounts for the entire response after each stimulus (see Fig. 9, A and B). The σ of the Gaussian fits to these spike rate tuning curves are invariant to contrast (Fig. 5, F and G).
Finally, the results from the two different calcium imaging experiments suggest that high-contrast drifting gratings should yield broader tuning curves than high-contrast flashed gratings. We tested this prediction by analyzing data from a population of neurons that were driven by both drifting and flashed gratings at full contrast (Fig. 7). The log of the bandwidth ratio for the two stimuli, log10(σdrifting/σflash), was used to quantify the dependence. Bandwidth for drifting gratings was 16% larger than for flashed gratings (geometric mean of σdrifting/σflash ∼ 1.16; P ∼ 0.0007). See Nishimoto et al. (2005) for a similar analysis showing that spike rate tuning curves of cat V1 neurons have similar bandwidths for the two stimulus types.
A Different Operating Point Within the Fluorescence Nonlinearity for Drifting Gratings and Noise Stimuli
In the previous section, we showed how noise stimuli can be used for calcium imaging to compute a kernel that is proportional to the receptive field. However, a remaining pitfall of the saturating nonlinearity is that it will potentially degrade signal-to-noise levels by attenuating the larger responses. In this section, we provide additional data from electrophysiology comparing spike responses to drifting gratings versus our noise stimulus. The data show that drifting gratings drive spiking with different temporal dynamics than the noise stimulus. These differences are likely to cause greater nonlinearity in the calcium-to-fluorescence transformation for the drifting gratings. In this experiment, we used metal electrodes to record the responses of V1 neurons to the same noise and drifting grating stimulus used in our two-photon imaging experiments. For each cell, we first identified the optimal grating. We then computed the firing rate (computed within 5-ms bins) to that stimulus as was done for the two-photon imaging data. That is, we took the response at the optimal spatial frequency and time delay for the flashed gratings, and we took the average within the stimulus window for the drifting gratings. This measure of firing rate response tends to be similar for the two different stimuli (see Fig. 8A; slope of linear fit to the data equals 0.93). We hypothesized that, despite this similarity in the expected spike rate for the optimal grating, noise stimuli might keep [Ca2+] at lower levels because of reduced temporal summation across overlapping calcium influx transients evoked by each spike. Lower [Ca2+] would get transformed by a Ca2+-to-fluorescence nonlinearity within a relatively linear regime. By convolving the observed firing rates with an exponential decay function (τ = 100 ms) we simulate the expected calcium signal that is the input to the fluorescence nonlinearity. We find that the modeled calcium transients accumulate to much higher values with drifting gratings than with noise stimuli (Fig. 8B), implying that there would be greater saturation with drifting gratings. The points form a fairly linear trend, so we quantified the difference in [Ca2+] between the two stimulus conditions using the slope (2.66) as a multiplicative factor. For example, the trend implies that the Ca2+ concentration reaches values that are 2.66-fold higher for the drifting gratings. The slope will change systematically with the time constant used, getting steeper for slower transients. This simple analysis shows that the firing patterns elicited by different stimuli have a dramatic impact on the operating point within the fluorescence saturation. The notion that noise stimuli keep [Ca2+] within a linear regime is also exemplified by the linearity of the scatterplot shown in Fig. 4C.
Justification for the Linear-Nonlinear Model of Fluorescence Generation
In the preceding analyses of the noise stimulus experiment, we had assumed that fluorescence responses in V1 can be modeled as the output of a linear-nonlinear cascade, which would allow us to recover the linear front end via cross-correlation with the response (Fig. 4). This rests on two key assumptions. The first is that the V1 receptive field can be modeled as purely linear (Fig. 4A), as opposed to a linear receptive field followed by a spiking nonlinearity. The second assumption is that the noise stimulus has spatiotemporal properties that are sufficient to recover both elements of the linear-nonlinear system. Below, we provide data to justify both assumptions. In short, the first assumption is justified by the short temporal integration time of the V1 receptive field, relative to the update period of the stimulus. The second assumption is justified by the long integration time of the calcium responses, relative to the update period of the stimulus.
To validate the first assumption, we characterize the spike rate response transients after the onset of the flashed grating noise stimulus. For each cell, we compute the average spike rate time course for three locations on the orientation tuning curve, at the optimal spatial frequency. Time courses were computed for the best orientation, for an orientation ±1σ away from the peak orientation, where σ is computed from the Gaussian fit, and for the orientation 90° from the optimal orientation. Figure 9, A and B, shows the resulting response transient for the single units, to demonstrate how the spike rates evolve within the time window of the stimulus period (333 ms). The average response decays to zero prior to the onset of the next stimulus. This lack of temporal integration of the stimulus within the receptive field warrants the purely linear description of the receptive field, as it allows one to consider the actual linear-nonlinear receptive field as a single linear component that has effectively assimilated the spike nonlinearity. To formulate, we first describe the transformation from visual stimulus to spike rate as
The spatiotemporal linear front end of the receptive field is defined as V(t,θ), where V is the generator potential, t is time, and θ is the stimulus parameter. For a linear system, the temporal transients add together across subsequent stimuli, which are indexed by i and are spaced in time by Δt. The spike rate is then generated by applying a spike nonlinearity, NL. In our case, Δt = 333 ms. When V(t), the transient response to a flashed grating, is much more narrow in time than Δt, there is effectively no summation of the transients across adjacent presentations, i. The lack of summation means that we can place the nonlinearity inside the summation term. The transformation from visual stimulus to spike rate can thus be approximated as Here, VNL = [V]NL. It is now a purely linear equation, where the new linear kernel is the original kernel, but processed by the nonlinearity. Tuning properties of VNL are theoretically the same as those obtained with a steady-state stimulus, such as drifting gratings. The entire cascade, from visual stimulus to fluorescence, can thus be approximated as a linear-nonlinear system as we have defined it in Fig. 4A. The linear stage of this system is the receptive field (VNL) convolved with the calcium transient, and the nonlinearity is the fluorescence saturation.
To test whether the noise stimulus indeed has spatiotemporal properties that are sufficient to recover both elements of the linear-nonlinear system, we performed the same time course analysis on data recorded with two-photon imaging. The time courses from the calcium responses were upsampled to the rate shown in Fig. 9, C and D (32 samples/s). This was done so that we could average time courses across multiple experiments that had slightly different acquisition rates, in addition to providing more resolution in our estimates of the dynamics. For both signals (Ca2+ and spikes), we computed the mean and standard deviation of time to peak (ttp) and full-width at half-maximum (FWHM) of the response to the optimal orientation. Table 1 summarizes these results. The relatively long integration time of the Ca2+ signal (FWHM = 741 ms for high contrast) should provide sampling density within the Ca2+ kernel such that corruption by the fluorescence saturation is minimized.
Reconstructing Spatiotemporal Receptive Fields with Noise Stimuli
In a previous section, we computed the “Ca2+ kernel,” which is the temporal convolution between the calcium transient and a cell's receptive field (Fig. 4A). In this section, we demonstrate how the receptive field alone can be estimated with noise stimuli. The first stage of this analysis is to simply compute the Ca2+ kernel as before by taking the expected response to each stimulus in the ensemble. Next, we compute the output nonlinearity by comparing the predicted output of this linear kernel to the actual output (Fig. 4, C and D) and fitting an exponential saturation function to the scatterplot. The output of the linear portion ([Ca2+]) can then be recovered by applying the inverse of the fitted function to the observed fluorescence. In most cases, we find that the computed nonlinearity is not very strong, such that the response appears to stay mostly within a linear operating regime. As shown in a previous section, the responses of a given cell to a noise stimulus are expected to stay more confined to the linear operating regime of a saturating nonlinearity than its responses to drifting gratings (Fig. 8). This made the inversion process quite stable. In some cases (Fig. 4C), the nonlinearity could effectively be ignored.
Once we have the input and output of the linear portion of the cascade model, estimating its parameters is relatively straightforward. We defined this purely linear cascade with the following equation: Here, si and Ci are the stimulus vector and calcium concentration at frame i and H is a matrix that defines the receptive field. The calcium exponential function is defined with a first-order difference equation. a is a scalar between 0 and 1 and is related to the exponential time constant of the calcium transient by τ = −T/ln(a), where τ is the time constant and T is the sample period of the two-photon data acquisition. To estimate the free parameters in the model, H and a, we used least-squares regression. This analysis is especially applicable when scan rates are sufficient to resolve the temporal dynamics of the receptive field. We demonstrate the results for a region of V1 that was scanned at 32 frames/s. The kernel, H, for each neuron, is a function of orientation, spatial frequency, and time. The orientation and spatial frequency preference maps given by these kernels, at the optimal time delay, are shown for the imaged region in Fig. 10, A and B. The kernels for three of the cells are shown in Fig. 10, C–F. For each cell, we took cross sections of the kernel at both the optimal orientation and the optimal spatial frequency to show the tuning dynamics for each parameter. The images show the response as a function of one stimulus parameter and time (Fig. 10, C and D). The curves show how the tuning changes before and after the time delay of maximum response (Fig. 10, E and F). The average tuning curves for the entire population (n = 18 cells), at three time delays, demonstrate how the average output for this small population evolves over time (Fig. 10, G and H). Consistent with previous studies in monkey V1, where receptive fields were measured based on spiking responses and reverse correlation (Ringach et al. 1997b), orientation tuning is sharper at the tail end of the response than near the peak. Estimating spatiotemporal receptive field properties with two-photon calcium imaging allows for measurements of the temporal dynamics of functional maps within unique functional compartments.
Exploiting Contrast Invariance to Measure the Impact of Fluorescence Saturation
We examined the influence of fluorescence saturation on functional measurements with two-photon calcium imaging. Our results show that in ferret V1 the nonlinearity in the spike-to-fluorescence transformation strongly distorts tuning curves measured in response to drifting gratings. Specifically, we find that with the use of drifting gratings, orientation tuning curves determined from calcium imaging data are not contrast invariant. We further demonstrate that with the use of flashed gratings and reverse correlation methods to extract receptive fields from two-photon calcium imaging data, orientation tuning is contrast invariant. Contrast invariance of orientation tuning proved to be particularly useful for demonstrating and quantifying nonlinearities in spike-to-fluorescence transformations. First, because orientation tuning has been a model for understanding how responses are shaped by the cortex, its properties are well described. To build on the decades of valuable research in V1, two-photon calcium imaging should at the very least meet fundamental principles established by studying spike rates. Second, there are clear predictions for how the fluorescence saturation encountered in calcium imaging will alter contrast invariance. These predictions allowed us to quantify the degree to which responses are distorted. Last, although simultaneously recording spikes may provide a more direct estimate of the nonlinearity for each cell, our measurements are taken from a large population of neurons and therefore yield strong statistics.
An alternative explanation for the lack of contrast invariance of the calcium imaging tuning curves measured in response to drifting gratings is that our recordings come from a population of neurons that are not contrast invariant. For example, two-photon imaging may be biased to sample from a different population of neurons than occurs from biases introduced by using extracellular electrodes, and the V1 population sampled by imaging may not be contrast invariant. However, if this were the case, we would expect a similar result when the cells were stimulated with the flashed gratings. Contrary to this prediction, the same neurons whose fluorescence responses to drifting gratings lacked contrast invariance displayed contrast-invariant fluorescence changes to flashed gratings.
Firing rates of neurons in ferret visual cortex are high, even under anesthesia. It is likely that the fluorescence saturation will have a smaller impact on two-photon experiments in species with lower firing rates such as rodents (Niell and Stryker 2008). However, we find that even in anesthetized mice OGB responses to drifting gratings are not contrast invariant, and can broaden receptive fields by as much as 16%.
Nonlinear Scaling Creates Biased Measurements
Two-photon calcium imaging is well suited for characterizing how functional properties differ within and between populations of neurons. Examples include different cell types, brain areas, or locations within functional maps. With more basic questions on functionality, such as how tuning preference (e.g., the location of the peak) differs between populations, a monotonic nonlinearity may be acceptable. However, nonlinearities create problems when measuring the “shape” of the response to a particular stimulus set, which is often more relevant when trying to understand how responses arise from the neural circuitry. For example, previous theoretical (Mechler and Ringach 2002) and experimental (Priebe et al. 2004) work has demonstrated that the bimodal distribution of the modulation ratio, F1/F0 (Skottun et al. 1991), used to classify cells as simple or complex, can be accounted for by the nonlinearity in spike generation. The theoretical results apply not only to F1/F0 but to vector-strength metrics in general. An appropriate example is circular variance, which is often used to quantify orientation selectivity. Nonlinearity in calcium imaging could therefore artificially enhance or diminish the observed functional differences between the populations under study, depending on how they might vary in other factors such as their spike rates and spiking dynamics (e.g., bursting) that can systematically alter the nonlinear transformations from spike rate to fluorescence. Such fundamental differences in the properties of each neuron can create biased results if the saturating nonlinearity with calcium imaging is not accounted for.
Maintaining a Linear Operating Regime
We demonstrated that calcium imaging can benefit from noise stimuli because the responses yield temporal spiking patterns that are expected to produce less temporal summation of calcium concentration and therefore less fluorescence saturation. We showed with extracellular spike data that drifting gratings may not necessarily drive V1 cells to greater maximum firing rates, but the combination of persistent spike rates and slow calcium dynamics will likely produce higher calcium concentration and thus more fluorescence saturation. The more a signal gets saturated, the more difficult it is to remove the saturation. For example, applying the inverse of a horizontal asymptote produces values at infinity. Reconstruction of precise spike times is much simpler if the nonlinearity is easily removed or nonexistent. Furthermore, having the signal in a more saturating dynamic range reduces the fluorescence change in response to a particular change in firing rate, which will naturally reduce the signal-to-noise ratio. Another important factor in establishing the dynamic range is the subspace of the chosen stimulus set. The noise stimulus in this study included a range of spatial frequencies that will drive most cells in ferret V1.
In addition to the design of stimuli for a given study with two-photon imaging, an important consideration is the calcium indicator used. Different indicators can yield unique signal properties related to temporal resolution, spatial resolution, and nonlinearity. OGB has been the most commonly employed indicator for in vivo studies such as this one. It is quite likely that the saturation observed here, for example, would be less dramatic if we had used a calcium indicator with lower affinity. However, there is a trade-off between saturation and signal-to-noise ratio that seems correlated with affinity. For example, the high affinity of OGB may provide the advantage of allowing for improved detection of single spikes with low firing rates (Grewe et al. 2010; Vogelstein et al. 2009) but result in stronger saturation at high firing rates. However, Wallace et al. (2008) demonstrated single spike detection with a genetic Ca2+ sensor. Another indicator-dependent variable to consider is the ability to scan neural processes, which is likely to require a genetically encodable indicator. All of these factors need to be considered in designing a two-photon experiment.
Benefits of Noise Stimuli
Although the fluorescence signal may stay within a more linear operating regime with noise stimuli than with drifting gratings, saturation is definitely still a possibility. However, noise stimuli have convenient properties for identifying and avoiding nonlinearities when the system can be modeled as a linear-nonlinear cascade (Bussgang 1952; Chichilnisky 2001; Korenberg and Hunter 1990). We modeled the receptive field and spike-to-Ca2+ transformation together as the linear front end (“Ca2+ kernel”) and the calcium-to-fluorescence transformation as the static output nonlinearity. Cross-correlation between the stimulus and fluorescence will theoretically produce the Ca2+ kernel, thus bypassing the nonlinear distortion of the spatial domain of the receptive fields. This means that measurements of functional properties such as tuning curve bandwidth will not be affected by the nonlinearity.
A few previous studies have applied noise stimuli with two-photon calcium imaging (Ramdya et al. 2006; Smith and Häusser 2010). In addition to the advantages described here that are particular to calcium imaging, noise stimuli have general advantages over steady-state stimuli for studying a population of receptive fields. For one, they make it possible to precisely map the layout of receptive fields in retinotopic space (Smith and Häusser 2010). The subspace stimulus in this study has the potential to make similar measurements if the cells have nonoverlapping on-off subregions (Ringach et al. 1997a). However, we were unable to reconstruct receptive fields in this way since most of the cells we recorded from have overlapping subregions, as expected for complex cells in superficial cortical layers. Noise stimuli also allow for measurement of the spatiotemporal dynamics of receptive fields. In addition to looking at the precise spatial tiling of receptive fields, this makes it possible to analyze the relative timing of the temporal dynamics within the population. The microarchitecture of the retinotopic layout, combined with the temporal dynamics of receptive fields, could provide key insights into how neural populations represent visual scenes.
In addition to the standard analysis of cross-correlation with a noise stimulus, we have implemented a simple method for recovering each cell's receptive field with two-photon calcium imaging. It should be noted that the method does not reconstruct precise spike times but simply the expected spike rate for a given stimulus and fluorescence response at particular times after the stimulus. Other studies have used relatively complex algorithms to extract spike times from the observed fluorescence. One study has done so by fitting a model that includes both a linear receptive field and an output nonlinearity to account for the fluorescence response (Vogelstein et al. 2009), very similar to the model we have described here. We believe that the simplicity of our approach, which exploits properties of noise stimuli, will provide a very general method for evaluating questions of receptive field properties with two-photon calcium imaging, not only in response to visual stimuli but also for other sensory modalities.
Our fastest image scanning was at 32 frames/s, which is at the lower limit for properly resolving the dynamic properties of V1 receptive fields (Ringach et al. 1997b). Furthermore, noise removal is much more difficult at this lower limit because the noise is aliased into the same frequency band as the signal. For example, any amount of low-pass filtering of the raw signal will also smooth the receptive field. In turn, reconstructing the spatiotemporal receptive field with calcium imaging would benefit from some of the recent advances in faster scan rates (see, e.g., Grewe et al. 2010; Reddy et al. 2008; Valmianski et al. 2010).
Summary and General Impact
In this study, we demonstrated that nonlinearities inherent in two-photon imaging can affect the measurement of tuning curves and need to be considered carefully when interpreting two-photon data. As the experimental scope of two-photon imaging continues to grow, we feel that the issues addressed here will become even more important. We expect that sensory stimulation of various cell types and brain areas, under different brain states, will elicit responses that are strong enough such that the nonlinearity becomes relevant in any species. For example, Niell and Stryker (2010) showed that V1 responses to visual stimuli in awake mice are more than twofold higher when the mouse is running than when it is stationary. The effect was similar for broad versus narrow spiking neurons, which also have very different firing rates.
A number of studies have examined the linearity of the spike-to-fluorescence transformation by patching on to cells, not all of which are consistent. Some of these are methods-driven and examine multiple indicators (Hendel et al. 2008; Pologruto et al. 2004; Yaksi and Friedrich 2006; Yasuda et al. 2004) while others are results-driven with OBG (Hofer et al. 2011; Kerlin et al. 2010; Kerr et al. 2005), with the latter leaning more toward a linear relationship between OBG and spikes. These patching studies provide data on a small number of cells, which may partially explain their inconsistency. One advantage our study provides is a data yield that would be extremely tedious to obtain by patching, thus providing strong statistics. Furthermore, the previous studies that identified the spike-to-fluorescence relationship as nonlinear did not provide any concrete quantification of the effect on specific functional measurements, such as orientation bandwidth. Perhaps most importantly, our results contribute to our understanding of how these distortions affect parameter distributions from the large populations provided by two-photon imaging. In addition, the present study provides the first results on how the saturation can change under different stimulus conditions. Finally, we showed that the use of contrast invariance gives investigators studying V1 with two-photon imaging a simple handle on the nonlinearity in their own measurements from a small number of animals.
This work was supported by National Eye Institute Grants EY-010742 to E. M. Callaway and EY-019821 to I. Nauhaus and an award from the Kavli Institute for the Brain and Mind (KIBM).
No conflicts of interest, financial or otherwise, are declared by the author(s).
Author contributions: I.N., K.J.N., and E.M.C. conception and design of research; I.N. and K.J.N. performed experiments; I.N. analyzed data; I.N., K.J.N., and E.M.C. interpreted results of experiments; I.N. prepared figures; I.N. drafted manuscript; I.N., K.J.N., and E.M.C. edited and revised manuscript; I.N., K.J.N., and E.M.C. approved final version of manuscript.
We thank Dario Ringach and E. J. Chichilnisky for helpful comments on the manuscript.
- Copyright © 2012 the American Physiological Society