Visual information processing in the retina has been actively investigated by exploring the relationship between light stimulation (input) and retinal ganglion cell (RGC) spiking (output). RGCs fire action potentials in response to a change in a specific visual stimulus (i.e., light contrast [1-3], direction [4-7], color [8, 9], etc.) in a spatially localized area called the receptive field (RF). Characterizing RGC RFs is a crucial step for understanding their functional roles in visual information processing [1, 3, 10, 11].

One of the widely used methods for estimating the RF of a spiking neuron is the spike-triggered average (STA), also called the reverse correlation [10, 12]. The STA is a linear estimate of a given neuron’s RF, is easily calculated, and has a simple and intuitive interpretation.

However, since the STA is based on a simple average, it fails to capture more complex RF structures. For instance, some RGCs respond to both ON and OFF stimulus patterns (ON-OFF cells) [12-14]. The physiological mechanism of such a bimodal response and its significance for visual information processing have been reported [5, 15-20]. Because ON and OFF stimulus patterns have similar spatiotemporal structures but opposite signs, they would cancel out in the STA analysis of an ON-OFF RGC. Thus, the STA provides little information about the actual RF of an ON-OFF RGC.

An alternative approach is based on the higher-order statistics of spike-triggered stimuli. For example, the spike-triggered covariance (STC) is the covariance matrix of spike-triggered stimuli [12]. However, the STC generally provides too much information, and therefore only a few significant STC eigenvectors are generally used to indirectly capture the structure of an RGC’s RF [12,14,21-26]. Such significant STC eigenvectors could be interpreted as additional RGC functional filters not found in the STA. Including such filters improved the accuracy of reconstruction of the presented stimuli from the RGC spikes [12, 13].

However, STC-based analyses are computationally cumbersome [12]. Another issue with STC is that these eigenvectors provide only stimulus directions, to which a given RGC’s response varies considerably. The potential problem lies in the interpretation of the eigenvector as a direction change for different RGC types. For example, for RGCs with unimodal RFs (such as ON and OFF cells), the eigenvectors provide the directions of the variability of the spike-triggered stimuli, providing a measure of RF specificity. In contrast, for RGCs with multi-modal RFs (such as ON-OFF cells), the first eigenvector indicates the direction along which the centers of ON and OFF responses are aligned. Such different interpretations of STC eigenvectors for different RGC types may hinder the adoption of STC-based analysis.

Despite the wealth of studies using STC eigenvectors, a consensus interpretation seems to be lacking [12,14,21-26]. For instance, Fairhall et al. [13] investigated the temporal structure of RGC RFs by STC analysis for rapidly changing but spatially uniform visual stimuli. They projected the spike-triggered stimuli onto a few significant eigenvectors and proposed five RGC sub-types based on the projection distribution. However, the subtype classification was based on manual inspection of spike-triggered stimuli projected onto the two eigenvectors, which is subjective and prone to classification errors. Cantrell et al. [22] investigated the development of ON and OFF retinal pathways by quantifying the bimodality of the non-centered STC, and calculated the STC without subtracting the average. They interpreted the eigenvector corresponding to the largest eigenvalue of the non-centered STC as the direction of ON or OFF responses for ON or OFF RGCs, respectively, and as the separation of ON and OFF responses for ON-OFF cells. Based on the projections of the spike-triggered stimuli onto the largest eigenvector, they calculated STC-NC bias to measure the relative response strength to ON or OFF stimulus patterns, and the STC-NC bias values were used to classify RGC subtypes (ON when STC-NC bias>0.6, OFF when STC-NC bias<-0.6, ON-OFF otherwise). However, they also admitted that this STC-NC bias was unreliable when the ON and OFF responses of an RGC are unbalanced (Fig. 4C in [22]).

Therefore, to overcome these limitations, we developed an alternative approach to analyze RGC RFs based on an unsupervised learning technique. Instead of providing a simple average such as the STA or a second-order moment such as the STC, we aimed to group spike-triggered stimuli into multiple clusters and provide information about cluster centers, especially with regard to the separation of ON- and OFF-centers in the RFs of ON-OFF cells.

The proposed method was motivated by the limitations of in vitro retinal experiments. To investigate retinal network functions, researchers typically aim to collect as many RGCs as possible from each retinal patch. However, ~70% of RGCs are identified by spike sorting as ON or OFF RGCs based on evident ON or OFF stimulus patterns, respectively, from STAs. The remaining ~30% of RGCs produce action potentials, each clearly identified by spike sorting, but the problem of determining the relationship between such spikes and the input stimuli is intractable by STA. Usually, researchers label such RGCs as unknown cells and exclude them from the analysis. Thus, in this study, we explored novel means to identify RGCs other than classification as ON and OFF types to reduce the number of unknown cells discarded in further analysis. Therefore, we developed a practical tool to analyze RGC RFs by overcoming limitations of the STA (too simple to capture complex RF structures) and STC (computationally cumbersome and hard to interpret).

Multi-electrode recordings of mouse RGCs

Retinae from male mice at postnatal day 56 and later (>P56) (C57BL/6J strain; The Jackson Laboratory, ME, USA) were used in this study (n=4). The experimental protocols were approved by the Institutional Animal Care Committee of Chungbuk National University (approval no. CBNUA-1172-18-02).

After a 20-min period of dark adaptation, visual stimuli were presented for 15 min as follows. A single visual stimulus frame consisted of 64 pixels (W=8 pixels for width and H=8 pixels for height) with a pixel width of 215 μm. Each pixel was independently turned on or off with equal probability, encoded as 0.5 and -0.5, respectively. The mean light intensity was 3.8 μW/cm². This random checkerboard stimulation was updated every 100 ms to match the highest temporal sensitivity (10 Hz) of mouse RGCs [27]. For a characterization of ON, OFF, or ON-OFF RGCs, we applied full-field illumination of 50 repeated trials of 4 s of a white screen followed by 4 s of a black screen. The mean light intensity was 40 μW/cm².

RGC spikes were recorded using a 60-channel multi-electrode array (MEA). The data acquisition system (MEA60 system; Multichannel Systems GmbH, Reutlingen, Germany) included a planar 60-channel, perforated MEA (60pMEA200/30iR), an amplifier (MEA1060), and temperature control units (TC01). The MEA contained 64 circular electrodes in an 8×8 grid with electrode diameters of 30 μm and inter-electrode distances of 200 μm. Multi-electrode recordings of retinal activity were obtained from 59 electrodes, excluding one reference electrode and four inactive electrodes, with a band pass from 1 to 3,000 Hz, a gain of 1,200, and a sampling rate of 25 kHz. This raw waveform was high-pass filtered with a 100-Hz cutoff, thresholded with a threshold four times the standard deviation of the background noise, and then spike-sorted using Offline Sorter^TM (Plexon Inc., TX, USA). Subsequent analysis was performed using custom-made Python codes and the Scikit-learn package [28].

Spike-triggered average and covariance

First, a spike-triggered stimulus was defined as the stimulus preceding a spike in a time window of T=700 ms. The time window included the boundary values (T=700 and 0 ms), resulting in eight visual stimulus frames. The time window size was chosen for the following two reasons. First, the length of STA waveforms in previous studies on mouse retinae was <1 s. Second, we calculated STAs with longer time windows and found that STA values beyond 700 ms were close to zero. Thus, the spike-triggered stimulus corresponding to a spike is a three-dimensional volume (W×H×T=8×8×8) of binary data. For notational convenience, each spike-triggered stimulus is concatenated into a 512-dimensional vector and annotated by S⁽ⁿ⁾, where n=1, 2, …, N indexes the corresponding spike and N is the total number of spikes.

The spike-triggered average (STA) is then the spike-triggered stimuli averaged over all spikes:

STA=1NΣn=1NS(n)

Intuitively, STA indicates the average stimulus pattern to which a given RGC is sensitive (blue square in Fig. 1B, inset). An STA has the same dimensionality (512) as a spike-triggered stimulus, which allows easier interpretation. A given STA was quantified using the peak-to-peak difference, defined as the difference between the maximum and minimum STA values out of the 512 (Fig. 2A, top left).

Similarly, the spike-triggered covariance (STC) is the covariance of spike-triggered stimuli:

STC=1N-1Σn=1N(S(n)-STA)(S(n)-STA) T,

where [·]^T represents the transpose of a given vector. Thus, each cell’s STC is a 512×512 matrix and conveys information about the ellipsoidal fit of spike-triggered stimuli (blue ellipse in Fig. 1B inset). STC eigendecomposition further reveals the major axes of this ellipse. The two largest eigenvalues (λ₁ and λ₂) of the STC and the corresponding eigenvectors (v₁ and v₂) were identified.

Using the STA, RGCs were manually classified into three types: ON, OFF, and unknown. Classification criteria were as follows: An RGC was classified as either ON or OFF if both temporal and spatial profiles showed a strong peak localized in time (typically about 200 ms before a spike) and space (close to the recorded electrode) with a preceding peak (typically about 300 ms before a spike) of the opposite sign. The sign of the former (closer in time to the spike) determines the cell type as ON or OFF. Otherwise, the RGC was labeled as of unknown type.

Spike-triggered clustering analysis

The proposed method, called spike-triggered clustering (STCL) analysis, comprises two stages of computations ? dimension reduction and clustering (Fig. 1B), details of which are as follows.

First, the dimension of the spike-triggered stimuli were reduced from 512 to two by projection onto the first and second STC eigenvectors (v₁ and v₂). Notationally, each S⁽ⁿ⁾ produced a dimension-reduced spike-triggered stimulus r⁽ⁿ⁾=(r₁⁽ⁿ⁾, r₂⁽ⁿ⁾), with each coordinate calculated by the projection (equivalent to the inner product):

ri(n)=S(n)·vi,

where i indexes the reduced dimension, and [·] represents the inner product of the two vectors. For visualization, scatterplots of r⁽ⁿ⁾ (Figs. 3~5A) and histograms of r₁⁽ⁿ⁾ (Figs. 3~5B) were calculated.

The dimension-reduced spike-triggered stimuli were then clustered into two groups using a Gaussian mixture model (GMM) [29]. The number of clusters is a free parameter. To identify the ON and OFF responses of ON-OFF RGCs, we set the number of clusters to two. More specifically, r⁽ⁿ⁾ for each RGC was fit to the GMM using the expectation-maximization (EM) algorithm, which considers an initial group index (g₀⁽ⁿ⁾=1 or 2) of each r⁽ⁿ⁾, iteratively updates the GMM to convergence, and produces the final estimated group index (g⁽ⁿ⁾=1 or 2) for each dimension-reduced spike-triggered stimulus.

To encourage the clustering algorithm to capture two stimulus patterns of opposite sign, we initialized the GMM’s initial group index based on the signs of the first coordinate of the dimension-reduced spike-triggered stimulus (r₁⁽ⁿ⁾):

g0(n)={2,  otherwise.1,  ifr1(n)>0,

After convergence, group indices (g⁽ⁿ⁾) were used to cluster the original spike-triggered stimuli into two clusters and calculate the two cluster centers c₁ and c₂ according to the predicted group index:

ck=Σg(n)=ks(n)Σg(n)=k1,

where k=1 and 2 index the groups. Note that each group center has the same dimension (512) and a similar interpretation as the STA, but with a determination of multiple centers suitable for analyzing RFs with complex structures.

These cluster centers were used to classify RGCs into four types: ON, OFF, ON-OFF, and unknown. If c₁ and c₂ of an RGC showed strong ON and OFF stimulus patterns with opposite signs, the RGC was classified as an ON-OFF cell. If either c₁ or c₂ of the RGC had a strong ON or OFF stimulus pattern, the RGC was classified as an ON or OFF RGC, respectively. Otherwise (no identifiable peaks in any of c₁ or c₂), the RGC was labeled as an unknown type.

Two measures were used to quantify cluster centers. First, similarly to the STA, the peak-to-peak difference of a cluster center was defined as the difference between a given group center's maximum and minimum values. Next, to quantify the arrangement of the two groups in the RF, we calculated the inner product (I) of the two cluster centers:

I=C1·C2

where [·] represents the inner product of the two vectors. The inner product would have a negative sign if the two cluster centers showed responses with opposite signs. The absolute value of the inner product is greater when the cluster centers are more dissimilar. Thus, a negative inner product with a large absolute value is a strong indicator of an ON-OFF RGC.

Criteria for choosing spike number for analysis

STC and spike-triggered clustering were performed only for RGCs with sufficient spikes (>4×stimulus dimension=2048), because the eigendecomposition becomes unstable and incorrect when sample size is less than 4 times of the data dimension. Thus, STC also requires far more spikes than does STA [12, 25]. Therefore, we used 171 RGCs with >2048 spikes to compare STA, STC, and STCL.

Availability of analysis software

Python code for STA and STC spike-triggered analysis and for the proposed clustering analysis is publicly available at https://github.com/ys7yoo/pysta2.

STA classification of RGC types

In four mouse retinal patches, a total of 389 RGCs responded to the light stimulus. Based on the STA, these RGCs were classified as ON (n=91, 23.4%), OFF (n=169, 43.4%), or unknown (n=129, 33.2%) (Fig. 2B). Consistent with a previous study [30], more OFF RGCs were collected than ON.

Fig. 2A shows the STAs of representative ON (top), OFF (middle), and unknown (bottom) RGCs, where the left and middle columns show the temporal and spatial STA profiles, respectively. The temporal STA profiles of ON RGCs (Fig. 2A, top left) showed two prominent peaks (a negative peak at 300 ms before a spike, followed by a prominent positive peak at 200 ms before a spike) consistent with the typical ON stimulus pattern. The STA of this RGC showed a large peak-to-peak difference of 0.89. These positive and negative peaks were spatially localized in the upper left corner (Fig. 2A, top middle), consistent with the location of the recorded electrode. The example OFF RGC had an STA (Fig. 2A, middle) similar to that of the example ON RGC, except for a reversed polarity (peak-to-peak difference=0.45). In contrast, the unknown RGC type did not show any noticeable peaks in its STA (Fig. 2A, bottom; peak-to-peak difference=0.09). In addition, with full-field illumination, functional types of RGCs were classified (Fig. 2A, right). Post-stimulus time histogram (PSTH) graphs indicated that the RGC types for each row were ON, OFF, and ON-OFF, respectively.

For classifying RGC types based on the STA, the peak-to-peak difference was a key feature. Fig. 2C shows the average peak-to-peak difference for each RGC type with the standard error indicated by the error bar. The peak-to-peak differences of ON and OFF RGCs were not significantly different (t-test, p>0.05), whereas those of the unknown RGCs were significantly less than those of ON and OFF RGCs (t-test, p<0.001).

STA, STC, and STCL comparisons

In the following section, cell type classification results by STA, STC, and STCL are compared for representative RGCs. The first case (Fig. 3) demonstrates an example where STA and STCL classification results are in agreement. The other two cases (Figs. 4, 5) illustrate cases where cell type classification results by STA and STCL differ, and the latter shows more detailed spatiotemporal RFs.

Fig. 3 shows an exemplar RGC classified as ON by both STA and STCL. Fig. 3A shows the scatterplot of the dimension-reduced spike-triggered stimuli (r⁽ⁿ⁾), whose mean (STA) and covariance (STC) are represented as blue squares and ellipses, respectively. Clustering analysis divided r⁽ⁿ⁾ into two groups, whose means and covariances are marked by red triangles and ellipses, respectively. One of the group centers is close to the STA with a covariance smaller than that of the STC. The other group is shifted toward the right. In Fig. 3B, the histogram of r₁⁽ⁿ⁾, which is the first dimension of r⁽ⁿ⁾, appears to be unimodally distributed, with a dominating projection from the first group identified by the STCL.

For such an ON RGC, the STA and STCL produced similar RFs. Based on the STA, this cell was classified as an ON RGC because temporal (Fig. 3C, left) and spatial (Fig. 3C, right) STA profiles show a strong ON response (peak-to-peak difference=0.62). However, STC analysis provided no further information. STC eigenvalues vary rather smoothly (Fig. 3D, left, inset), and the eigenvector corresponding to the largest eigenvalue does not show any noticeable peaks (peak-to-peak difference=0.24). In contrast, the first group center of the STCL shows a strong ON response (Fig. 3E, top), with an even greater peak-to-peak difference (0.79) than the STA (0.62). The other group center does not show any noticeable peaks (Fig. 3E, bottom; peak-to-peak difference=0.24). Thus, in this case, the STA and STCL classified the same RGC as an ON type, while the STCL further separated the stimulus-induced spikes (the first group) from the others (the second group). In addition, with full-field illumination, a functional type of the RGC was confirmed (Fig. 3F). The PSTH graph indicated that the RGC type was indeed ON.

Next, Fig. 4 shows an example RGC classified as ON by STA but as ON-OFF by STCL. Fig. 4A shows the scatterplot of the dimension-reduced spike-triggered stimuli (r⁽ⁿ⁾), which appears to include two overlapping groups. The STA (blue square) and STC (blue ellipse) considered r⁽ⁿ⁾ as a whole and failed to capture the finer structure of r⁽ⁿ⁾. In contrast, STCL correctly identified the centers of the two groups (red triangles) with smaller covariances (red ellipses) than STC. In Fig. 4B, the histogram of r₁⁽ⁿ⁾ appears to be bimodal with two overlapping modes with different strengths, consistent with the covariances’ relative sizes provided by STCL. For an RGC with an ambiguous STA, STCL recovered a more refined RF structure, while STC provided no additional information. The temporal STA profile (Fig. 4C, left) showed three peaks at 500 ms, 300 ms, and 200 ms before a spike. While the second and the third peaks agree with a typical ON response, the weaker peak at 500 ms before a spike is somewhat unusual in RGC STAs. Nevertheless, because of the high peak-to-peak difference (0.62) and localized spatial profile (Fig. 4C, right), this RGC was classified as an ON type by STA criteria alone. STC analysis failed to capture any meaningful structure in the variance. One STC eigenvalue stands out (Fig. 4D, inset), but the corresponding eigenvector does not contain any noticeable peaks (Fig. 4D). In contrast, STCL identified strong ON and weak OFF responses in both groups (Fig. 4E).

The unusual STA of this RGC is well interpreted by STCL as follows: First, the average of the two group responses identified by STCL would be a weak ON response, which agrees with the STA’s polarity. In addition, the waveform of this RGC’s OFF response (Fig. 4E, bottom), with peaks at 400 and 200 ms before a spike, changes more slowly than does that of the ON response (Fig. 4E, top), which has peaks at 300 and 200 ms before a spike, implying that this RGC responds to a typical ON-type stimulus as well as to a more slowly changing OFF-type stimulus. Because of this temporal mismatch between the ON and OFF components, the positive peak of the OFF response at around 400 ms (Figottom) is not fully canceled by the negative peak of the ON response at around 300 ms (Figop), which results in a weak positive peak at around 500 ms of the STA. Thus, this case demonstrates that some RGCs classified as ON or OFF may be ON-OFF RGCs, and that clustering analysis can reveal unbalanced bipolar properties of ON-OFF cells. PSTH analysis showed that the RGC type was ON-OFF (Fig. 4F). However, the strength of the ON response was much greater than that of the OFF response, indicating an unbalanced response.

Finally, Fig. 5 shows an example RGC classified as an unknown type by STA but as an ON-OFF type by STCL. Fig. 5A shows the scatterplot of the dimension-reduced spike-triggered stimuli (r⁽ⁿ⁾), which clearly contains two ellipsoidal blobs of approximately equal size. The blobs are separated into two groups by spike-triggered clustering analysis. Cluster centers (red triangles in Fig. 5A) and individual covariances (red ellipses in Fig. 5A) are reliably estimated despite some overlap between the two groups.

Clustering analysis identified two groups in r⁽ⁿ⁾ as indicated by their means (red triangles in Fig. 5A) and covariances (red ellipses in Fig. 5A). In contrast, the STA (blue square in Fig. 5A) is located close to the middle of the cluster centers, and to the origin. However, few samples are close to the STA, and thus the STA itself fails to capture the bimodal RF structure. The STC (blue ellipse in Fig. 5A) includes the two blobs, but its interpretation requires further analysis. In Fig. 5B, the histogram of r₁⁽ⁿ⁾ has two distinctive peaks with opposite signs, an indicator of ON-OFF?type RGCs.

Both the STA and STC failed to capture the RF of the ON-OFF RGCs. Based on STA analysis only, this cell was classified as an unknown type because of the weakness of the temporal and spatial STA profiles (peak-to-peak difference=0.14; Fig. 5C), in agreement with the previous observation in Fig. 5A that the STA, being a simple average of spike-triggered stimuli, is close to the origin. On the other hand, one STC eigenvalue (Fig. 5D, left, inset) stands out from the rest, and the corresponding eigenvector appears to be temporally (Fig. 5D, left) and spatially (Fig. 5D, right) an ON response. However, the RGC polarity (ON or OFF) cannot be determined from the STC eigenvector because of its scale ambiguity: any scalar (including a negative value) multiple of an eigenvector is, in principle, another valid eigenvector corresponding to the same eigenvalue. The only valid observation for eigenvector profiles is that spike-triggered stimuli vary so that stimuli at 300 ms and 200 ms before a spike are anti-correlated. Therefore, ON-OFF RGCs cannot be identified by STA or STC. These limitations motivated us to develop STCL to recover the actual centers and covariances of individual blobs in the RFs of ON-OFF RGCs.

In contrast to STA and STC analysis, spike-triggered clustering analysis identified the bipolar RFs of ON-OFF RGCs. The cluster centers identified by STCL correspond to strong ON and OFF stimulus patterns (Fig. 5E; peak-to-peak differences=0.77 and 0.65, respectively). Because these two responses with opposite signs vary in tandem, their average is close to nil, consistent with the flat STA shown in Fig. 5C. In addition, PSTH analysis confirmed that this RGC is an ON-OFF-type RGC (Fig. 5F). Therefore, STCL analysis reveals bipolar properties in the RFs of ON-OFF cells that cannot be obtained from the STA.

Comparison of STCL and STC-NC

The proposed method (STCL) is compared to the STC-NC method by Cantrell et al. [22] using the experimental data and analysis code released by Cantrell et al. [22], available at code.google.com/archive/p/non-centered-spike-triggered-covariance. The provided dataset contains 12 RGCs, which were classified as 4 ON, 3 OFF, and 5 ON-OFF RGCs by STC-NC based on the STC-NC bias values calculated by the accompanied Matlab code (Fig. 6A, inner pie chart). In contrast, the STCL analysis classified the 12 RGCs as 4 ON, 1 OFF, and 7 ON-OFF RGCs (Fig. 6A, outer pie chart). The STCL analysis identified two more ON-OFF RGCs, which were labeled as OFF RGCs by STC-NC.

The panel B of Fig. 6 shows a comparison of STC-NC bias values with STCL inner product values. The STC-NC bias and STCL inner product were in agreement for ON and ON-OFF RGCs. More specifically, ON RGCs had higher STC-NC bias and STCL inner product values (red dots in Fig. 6B) while ON-OFF RGCs had lower STC-NC bias and STCL inner product values (green dots in Fig. 6B), indicating bimodality of the RFs. However, two RGCs (black squares in Fig. 6B) labeled as OFF by STC-NC had low STC-NC bias values (<-0.6), indicating strong unimodality with negative sign, but low STCL inner product, indicating bimodality of the RFs.

The panel C of Fig. 6 shows the histogram of the projections to the first eigenvector of the non-centered STC, used for calculating the STC-NC bias value (-0.91), of the misclassified RGC (ch 37 unit 5). The projection histogram contained a prominent peak on the right and a much weaker peak on the left. Thus, the STC-NC bias value (-0.91) was much smaller than the threshold value (-0.6), and this RGC was classified as OFF RGC by STC-NC. This demonstrates that the STC-NC bias, based on the projection to the first eigenvector, fails to identify ON-OFF RGCs with unbalanced bimodality.

In contrast, this unbalanced bimodality is clearly shown in the two-dimensional projection in the panel D of Fig. 6. The centers of the two groups were accurately identified by the STCL centers (red triangles in Fig. 6D).

One might suspect that the less prominent group of spike-triggered stimuli, additionally identified by the STCL analysis, might be due to some random noise. To test this, we have plotted the spatio-temporal profiles of the STCL centers in Fig. 6E. Even though Cluster center 1 (Fig. 6E top) is somewhat noisier than Cluster center 2 (Fig. 6E bottom), the former has clear ON response and shows even higher peak-to-peak difference (0.22) than that of the latter (0.21). Thus, we conclude that this RGC responded to both ON and OFF stimulus patterns and should have been classified as an ON-OFF RGC.

In summary, STCL identified more ON-OFF RGCs than STC-NC did. The additionally found ON-OF RGCs had unbalanced bimodality in RFs, which STC-NC failed to capture.

Cell type classification using STA versus STCL

The cell type classifications of the 171 RGCs with sufficient spikes were compared for STA and STCL. The nested pie chart in Fig. 7A shows the cell counts of each type as classified using STA (inner pie) and STCL (outer pie). Using STA analysis, the 171 RGCs were classified as ON (n=48, 28%), OFF (n=84, 49%), and unknown (n=39, 23%) (Fig. 7A, inner pie). Using STCL analysis, the same 171 RGCs were classified as ON (n=35, 20.4%), OFF (n=79, 46.2%), ON-OFF (n=23, 13.4%), and unknown (n=34, 20%) (Fig. 7A, outer pie). Among RGCs labeled as ON (n=48) or OFF (n=84) by STA, 13 RGCs labeled as ON and five labeled as OFF were reclassified as ON-OFF cells (Fig. 7A, left, in green). On the other hand, among the 39 RGCs labeled as unknown by STA analysis, five were reclassified as ON-OFF cells (Fig. 7A, right, in green). Therefore, using STCL, a total of 23 (13 %) ON-OFF RGCs were identified among the 171 RGCs. Interestingly, more ON-OFF RGCs were found among RGCs previously labeled as ON or OFF than among the unknown types.

Characteristics of ON-OFF RGCs identified by STCL

With STCL, only the ON-OFF RGCs had peak-to-peak differences in the cluster centers significantly greater than those of STAs (Fig. 7B). For RGCs consistently labeled as ON or OFF by STA and STCL, the peak-to-peak STCL differences did not differ from those of the STA (one-sample t-test, p>0.05). In contrast, ON-OFF RGCs identified by STCL had significantly larger peak-to-peak differences in the cluster centers than did those of STAs (one-sample t-test, p<0.001). This finding indicates that the cancellation of ON and OFF responses in the STAs of ON-OFF RGCs was successfully recovered by STCL analysis.

A more negative inner product (I) between clustering centers is a strong indicator of ON-OFF RGCs. We observed no significant difference between the inner product values of ON and OFF RGCs (Fig. 7C, t-test, p>0.05). The inner products between the cluster centers were significantly less for ON-OFF RGCs than for ON or OFF RGCs (Fig. 7C, t-test, p<0.001), which indicates that the cluster centers of ON-OFF RGCs were separated further toward opposite polarities.

Accurate RGC RF mosaics with STCL

Fig. 7D shows STA-based mosaics of ON and OFF RGCs for a retinal patch. Cells for which STA and STCL agree as to type are shown as solid ellipses (ON RGCs in red, left; OFF RGCs in blue, right). The RFs of ON and OFF RGCs identified as ON-OFF by STCL are shown as dashed ellipses in green. RF ellipses were mainly harvested in the top left but not in the bottom right region might be due to the weak attachment between MEA and retinal patch in the bottom right region.

These mosaics are consistent with previous studies where ON and OFF RGCs tile the retina [31, 32]. However, some ON and OFF RFs (dashed ellipses in green) were incorrectly classified here by conventional analysis, and information about ON-OFF RGCs in previous studies is lacking.

In contrast, STCL provided more accurate mosaics that included ON-OFF RGCs. Fig. 7E shows mosaics of ON, OFF, and ON-OFF RGCs identified by STCL for the same retinal patch shown in Fig. 7D. The RFs of RGCs identified as ON (1st, red) and OFF (2nd, blue) are shown as ellipses. In the third panel, the RFs of ON-OFF RGCs are shown, with ON responses in red and OFF responses in blue. The RFs of ON-OFF RGCs (Fig. 7E, 3rd) overlap with those of ON and OFF RGCs (Fig. 7E, 1st and 2nd). The ON and OFF RFs of each ON-OFF RGC significantly overlapped. However, for cell 61a ON-OFF RGCs, the ellipsoidal fits of ON and OFF responses of the RGC appear to be oriented toward different angles. Therefore, we double checked the spatial profile using STA and STCL. The spatial profiles of the STA (left) and STCL centers (middle and right) of the RGC in inset figure E. Solely based on the STA, this RGC appeared to be an ON RGC with a weak light modulation and RF size was about one pixel. However, the STCL analysis showed that this RGC was actually ON-OFF RGC with stronger ON and slightly weaker OFF responses. Although the ellipsoidal fits of ON and OFF responses of the RGC appear to be oriented toward different angles, the RFs of STCL centers overlap, which also coincide with the RF derived from STA. To confirm the overlapping, we also have measured the offset of centers in the ellipses between ON and OFF responses for each ON-OFF RGC (Fig. 7E, 4th). Both x- and y-offsets are not significantly different from 0 (one-sample t-test, p>0.05).

Analysis of RGC RFs: beyond the STA and STC

We first classified RGCs as to cell type based on the STA. For this, the peak-to-peak difference was a powerful feature. Specifically, RGCs with large peak-to-peak differences were unambiguously identified as either ON or OFF types. However, those with smaller peak-to-peak differences were not clearly determined. Consequently, RGCs classified as unknown type accounted for about 30% of the total RGCs recorded (129 of 389).

Spike-triggered clustering analysis overcomes this limitation of the STA and reveals multi-modal RFs. By definition, the STA provides information about the centers of spike-triggered stimulus distributions, which works well for an RGC with a simple unimodal structure. However, the STA fails to capture more complex RFs with multiple modes or groups. As a result, ON-OFF RGCs are undetected by STAs and are classified as unknown types or misclassified as ON or OFF types. Spike-triggered clustering analysis correctly identifies ON-OFF RGCs, which otherwise may be incorrectly classified as ON or OFF by simple STA analysis or left unknown. Therefore, with clustering analysis, it is possible to include ON-OFF RGCs and thereby obtain a more holistic understanding of the retinal network's information processing.

Compared to STC analysis, spike-triggered clustering analysis is more intuitive and bias-free. STC captures the second-order RF moment. In contrast to STC, the proposed method provides the actual centers of RF clusters. For ON-OFF RFs, spike-triggered clustering analysis separated spike-triggered stimuli into two groups, which correspond to typical ON and OFF stimulus patterns. Each group center has the same interpretation as in STA and is thus easier to interpret than the STC. If needed, the covariance of each group center may also be used for STC-based analysis. We also provide a practical guideline for using spike-triggered clustering analysis for analyzing RGC RFs. First, the number of spikes should be at least 4 times the stimulus dimension. For RGCs with few spikes, the RF should be analyzed only with the STA. Second, clustering analysis is recommended for apparent STA ON or OFF cells, even with large peak-to-peak differences, because such RGCs may have responded to both ON and OFF stimulus patterns with different sensitivities, and in the STA, the weaker response may have been masked by the stronger response. Third, an increase in the peak-to-peak differences going from STA to STCL and a negative inner product of cluster centers are strong indicators of ON-OFF RGCs.

Spike-triggered clustering analysis may be tailored to the demands of specific experiments by substituting other computational algorithms into each step of the clustering analysis. The first step of clustering analysis is dimension reduction. In this study, we used the STC eigenvectors for easier comparison with STC, but any other dimension-reduction algorithm would be suitable. For instance, we are currently exploring dimension-reduction algorithms that are robust to noise and outliers, which would improve the clustering accuracy and stability. Second, numerous clustering algorithms have been developed for different requirements [33-36]. We here chose the GMM, which is ideal for efficiently finding clusters with spherical shapes. We believe that this simple algorithm should work well for analyzing the RFs of most RGCs. For investigating RGCs with more complex non-spherical components, other clustering algorithms may be needed. However, it should be noted that more flexible and advanced clustering algorithms would require more spikes from each RGC and may suffer from overfitting problems [37].

In this study, we intentionally fixed the number of clusters at two for identifying ON-OFF cells; the focus was on recovering ON and OFF responses from RGCs that could not be detected by STA. Clustering analysis could be used to investigate RFs at finer scales. Recently, RGC RF subunits were shown to co-localize with the RFs of simultaneously recorded bipolar cells [14]. Therefore, RGC RF subunits may provide information about presynaptic bipolar-cell inputs and thus the functional connectivity of the retinal network. The number of clusters could be increased to identify further RGC RF subunits.

Physiological implications of ON-OFF RGCs

ON-OFF RGCs that respond to both light-on and light-off are known as bistratified cells [5, 38, 39]. Their dendrites ramify in both sublamina a and b of the inner plexiform layer, where one arborizes with the axon terminations of ON bipolar cells and the other with the axon terminations of OFF bipolar cells. Further, ON-OFF RGCs comprise 30% of the RGC population in the mouse retina [40]. In the mouse, their axons project mostly to the dorsal lateral geniculate nucleus (dLGN), a specialized relay that carries information about stimulus motion and eye movement [41]. In contrast, in the primate retina, only 10% of RGCs are ON-OFF type [42, 43]. Their destination is the LGN’s koniocellular layer, involved in a variety of visual functions, including color vision, eye movement, and motion detection [44, 45].

In visual processing, ON-OFF RGCs are known to serve four main functions. The first is movement detection: frog ON-OFF cells are considered fly detectors or moving-edge detectors, primarily based on their high sensitivity to motion [46]. Some pigeon ON-OFF cells are classified as motion-sensing units [47]. Moreover, in most functional studies, ON-OFF cells are known as direction-selective RGCs, sensitive to the direction of moving objects [48-50]. Second, some ON-OFF RGCs have been described as local edge detectors [6, 17, 51]. They are non-direction selective and are involved in encoding edge positions. The third is encoding conditional stimulus sequences. ON-OFF cells have been reported to respond to specific sequential changes. For instance, Sakai et al. [52] reported that spike firing is enhanced by specific green/red sequential changes in the gourami retina, while Uchiyama et al. [20] showed that sequential OFF/ON changes in light intensity were encoded by spike pairs with a 20-ms interspike interval. In an in vitro study using salamander RGCs, Geffen et al. [16] showed that some ON-OFF RGCs, usually dominant for OFF responses, change RF polarity from OFF to ON when an image shift occurs in the surrounding region of the RGC RF. These results suggest that ON-OFF cells go beyond motion detection to encode more specific sequential changes. The fourth is color coding. Particularly in the primate retina, small bistratified cells (SBCs), a subgroup of ON-OFF RGCs, are involved in color vision, always responding ON to blue and OFF to yellow [19, 53, 54]. SBCs exhibit two sub-receptive fields, with an S-ON component (S-cone ON input dominant) and an (L+M)-OFF component (L and M-cone OFF input dominant). As such, ON-OFF RGCs play a variety of functional roles in the retina, an early processing stage of the visual system.

For this reason, it is necessary to reliably classify ON-OFF RGCs for advanced studies of visual processing in the retina. In this regard, STCL provides a useful tool as an RGC classifier.

In our current study, we showed diversity in ON-OFF cells depending on symmetry (balance) of bimodality. RGCs classified as ON (Fig. 4) or unknown (Fig. 5) turned out to be ON-OFF RGCs by clustering analysis.

The ON-OFF cell responds to both light-on and light-off, originating from both inputs of on-bipolar cell and off-bipolar cell. Therefore, the response bias of ON-OFF is determined by the response preference for the bipolar cell input. This bias in the ON-OFF cell has been reported in several studies. Tikidji-Hamburyan et al. [55] showed the on-off response of mouse ON-OFF RGCs varied with light intensity. As mentioned already, Geffen et al. [16] showed that some ON-OFF RGCs, usually dominant for OFF responses, change RF polarity from OFF to ON in response to peripheral image shift. These studies suggest that the response preference of ON-OFF cells can be altered depending on the visual environment.

Thus, STCL analysis can be applied to reveal the response preference of ON-OFF RGCs, where it is possible to investigate the dominant input of on- and off-bipolar cells for specific visual stimuli. If the two opposite bipolar inputs are the same, the RGC can be classified as a balanced ON-OFF cell, whereas, if it is asymmetric, the RGC can be an unbalanced ON-OFF cell.

Further investigation of functional roles of ON-OFF RGCs would be another direction for future work. Information analysis in retinal networks have been studied usually including ON and OFF RGC only [24, 56]. Including ON-OFF RGCs in the retinal population pool could provide more accurate analysis.

Limitations and future direction

This study has a few limitations. First, since we did not use a moving-bar stimulus in our experiments to properly isolate direction-selective RGCs (DS RGCs), there is a possibility that ON-OFF RGCs may include DS RGCs. In future studies, we plan to include moving-bar stimuli to discriminate DS RGCs.

Second, different RGC types have different RF sizes (50~400 μm) in the mouse retina [40, 57]. For instance, it is known that the RF size of transient cells is much larger than the total average RGC (RF diameter ~200 μm), while sustained cells have RF size similar to or slightly smaller than average. Since we used a 215-μm pixel size for the stimulus, there is a possibility of misclassification for RGCs with smaller RF size, including sustained RGCs. With higher spatiotemporal resolution of the visual stimuli, spike-triggered clustering analysis could be implemented to identify and quantify the RGC RFs in fine scale. Our proposed clustering method could be also implemented to more complex receptive fields in higher visual pathways and to cell-type classification in general.

This research was supported in part by the Bio & Medical Technology Development Program (NRF-2017M3A9E2056460), the Basic Science Research Program (NRF-2018R1A2B6003917) of the National Research Foundation of Korea (NRF) funded by the Ministry of Science and ICT (MSIP), and an Incheon National University research grant (2017).